¿Cómo extraigo los efectos fijos de lmer por observación?

Question

I tienen un objeto LME, construido a partir de algunas medidas repetidas de datos de ingesta de nutrientes (dos períodos de admisión 24-horas por IdPersona):

Male.lme2 <- lmer(BoxCoxXY ~ -1 + AgeFactor + IntakeDay + (1|RespondentID), 
    data = Male.Data, 
    weights = SampleWeight) 

y puedo recuperar correctamente los efectos aleatorios por RespondentID utilizando ranef(Male.lme1). También me gustaría recopilar el resultado de los efectos fijos en RespondentID. coef(Male.lme1) no proporciona exactamente lo que necesito, como muestro a continuación.

> summary(Male.lme1) 
Linear mixed model fit by REML 
Formula: BoxCoxXY ~ AgeFactor + IntakeDay + (1 | RespondentID) 
    Data: Male.Data 
    AIC BIC logLik deviance REMLdev 
    9994 10039 -4990  9952 9980 
Random effects: 
Groups  Name  Variance Std.Dev. 
RespondentID (Intercept) 0.19408 0.44055 
Residual     0.37491 0.61230 
Number of obs: 4498, groups: RespondentID, 2249 

Fixed effects: 
        Estimate Std. Error t value 
(Intercept)   13.98016 0.03405 410.6 
AgeFactor4to8  0.50572 0.04084 12.4 
AgeFactor9to13  0.94329 0.04159 22.7 
AgeFactor14to18  1.30654 0.04312 30.3 
IntakeDayDay2Intake -0.13871 0.01809 -7.7 

Correlation of Fixed Effects: 
      (Intr) AgFc48 AgF913 AF1418 
AgeFactr4t8 -0.775      
AgeFctr9t13 -0.761 0.634    
AgFctr14t18 -0.734 0.612 0.601  
IntkDyDy2In -0.266 0.000 0.000 0.000 

He anexado los resultados ajustados a mis datos, head(Male.Data) muestra

NutrientID RespondentID Gender Age SampleWeight IntakeDay IntakeAmt AgeFactor BoxCoxXY lmefits 
2   267  100020  1 12 0.4952835 Day1Intake 12145.852  9to13 15.61196 15.22633 
7   267  100419  1 14 0.3632839 Day1Intake 9591.953 14to18 15.01444 15.31373 
8   267  100459  1 11 0.4952835 Day1Intake 7838.713  9to13 14.51458 15.00062 
12  267  101138  1 15 1.3258785 Day1Intake 11113.266 14to18 15.38541 15.75337 
14  267  101214  1 6 2.1198688 Day1Intake 7150.133  4to8 14.29022 14.32658 
18  267  101389  1 5 2.1198688 Day1Intake 5091.528  4to8 13.47928 14.58117 

El primer par de líneas de coef(Male.lme1) son:

$RespondentID 
     (Intercept) AgeFactor4to8 AgeFactor9to13 AgeFactor14to18 IntakeDayDay2Intake 
100020 14.28304  0.5057221  0.9432941  1.306542   -0.1387098 
100419 14.00719  0.5057221  0.9432941  1.306542   -0.1387098 
100459 14.05732  0.5057221  0.9432941  1.306542   -0.1387098 
101138 14.44682  0.5057221  0.9432941  1.306542   -0.1387098 
101214 13.82086  0.5057221  0.9432941  1.306542   -0.1387098 
101389 14.07545  0.5057221  0.9432941  1.306542   -0.1387098 

Para demostrar si los coef resultados se refieren a las estimaciones ajustadas en Male.Data (que fueron tomadas usando Male.Data$lmefits <- fitted(Male.lme1), para el primer RespondentID, que tiene la edad) Nivel de factor 9-13: - el valor ajustado es 15.22633, que es igual - a partir de los coeffs - (Intercept) + (AgeFactor9-13) = 14.28304 + 0.9432941

¿Hay un comando inteligente para mí utilizar que hará quiero Quiero automáticamente, lo que es extraer el efecto fijo estimar para cada sujeto, o ¿me enfrento a una serie de declaraciones if tratando de aplicar el nivel correcto de AgeFactor a cada sujeto para obtener el estimado de efecto fijo correcto, después de deducir la contribución del efecto aleatorio de Interceptar?

Actualización, disculpas, estaba tratando de reducir la salida que estaba proporcionando y se olvidó de str(). La salida es: no se están utilizando

>str(Male.Data) 
'data.frame': 4498 obs. of 11 variables: 
$ NutrientID : int 267 267 267 267 267 267 267 267 267 267 ... 
$ RespondentID: Factor w/ 2249 levels "100020","100419",..: 1 2 3 4 5 6 7 8 9 10 ... 
$ Gender  : int 1 1 1 1 1 1 1 1 1 1 ... 
$ Age   : int 12 14 11 15 6 5 10 2 2 9 ... 
$ BodyWeight : num 51.6 46.3 46.1 63.2 28.4 18 38.2 14.4 14.6 32.1 ... 
$ SampleWeight: num 0.495 0.363 0.495 1.326 2.12 ... 
$ IntakeDay : Factor w/ 2 levels "Day1Intake","Day2Intake": 1 1 1 1 1 1 1 1 1 1 ... 
$ IntakeAmt : num 12146 9592 7839 11113 7150 ... 
$ AgeFactor : Factor w/ 4 levels "1to3","4to8",..: 3 4 3 4 2 2 3 1 1 3 ... 
$ BoxCoxXY : num 15.6 15 14.5 15.4 14.3 ... 
$ lmefits  : num 15.2 15.3 15 15.8 14.3 ... 

El peso corporal y sexo (estos son los datos machos, por lo que todos los valores de género son los mismos) y el NutrientID se fija de manera similar para los datos.

He estado haciendo horribles declaraciones de ifelse que publiqué, por lo que intentaré su sugerencia inmediatamente. :)

Update2: esto funciona perfectamente con mis datos actuales y debe ser a prueba de futuro para nuevos datos, gracias a DWin por la ayuda adicional en el comentario para esto. :)

AgeLevels <- length(unique(Male.Data$AgeFactor)) 
Temp <- as.data.frame(fixef(Male.lme1)['(Intercept)'] + 
c(0,fixef(Male.lme1)[2:AgeLevels])[ 
     match(Male.Data$AgeFactor, c("1to3", "4to8", "9to13","14to18", "19to30","31to50","51to70","71Plus"))] + 
c(0,fixef(Male.lme1)[(AgeLevels+1)])[ 
     match(Male.Data$IntakeDay, c("Day1Intake","Day2Intake"))]) 
names(Temp) <- c("FxdEffct") 

Answer 1

Va a ser algo como esto (aunque realmente debe haber dado los resultados de str (Male.Data) debido a los resultados del modelo no nos dice los niveles de los factores de los valores basales :)

#First look at the coefficients 
fixef(Male.lme2) 

#Then do the calculations 
fixef(Male.lme2)[`(Intercept)`] + 
c(0,fixef(Male.lme2)[2:4])[ 
      match(Male.Data$AgeFactor, c("1to3", "4to8", "9to13","14to18"))] + 
c(0,fixef(Male.lme2)[5])[ 
      match(Male.Data$IntakeDay, c("Day1Intake","Day2Intake"))] 

básicamente estás ejecutando los datos originales a través de una función match para recoger el coeficiente (s) correcta para agregar a la intersección ... que será 0 si los datos son la base del nivel del factor (cuyo deletreo estoy adivinando.)

EDITAR: Me acabo de dar cuenta de que pone un "-1" en la fórmula así que quizás todos sus términos AgeFactor se enumeran en la salida y puede descubrir el 0 en el vector de coeficientes y el nivel inventado AgeFactor en la tabla de coincidencias vector.

Answer 2

A continuación es cómo siempre me ha parecido más fácil extraer los efectos fijos de los individuos y los componentes de efectos aleatorios en el lme4 -package. En realidad, extrae el ajuste correspondiente a cada observación. Suponiendo que tenemos un modelo de efectos mixtos de la forma:

y = Xb + Zu + e 

donde Xb son los efectos fijos y Zu son los efectos aleatorios, podemos extraer los componentes (usando lme4 's sleepstudy como ejemplo):

library(lme4) 
fm1 <- lmer(Reaction ~ Days + (Days|Subject), sleepstudy) 

# Xb 
fix <- getME(fm1,'X') %*% fixef(fm1) 
# Zu 
ran <- t(as.matrix(getME(fm1,'Zt'))) %*% unlist(ranef(fm1)) 
# Xb + Zu 
fixran <- fix + ran 

Sé que esto funciona como un enfoque general para extraer componentes de modelos lineales de efectos mixtos. Para modelos no lineales, la matriz modelo X contiene repeticiones y es posible que tenga que adaptar un poco el código anterior. Aquí hay alguna salida validación, así como el uso de la visualización de celosía:

> head(cbind(fix, ran, fixran, fitted(fm1))) 
     [,1]  [,2]  [,3]  [,4] 
[1,] 251.4051 2.257187 253.6623 253.6623 
[2,] 261.8724 11.456439 273.3288 273.3288 
[3,] 272.3397 20.655691 292.9954 292.9954 
[4,] 282.8070 29.854944 312.6619 312.6619 
[5,] 293.2742 39.054196 332.3284 332.3284 
[6,] 303.7415 48.253449 351.9950 351.9950 

# Xb + Zu 
> all(round((fixran),6) == round(fitted(fm1),6)) 
[1] TRUE 

# e = y - (Xb + Zu) 
> all(round(resid(fm1),6) == round(sleepstudy[,"Reaction"]-(fixran),6)) 
[1] TRUE 

nobs <- 10 # 10 observations per subject 
legend = list(text=list(c("y", "Xb + Zu", "Xb")), lines = list(col=c("blue", "red", "black"), pch=c(1,1,1), lwd=c(1,1,1), type=c("b","b","b"))) 
require(lattice) 
xyplot(
    Reaction ~ Days | Subject, data = sleepstudy, 
    panel = function(x, y, ...){ 
     panel.points(x, y, type='b', col='blue') 
     panel.points(x, fix[(1+nobs*(panel.number()-1)):(nobs*(panel.number()))], type='b', col='black') 
     panel.points(x, fixran[(1+nobs*(panel.number()-1)):(nobs*(panel.number()))], type='b', col='red') 
    }, 
    key = legend 
)