In phenological applications observations should be the julian day
of observation of a certain phase, factor 1 should be the observation year
and factor 2 should be a station-id.
Usually this is much easier created by:
y <- factor(f1),
s <- factor(f2),
ddm <- as.matrix.csr(model.matrix(~ y + s -1, contrasts=list(s=("contr.sum"))))
.
However, this procedure can be quite memory demanding and might exceed storage
capacity for large problems.
This procedure here is much less memory comsuming.
Moreover, in order to get direct estimates for all coefficients,
an additional row is appended to the matrix, where the columns for the second factor are set to 1.
Therefore, dimensions of ddm
are (nlevels(factor1)+1)x(nlevels(factor2)).