[R] analysing mixed effects/poisson/correlated data

From: Alexandra Bremner <alexandra.bremner_at_uwa.edu.au>
Date: Sat, 08 Mar 2008 17:57:02 +0900


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_at_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.
Received on Sat 08 Mar 2008 - 09:13:27 GMT

Archive maintained by Robert King, hosted by the discipline of statistics at the University of Newcastle, Australia.
Archive generated by hypermail 2.2.0, at Sat 08 Mar 2008 - 14:30:19 GMT.

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

list of date sections of archive