I'm examining trends (between 1998 and 2011) in mortality rates among patients with Crohn's disease. Each patient (case) have been included during 1998 to 2011. At inclusion, each patient have been matched to a healthy control with the same age and sex. I'm analysing trends in mortality rates. When doing this directly, without any adjustments, I obtain fluctuating mortality rates over time, which is presumably due to the fact that individuals included a given year wont be comparable to those included another year. Thus I aim to adjust the mortality rates. I expect that mortality rates in both groups (cases and controls) will decline over time and the gap between cases and controls will narrow successively.
My idea is to do the adjustment by means of Poisson regression. My data is on individual level. I wish to obtain one estimate on the incidence rate (per 1000 person-years) for cases and controls each year from 1998 to 2011. Survival time would be included as the offset in the model. Something similar has been done here.
I've attached the 200 first rows from my data set, which consists of 1500 individuals. Here is the data. Variable explanation:
- dead = if patient died or not during follow-up
- surv = survival time in days
- agegroup = categorized age group (4 groups)
- gender = male/female
- diagnosis = 0 for healthy control, 1 for Crohns disease
- age = age in years
- inclusion_year = year of inclusion in the study
What did I try so far? I've tried to fit Poisson models with the glm() function in R, using individual observations (log(surv) as the offset), but I either received an error or could not figure out how to use the fits. I've also aggregated the data into groups, and then analysed the death counts in glm(); when I used the fit to obtain incidence rates I could only obtain rates for a specific age/agegroup and sex (as needed to be specified in the predict() function).
I'd really appreciate some statistical advice and coding examples, which can be done on the attached data set.
contrasts<-
(*tmp*
, value = contr.funs[1 + isOF[nn]]) : contrasts can be applied only to factors with 2 or more levels $\endgroup$diagnosis*inclusion_year
interaction terms. If you just use the current model, the case num will only differ by the beta ofdiagnosis
, constant along years because it's not allowed to interact. Afterwards, the predicts will be just substitution. I am not too picky so I'd just subs the mean age and mean male percentage. $\endgroup$