Is there a difference whether to fit single subjects and use the mean out of them or to fit all together directly?

Question

I have 10 devices which I fit with a linear model in the form of

log(Y) ~ log(X)

then I get a slope for each device of which I use the mean to provide their "general" behaviour.

Could I, instead, fit just all devices all together by having only a single slope at the end so I don't have to apply a mean? Will there be a difference?

Answer 1

You are looking for a mixed model (i.e. model which will allow you to use random effects along with the fixed effects you already use). You want to include device-specific random effect for the slope of log(X). Say you have a data.frame with columns logY, logX, deviceNumber. Than you get your desired result by:

require(lme4)
m <- lmer(logY ~ logX + (1 + logX|deviceNumber))

This will give you a general mean slope (over all devices), which will be something like an "average over device specific slopes". Along with proper standard errors and p-values.

It is better to fit it like this in a single model, than fitting in separate models and then use the mean slope. It is because the mean slopes are uncertain values with given standard error, and these standard errors wouldn't be involved in estimation of the standard error of the mean slope.

It would also be a mistake to fit them in a single model without random effects, like lm(logY ~ logX). Because the linear model assumes independence of the data points, but the data points of the same device are likely to be more similar or have more similar slope, so the assumption of the independence would be violated - which would result in incorrect standard errors.

Answer 2

I would suggest trying some market mix models which does randomization of subjects and provide estimates at subject level.

Answer 3

I will give you an example in R, using two contrast type, the default one and one that could be useful to you. The default contrast in R is to compare each level of a factor to the "reference", i.e. a baseline, the effect is estimated in the intercept term. Another set of contrasts is called "contr.sum", the intercept term now represents the mean intercept of all factor levels.

An example with mtcars toy dataset, predicting mpg from gear, which is a factor with 3 levels.

mtcars$gear=factor(mtcars$gear)
summary(lm(mpg~gear,data=mtcars))

which results in

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)   16.107      1.216  13.250 7.87e-14 ***
gear4          8.427      1.823   4.621 7.26e-05 ***
gear5          5.273      2.431   2.169   0.0384 *  

the intercept term has the mean score of the first level of gear (gear==3), we can check this by manually calculating the mean of mpg for gear3

mean(mtcars$mpg[mtcars$gear==3])
[1] 16.10667

Now we use the "contr.sum" contrast

summary(lm(mpg~gear,data=mtcars,contrasts=list(gear="contr.sum")))

which results in

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  20.6733     0.9284  22.267  < 2e-16 ***
gear1        -4.5667     1.1639  -3.924 0.000492 ***
gear2         3.8600     1.2156   3.175 0.003534 ** 

now the intercept term has the overall mean of each factor, we can check manually

tapply(mtcars$mpg,mtcars$gear,mean)
       3        4        5 
16.10667 24.53333 21.38000

and the overall mean

mean(tapply(mtcars$mpg,mtcars$gear,mean))
[1] 20.67333

The other terms are interpreted as deviations from the overall mean.