Brian, many thanks for these helpful pointers:
On 23 Feb 2005, at 15:45, Prof Brian Ripley wrote:
MASS4 p.150 White Book p.38
Those are the only two reasonably comprehensive accounts that I am aware of (and they have only partial overlap).
The underlying motivation is to span the _additional_ vector space covered by the term, the complement to what has gone before. Put another way, as each term is added, only enough columns are added to the model matrix to span the same space as if dummy coding had been used for that term and its predecessors. So think of this as a way to produce a parsimonious (usually full-rank) basis for the model space.
Yes, indeed. My surprise was that *this* particular basis (dummy coding) was the one used here. I should have got a clue from what contrasts() does, and is documented to do,
> options("contrasts")
$contrasts
[1] "contr.treatment" "contr.poly"
> contrasts(a)
.L .Q
[1,] -7.071068e-01 0.4082483
[2,] -9.073800e-17 -0.8164966
[3,] 7.071068e-01 0.4082483but when there's no intercept in the model the contrasts used appear to be
> contrasts(a, contrasts = FALSE) -1 0 1 -1 1 0 0 0 0 1 0 1 0 0 1
which are not the same as
> contr.poly(a, contrasts = FALSE)
^0 ^1 ^2
[1,] 1 -7.071068e-01 0.4082483
[2,] 1 -9.073800e-17 -0.8164966
[3,] 1 7.071068e-01 0.4082483-- which is what I had naively expected to get in my model matrix.
This is of course all a matter of convention. The present convention does seem a touch confusing though: the basis for the space spanned by a factor is determined by options("contrasts") or by the contrasts attribute of the factor or by the contrasts argument in the call, *except* when there's no intercept or other factor earlier in the model in which case all such settings are ignored (well, not *quite* ignored: they do get put in the contrasts attribute of the resultant model matrix).
On the other hand, one good reason to use dummy coding for the first-encountered factor when there's no intercept is that the associated parameters are then often interpretable as group-specific intercepts.
Would it be an improvement, though, if the "contrasts" attribute of the resultant model matrix contained "contr.treatment" in such cases instead of the name of a contrast function that was not actually used?
Best wishes, David
On Wed, 23 Feb 2005, David Firth wrote:
I was surprised by this (in R 2.0.1):
a <- ordered(-1:1) a[1] -1 0 1 Levels: -1 < 0 < 1
model.matrix(~ a)(Intercept) a.L a.Q 1 1 -7.071068e-01 0.4082483 2 1 -9.073800e-17 -0.8164966 3 1 7.071068e-01 0.4082483 attr(,"assign") [1] 0 1 1 attr(,"contrasts") attr(,"contrasts")$a [1] "contr.poly"
model.matrix(~ -1 + a)a-1 a0 a1 1 1 0 0 2 0 1 0 3 0 0 1 attr(,"assign") [1] 1 1 1 attr(,"contrasts") attr(,"contrasts")$a [1] "contr.poly"
Without the intercept, treatment contrasts seem to have been used (this despite the "contr.poly" in the "contrasts" attribute).
It's not restricted to ordered factors. For example, if Helmert contrasts are used for nominal factors, the same sort of thing happens.
I suppose it is a deliberate feature (perhaps to protect the user from accidentally fitting models that make no sense? or maybe some better reason?) -- is it explained somewhere?
-- Brian D. Ripley, [EMAIL PROTECTED] Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/ University of Oxford, Tel: +44 1865 272861 (self) 1 South Parks Road, +44 1865 272866 (PA) Oxford OX1 3TG, UK Fax: +44 1865 272595
______________________________________________ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
