3
$\begingroup$

I'm studying the difference of feed intake between more than 150 horses. From each horse we have their feed intake at different week points.

My data is not normal-distributed. Therefore, in order to construct a linear mixed model (fixed effect: weeks and random effect: horse) I need to go through glmer() instead of lmer() in package lme4().

In glmer, the lower AIC obtained is through Gamma distribution.

Nevertheless I cannot understand that in glmer()model the variance of my random effect (horse) is only 6%(0.0646), in contrast if we considering normality and performing lmer(), then this 6% increase up to 33% (0.33159).

I'm more agree with the variance of my random effect stated by lmer() model instead of glmer() because a Horse will be an important variable to explain the model.

But, how I should interpret my results in order to discuss it properly?


# ASSUMING NORMALITY OF MY RESPONSE VARIABLE (Feed_kg_DM_day)

m_avg=lmer (Feed_kg_DM_day ~ factor(Week) + (1|Horse), data=dietdef)
summary(m_avg)

# Linear mixed model fit by REML ['lmerMod']
# Formula: Feed_kg_DM_day ~ factor(Week) + (1 | Horse)
# Data: dietdef
# 
# REML criterion at convergence: 551.5
# 
# Scaled residuals: 
#   Min      1Q  Median      3Q     Max 
# -8.6200 -0.1579 -0.0242  0.1094  7.6845 
# 
# Random effects:
#   Groups   Name        Variance Std.Dev.
# Horse    (Intercept) 0.33159  0.5758  
# Residual             0.05074  0.2253  
# Number of obs: 2909, groups:  Horse, 194
# 
# Fixed effects:
#   Estimate Std. Error t value
# (Intercept)     3.162e+00  4.439e-02  71.222
# factor(Week)13 -4.479e-13  2.287e-02   0.000
# factor(Week)15 -6.610e-03  2.290e-02  -0.289
# factor(Week)17  2.000e-03  2.287e-02   0.087
# factor(Week)19  2.000e-03  2.287e-02   0.087
# factor(Week)2  -1.233e-02  2.287e-02  -0.539
# factor(Week)22  2.000e-03  2.287e-02   0.087
# factor(Week)24  2.000e-03  2.287e-02   0.087
# factor(Week)4  -1.233e-02  2.287e-02  -0.539
# factor(Week)43 -8.332e-02  2.287e-02  -3.643
# factor(Week)45  3.582e-02  2.287e-02   1.566
# factor(Week)47 -1.066e-03  2.287e-02  -0.047
# factor(Week)48 -1.852e-03  2.287e-02  -0.081
# factor(Week)7  -8.765e-03  2.287e-02  -0.383
# factor(Week)9  -5.115e-03  2.287e-02  -0.224


# ASSUMING NOT NORMALITY

modelg2<-glmer(Feed_kg_DM_day~ factor(Week)+(1|Horse), data=dietdef, family=Gamma(link=identity))
summary(modelg2)

# Generalized linear mixed model fit by maximum likelihood (Laplace Approximation) ['glmerMod']
# Family: Gamma  ( identity )
# Formula: Feed_kg_DM_day ~ factor(Week) + (1 | Horse)
# Data: dietdef
# 
# AIC      BIC   logLik deviance df.resid 
# 92.9    194.5    -29.4     58.9     2892 
# 
# Scaled residuals: 
#   Min      1Q  Median      3Q     Max 
# -6.7920 -0.1935 -0.0464  0.1683  6.1595 
# 
# Random effects:
#   Groups   Name        Variance Std.Dev.
# Horse    (Intercept) 0.064600 0.2542  
# Residual             0.006972 0.0835  
# Number of obs: 2909, groups:  Horse, 194
# 
# Fixed effects:
#   Estimate Std. Error t value Pr(>|z|)    
# (Intercept)     3.304e+00  6.262e-02  52.753  < 2e-16 ***
#   factor(Week)13 -6.481e-06  2.196e-02   0.000  0.99976    
# factor(Week)15 -8.264e-03  2.197e-02  -0.376  0.70678    
# factor(Week)17 -5.430e-03  2.195e-02  -0.247  0.80466    
# factor(Week)19 -5.430e-03  2.195e-02  -0.247  0.80465    
# factor(Week)2  -1.186e-02  2.191e-02  -0.541  0.58832    
# factor(Week)22 -5.429e-03  2.195e-02  -0.247  0.80469    
# factor(Week)24 -5.434e-03  2.195e-02  -0.248  0.80451    
# factor(Week)4  -1.186e-02  2.191e-02  -0.541  0.58846    
# factor(Week)43 -3.059e-02  2.181e-02  -1.403  0.16067    
# factor(Week)45  7.184e-02  2.222e-02   3.233  0.00123 ** 
#   factor(Week)47  3.180e-02  2.206e-02   1.442  0.14942    
# factor(Week)48  4.943e-02  2.206e-02   2.241  0.02505 *  
#   factor(Week)7  -9.628e-03  2.192e-02  -0.439  0.66052    
# factor(Week)9  -7.088e-03  2.193e-02  -0.323  0.74657    
# ---
#   Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
``` 
$\endgroup$
  • 1
    $\begingroup$ Not sure why people are voting to close this since there is a statistical issue here (I think, or at least an issue with the OP understanding of the random effects) which I have tried to answer. $\endgroup$ – mdewey Apr 26 '19 at 11:47
3
$\begingroup$

I do not think you have a problem here as the scales are different. Notice in particular the ratio between the random effects variances for horse and the residual. In both cases the majority of the variance is between horses which is what you say you expected on theoretical grounds rather than within horses. For the gamma model the ratio is about 9 for the normal model about 6.

Note that the (estimated) variance for horse is the variance of the distribution of the random intercepts for horse. The residual variance is the remaining variance. Each is expressed on the scale of the model so in one case on the log scale on the other the linear scale. You cannot convert from one to the other and I do not think your attempt to convert to a percentage is correct and even if it were I do not think it would be helpful

$\endgroup$
  • $\begingroup$ Thank you mdewey for your response and solving my doubt! Nevertheless, it exists a procedure to rescale my glmer() variance results in order to obtain the lmer() variance values? Therefore, try to use the same scale as lmer() model to express variance of random effect $\endgroup$ – ribelles Apr 26 '19 at 14:12
  • $\begingroup$ As far as I know there is no way of converting variances on the log scale to the linear scale. $\endgroup$ – mdewey Apr 26 '19 at 14:33
  • $\begingroup$ I'm newbie in StackExchange and the previous comment it was not my definitive doubt. Sorry for my disorder. By the lmer model I can interpret my variance as percentage(33% of variance due to horse as random effect). But in glmer () model because it uses another scale, then I understand that I cannot transform this 0.06 (variance of random effect of horse) to percentage. But then, how I can discuss the results of the model used. In which units is expressed the variance of random effects in glmer() models? Thanks another time. $\endgroup$ – ribelles Apr 26 '19 at 15:04

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.