I have data from patients treated with 2 different kinds of treatments during surgery.
I need to analyze its effect on heart rate.
The heart rate measurement is taken every 15 minutes.
Given that the surgery length can be different for each patient, each patient can have between 7 and 10 heart rate measurements.
So an unbalanced design should be used.
I’m doing my analysis using R. And have been using the ez package to do repeated measure mixed effect ANOVA. But I do not know how to analyse unbalanced data. Can anyone help?
Suggestions on how to analyze the data are also welcomed.
As suggested, I fitted the data using the
lmerfunction and found that the best model is:
heart.rate~ time + treatment + (1|id) + (0+time|id) + (0+treatment|time)
with the following result:
Random effects: Groups Name Variance Std.Dev. Corr id time 0.00037139 0.019271 id (Intercept) 9.77814104 3.127002 time treat0 0.09981062 0.315928 treat1 1.82667634 1.351546 -0.504 Residual 2.70163305 1.643665 Number of obs: 378, groups: subj, 60; time, 9 Fixed effects: Estimate Std. Error t value (Intercept) 72.786396 0.649285 112.10 time 0.040714 0.005378 7.57 treat1 2.209312 1.040471 2.12 Correlation of Fixed Effects: (Intr) time time -0.302 treat1 -0.575 -0.121
Now I’m lost at interpreting the result.
Am I right in concluding that the two treatments differed in affecting heart rate? What does the correlation of -504 between treat0 and treat1 means?
The lme/lmer functions from the nlme/lme4 packages are able to deal with unbalanced designs. You should make sure that time is a numeric variable. You would also probably want to test for different types of curves as well. The code will look something like this:
library(lme4) #plot data with a plot per person including a regression line for each xyplot(heart.rate ~ time|id, groups=treatment, type= c("p", "r"), data=heart) #Mixed effects modelling #variation in intercept by participant lmera.1 <- lmer(heart.rate ~ treatment * time + (1|id), data=heart) #variation in intercept and slope without correlation between the two lmera.2 <- lmer(heart.rate ~ treatment * time + (1|id) + (0+time|id), data=heart) #As lmera.1 but with correlation between slope and intercept lmera.3 <- lmer(heart.rate ~ treatment * time + (1+time|id), data=heart) #Determine which random effects structure fits the data best anova(lmera.1, lmera.2, lmera.3)
To get quadratic models use the formula “heart.rate ~ treatment * time * I(time^2) + (random effects)”.
In this case where treatment is a between-subjects factor, I would stick with the model specifications above. I don’t think the term (0+treatment|time) is one that you want included in the model, to me it doesn’t make sense in this instance to treat time as a random-effects grouping variable.
But to answer your question of “what does the correlation -0.504 mean between treat0 and treat1” this is the correlation coefficient between the two treatments where each time grouping is one pair of values. This makes more sense if id is the grouping factor and treatment is a within-subjects variable. Then you have an estimate of the correlation between the intercepts of the two conditions.
Before making any conclusions about the model, refit it with lmera.2 and include REML=F. Then load the “languageR” package and run:
Then you can get p-values, but by the looks of it, there is probably a significant effect of time and a significant effect of treatment.