I am attempting to model data with the following variables:
timepoint - n=48, monthly over 4 years hospital - n=3 opsn1 - no of outcomes total.patients skillmixpc - skill mix percentage nurse.hours.per.day Aims To determine if skillmix affects rate (i.e. no.of.outcomes/total.patients). To determine if nurse.hours.per.day affects rate. To determine if rates vary between hospitals. There is first order autoregression in the data. I have attempted to use the lmer function (and lmer2) with correlation structure and weights: test1 <-lmer(opsn1~timepoint+as.factor(hospital)+ skillmixpc + nursehrsperpatday +(timepoint|hospital) +offset(log(totalpats)),family=poisson, data=opsn.totals) test2 <-lmer(opsn1~timepoint+as.factor(hospital)+ skillmixpc + nursehrsperpatday +(timepoint|hospital)+offset(log(totalpats)),family=poisson, data=opsn.totals, correlation=corAR1(form=~1|hospital)) test3 <-lmer(opsn1~timepoint+as.factor(hospital)+ skillmixpc + nursehrsperpatday +(timepoint|hospital)+offset(log(totalpats)),family=poisson, data=opsn.totals, correlation=corAR1(form=~1|hospital),weights=varIdent(form=~1|hospital)) Test1 & test2 give the same output (below). Does this mean that the correlation structure is not being used? Test3 produces the following error message (I notice there are others who have had problems with weights). Error in model.frame(formula, rownames, variables, varnames, extras, extranames, : variable lengths differ (found for '(weights)') > summary(test1) Generalized linear mixed model fit using Laplace Formula: opsn1 ~ timepoint + as.factor(hospital) + skillmixpc + nursehrsperpatday + (timepoint | hospital) + offset(log(totalpats)) Data: opsn.totals Family: poisson(log link) AIC BIC logLik deviance 196.2 223.0 -89.12 178.2 Random effects: Groups Name Variance Std.Dev. Corr hospital (Intercept) 3.9993e-03 6.3240e-02 timepoint 5.0000e-10 2.2361e-05 0.000 number of obs: 144, groups: hospital, 3 Estimated scale (compare to 1 ) 1.057574 Fixed effects: Estimate Std. Error z value Pr(>|z|) (Intercept) -2.784857 1.437283 -1.9376 0.0527 . timepoint -0.002806 0.002709 -1.0358 0.3003 as.factor(hospital)2 -0.030277 0.120896 -0.2504 0.8022 as.factor(hospital)3 -0.349763 0.171950 -2.0341 0.0419 * skillmixpc -0.026856 0.015671 -1.7137 0.0866 . nursehrsperpatday -0.025376 0.043556 -0.5826 0.5602 --- Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 Correlation of Fixed Effects: (Intr) timpnt as.()2 as.()3 skllmx timepoint -0.606 as.fctr(h)2 -0.382 0.132 as.fctr(h)3 -0.734 0.453 0.526 skillmixpc -0.996 0.596 0.335 0.715 nrshrsprptd 0.001 -0.297 0.204 -0.130 -0.056 Can anyone suggest another way that I can do this analysis please? Or a work around for this method? Any help will be gratefully received as I have to model 48 different opsn outcomes. Alex ______________________________________________ R-help@r-project.org 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.