# Is this longitudinal data too complicated for GLMM or GEE?

After writing this post, I've realized that I am running around in circles, chasing my tail. Any help approaching this problem would be greatly appreciated, as I think I just need to bounce ideas around and don't have colleagues that can help with statistical methods at my new position.

I am working on an incomplete longitudinal/repeated measure dataset, which is very new to me. After spending a week working with the data, I have come to the (possibly incorrect) conclusion that I should be approaching the data with a GLMM or GEE. I say incorrect, because I wonder if I have too many dimensions to model through these approaches (potentially leading to a PCA to reduce my dimensions?). I've been attempting to approach the data with lme4() in Program R, but am having a major issue in wrapping my head around how to account for the various variables of interest. I believe my data has a layer or two of complexity beyond the 'sleepstudy' data in Bates Chapter 3 (on the use of Mixed Effects Models in R), so would appreciate any advice on what I am missing and what direction I may pursue, as my brain has reached it's carrying capacity for selflearning at the moment

My simulated dataset can be found here. A treatment (drug vs control) was administered to a set of subjects at TimePeriod 1. A Measurement of interest was taken at TimePeriod 0 (before treatment), and at TimePeriods 1, 2, 3, 4, 5 (spaced each 12hrs). These treatments are conducted on newborn babies born at different gestational ages (GA), but are conducted at similar (but not identical) HOL (hour of life). The dataset includes:

• 'GA' (gestational age at birth)
• the biological 'Measurement' of interest
• 'TimePeriod' (which timeperiod the measurement was taken at. 0=predrug/control administration, 1=post, 2=12hrs after, 3=24hrs after, etc)
• HOL (hour of life at eachmeasurement)
• 'Group' (term vs preterm birth, which is a categorization of GA that I don't think is needed)
• 'Treatment'(drug/control).

There are 6 measurements per subject each taken at a similar HOL between subjects. The main question is if treatment (drug/control) affects the measurement of interest. Another question is how GA and HOL relate to the measurement of interest.

We can open this in R via

newData<-"[paste linked-data here]"


We then can use ggplot2 to visualize the data, using my 4-dimensions of interest (GA, HOL, Treatment, and Measurement)

library(ggplot2)
library(grid)
s <- ggplot(Data, aes(x=HOL,y=GA, shape=Treatment,fill=Subject,color=Subject))
s <- s+geom_point(aes(size=Measurement))
s + opts(legend.key.height=unit(.3, "cm"))


. What we see with this plot is GA (gestational age) on the y-axis, and HOL (hour of life on the xaxis). The size of the points give information on the Measurement of interest, while shape gives Treatment group. Clearly, the Measurement increases with increasing GA and HOL, and each subject has a unique GA with a ~15week spread of values, while the 6 measurements per individual are at similar HOLs. It is important to note that this is based on Biological Relevance. In our real data, we expect a strong effect of increasing HOL and GA on Measurement

1) The first question to analyze. Does Treatment (control/drug) affect the measurement of interest vs 'TimePeriod' (which tracks time of dosing, where solid black line is the dose/treatment. X-axis has been 'jittered' for clarity).

I envisioned a GLMM where TimePeriod is X, Measurement is Y (graph above), and then HOL and GA are the Random Effects. But then, how do I test for a Treatment effect between Drug/Control? I have created a plot of mean values and confidence intervals (normalized for within-subject variability via Morey 2008), which shows no effect of Treatment, but I wonder how to approach this with a GLMM (if that's appropriate), instead of just visualizing group means and CIs. I believe I understand how to make a GLMM for the Drug or Control group separately, but do not see how to compare the two through a statistical test (beyond plotting them on the same graph). I am clearly overlooking something in my brain-fog.

2) The first question is really my main interest. The second question is then, assuming there is no Treatment effect, can I model the relative importance of HOL and GA on the measurement of interest (again, we believe we will find no Treatment effect, but will find a HOL and GA effect)? This relationship is seen in the first graph I provide (measurement increases as you move up on both axises), but for some reason, I am having a complete mental block on approaching this question statistically.

EDIT (6/26/14): I've been working through this, and have come to this approach.

 lmer.full<-lmer(log(Measurement)~HOL*Treatment+GA+(1|Subject))


This allows a treatment of longitudinal data through the (1|Subject) random effect. It also simplifies the data since 'TimePeriod' and 'HOL' are so heavily correlated. Instead of asking how the Measurement changes per time-period, I can explore the interaction between Hour-of-Life and Treatment.

Here is my new issue. Gestational Age is heavily correlated with the Subject effect (which makes sense since each subject has it's own Gestational Age). My questionis, if due to this high correlation, I should leave out GA, and just focus on the (1|Subject) random effect to control for this variability.

EDIT 2: I have improved the model on advice of Bodo Winter, who has a great personal website with some nice beginner tutorials I highly recommend www.bodo-winter.net

First, I zero-centered my Time parameter

 Data$HOL<-Data$HOL-mean(Data\$HOL,na.rm=TRUE)


I also removed all initial measurements before treatment administration

 Data2<-subset(Data,Data[,6]>0)


This removed convergence issues I was having when adding a correlated random slope/intercept as such

 lmer.full.slope<-lmer(log(Measurement)~HOL*Treatment+GA+(1+HOL|Subject))

• If you're still working on this thing, you might have some luck with a Bayesian approach in Stan. Commented Jun 5, 2015 at 15:10