Re: [R] Parameterization puzzle

From: Berwin A Turlach <berwin_at_maths.uwa.edu.au>
Date: Fri 21 Jul 2006 - 18:16:32 EST

G'day all,

>>>>> "BDR" == Prof Brian Ripley <ripley@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


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 and provide commented, minimal, self-contained, reproducible code. Received on Fri Jul 21 18:23:21 2006

Archive maintained by Robert King, hosted by the discipline of statistics at the University of Newcastle, Australia.
Archive generated by hypermail 2.1.8, at Fri 21 Jul 2006 - 20:16:50 EST.

Mailing list information is available at https://stat.ethz.ch/mailman/listinfo/r-help. Please read the posting guide before posting to the list.