Averaging beta and slope from several regressions, spurious regression slopes by group

Question

This question is partly a question from ignorance of the statistics, and part an R question. I've moved this question over from Stack Overflow since it's more about the theory.

I have growth data for multiple animals, in different habitats, across a time-span of their life. I have fitted linear regression models to each one, and as expected, there is quite a bit of variation between their growth trajectory over time. Some grow faster through time, some grow slower. What I want to calculate is the average growth trajectory among animals in each location, so that I can use that information to inform a larger model...essentially, I care about the slope of their growth/day as a predictor for a larger model. I'm not all that worried about the intercept.

My data has much more variation than the example dataset below. But, like this example, there is a mix of negative and positive slopes and variation in the length and starting points of the growth data for each individual. The data was opportunistically available in each location, so it's not like I'm working from a controlled experiment.

I'm re-writing the question with a different example dataset to better get at what I'm trying to determine, which is:

"can I just fit a regression to each animal, and average the slope coefficients within each location? Or is there a more statistically valid way to do this?"

---

dput(example)
structure(list(Day = c(1, 2, 3, 1, 2, 3, 3, 4, 5, 1, 2, 3, 3.5, 
4.5, 5.5, 1, 2, 3, 1, 2, 3, 1, 2, 3, 4, 5, 6), Animal = structure(c(1L, 
1L, 1L, 2L, 2L, 2L, 3L, 3L, 3L, 4L, 4L, 4L, 5L, 5L, 5L, 6L, 6L, 
6L, 7L, 7L, 7L, 8L, 8L, 8L, 9L, 9L, 9L), .Label = c("a", "b", 
"c", "d", "e", "f", "g", "h", "i"), class = "factor"), Growth = c(5, 
4, 3, 6, 5, 4, 7, 9, 11, 1, 2, 3, 3, 4, 5, 6, 7, 8, 3, 2, 1, 
5, 4, 3, 7, 6, 5), Location = c("X", "X", "X", "X", "X", "X", 
"X", "X", "X", "Y", "Y", "Y", "Y", "Y", "Y", "Y", "Y", "Y", "Z", 
"Z", "Z", "Z", "Z", "Z", "Z", "Z", "Z"), Var1 = c(NA, NA, NA, 
NA, NA, NA, NA, NA, NA, NA, NA, NA, 3.5, 4.5, 5.5, NA, NA, NA, 
NA, NA, NA, NA, NA, NA, NA, NA, NA)), row.names = c(NA, -27L), .Names = c("Day", 
"Animal", "Growth", "Location", "Var1"), class = "data.frame")

ggplot(example, aes(x=Day, y=Growth, color = Animal, group = Location)) + geom_point() + geom_smooth(method = "lm", se = FALSE) + geom_smooth(method = "lm", aes(group = Animal)) + facet_grid(.~Location)

Clearly, there is a problem here with fitting an overall linear regression to each location. Because of the different starting points and lengths of each regression, different animals exert undue influence over the fit of the line (namely animals "c", "e", and "I"). For example, the actual slopes in X are -1, -1, 2. If animal "c" started at point 1 the slope would be zero.

When I asked this before I got quite a few answers telling me to fit a mixed effects model with Location and Individual as random effects, then to use the coefficients for Location within that model. But this does not escape the problem of spurious slopes no matter how you set up the model:

mem = lmer(Growth ~ (Day|Animal/Location), data = example)
sjp.lmer(mem, type = "pred", vars = c("Day", "Location"))
#ignore the fact that the model doesn't converge...

So...is it kosher to simply take an average of the slope coefficients for the linear regression of each individual, or is there a more rigorous way to get at this?

Answer 1

You can fit this in a single regression. The issue, to me, that it seems that you are dealing with is independence of observation. The example you gave is a really good example of the difference between within and between variance.

You can either do it as you see there, or you can toss everything within a regression model.

For example, each of those factors can be treated as binary variables (aka Boolean). Either the animal was in the shelter or it was not. So 1 if it was there, 0 if it was not. Doing this gives you an X by X dummy variable matrix. So for animal n, it will have a 1 in a single column and then its remaining row will be filled with zeros. This lets you put it into a regression. Then you include each of those variables (save one, which will be your baseline) within your regression model.

You regression ends up looking like this-

$growth = time + habitat_1 + habitat_2 + habitat_3 + error$.

I this case, $habitat_4$ is my baseline. The coefficients for $habitat_{1-3}$ will tell me the mean difference in growth between that habitat and $habitat_4$.