[R] Parameterization puzzle

Berwin A Turlach berwin at maths.uwa.edu.au
Fri Jul 21 10:16:32 CEST 2006


G'day all,

>>>>> "BDR" == Prof Brian Ripley <ripley at stats.ox.ac.uk> writes:

    BDR> R does not know that poly(age,2) and poly(age,1) are linearly
    BDR> dependent.
Indeed, I also thought that this is the reason of the problem.

    BDR> (And indeed they only are for some functions 'poly'.)
I am surprised about this.  Should probably read the help page of
'poly' once more and more carefully.

    BDR> I cannot reproduce your example ('l' is missing), [...]
My guess is that 'l' is 'pyears'.  At least, I worked under that
assumption.  

Interestingly, on my machine (using R 2.3.1, 2.2.1 and 2.1.1) I cannot
fit any of the Poisson GLM that Murray tried.  I always get the error
message:

Error: no valid set of coefficients has been found: please supply starting values

But I have to investigate this further.  I can fit binomial models
that give me similar answers.

    BDR> [...] but perhaps
    BDR> glm(deaths ~ poly(age,2) + poly(age,1)*Smoke + offset(l),
    BDR> poisson)
    BDR> was your intention?
In this parameterisation a 'poly(age,1)' term will appear among the
coefficients with an estimated value of NA since it is aliased with
'poly(age, 2)1'.  So I don't believe that this was Murray's intention.

The only suggestion I can come up with is:

> summary(glm(cbind(deaths, l-deaths) ~ age*Smoke+I(age^2), family=binomial))

[...]

Coefficients:
              Estimate Std. Error z value Pr(>|z|)    
(Intercept)  -10.79895    0.45149 -23.918  < 2e-16 ***
age            2.37892    0.20877  11.395  < 2e-16 ***
SmokeYes       1.44573    0.37347   3.871 0.000108 ***
I(age^2)      -0.19706    0.02749  -7.168  7.6e-13 ***
age:SmokeYes  -0.30850    0.09756  -3.162 0.001566 ** 

[...]

Which doesn't use orthogonal polynomials anymore.  But I don't see how
you can fit the model that Murray want to fit using orthogonal
polynomials given the way R's model language operates.

So I guess the Poisson GLM that Murray wants to fit is:

        glm(deaths~ age*Smoke+I(age^2)+offset(l), family=poisson)

Cheers,

        Berwin

========================== Full address ============================
Berwin A Turlach                      Tel.: +61 (8) 6488 3338 (secr)   
School of Mathematics and Statistics        +61 (8) 6488 3383 (self)      
The University of Western Australia   FAX : +61 (8) 6488 1028
35 Stirling Highway                   
Crawley WA 6009                e-mail: berwin at maths.uwa.edu.au
Australia                        http://www.maths.uwa.edu.au/~berwin



More information about the R-help mailing list