[R] APC Modelling and the GLM function

Sue_Paul@moh.govt.nz Sue_Paul at moh.govt.nz
Mon Mar 24 02:14:22 CET 2003


Hi all
Apologies for any cross posting.
I have encountered a rather bizarre "problem" in Splus and R.  I am using Age-Period-Cohort models to model cervical cancer and have run the same data
on both R (v.1.4.1 & v1.6.2) and Splus (version 6.0).  I used the same command line in both Splus and R:  glm(cases~-1+as.factor(age)
+as.factor(period)+as.factor(cohort)+offset(log(person.years)),family="poisson",data=mydata)
While Splus and R fit APC models using different constraints, the fitted values should be identical.  However, I have found the following:
   Both Splus and R models give you the same value for the residual deviance;
   If I use the function fitted.values on the glm object then both Splus AND R models returns the same number of cases (and hence the same incidence
   rate once you have divided by person years at risk);
   However, if I try to derive the fitted values "manually":  i.e. fitted incidence rate = exp{ age.effect+period.effect+cohort.effect} then I get a
   completely different set of fitted incidence rates.
To do a quick check I also looked at second differences to see if these were identifiable, and found that the second differences for the age effects
are consistent in both R and Splus.  The period and cohort effects however, yield completely different second differences (in R & Splus).  I guess
this kind of narrows down the problem to the age and period effects, although I still cannot understand why glm would return the same deviance and
fitted number of cases, if all the second differences and fitted rates were not identifiable.
I am quite puzzled by this and can't seem to figure out what is going wrong.
I would really appreciate any help that anyone can give me.
Thanking you in advance

Kind Regards
Sue Paul
Advisor (Statistics)
Public Health Intelligence
Ministry of Health
DDI: 04 460 4926
Mobile: 021 100 3340
Fax: 04 495 4401

http://www.moh.govt.nz/PHI
mailto:sue_paul at moh.govt.nz



More information about the R-help mailing list