Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Discrepancy between plot_model output and estimate from lmer summary #424 #934

Open
tfreitas-kirk opened this issue Apr 16, 2024 · 2 comments

Comments

@tfreitas-kirk
Copy link

tfreitas-kirk commented Apr 16, 2024

Hello all.
I am running a mixed linear models, but I have been using the plot_model to generate the regression model graph
When I run linear mixed effect model (lmer function) the out put doesn't match with the intercept and slope for the 2 periods.Which function does plot_model use to generate the graph? Can you get the equation fro plot_model?

summary(lmer(log_Grazing_Percentage ~ Period + Count_RSSI_GE_Minus50 + Period:Count_RSSI_GE_Minus50 +
       (1 | Animal_id:Period),
     data = chap2))

REML criterion at convergence: -511.7

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-3.7017 -0.4717 -0.0746  0.3289  4.0732 

Random effects:
Groups           Name        Variance Std.Dev.
Animal_id:Period (Intercept) 0.005753 0.07585 
 Residual                     0.017927 0.13389 
Number of obs: 503, groups:  Animal_id:Period, 28

Fixed effects:
                               Estimate Std. Error         df t value Pr(>|t|)    
(Intercept)                    1.274435   0.058232  44.430851  21.885   <2e-16 ***
Period                        -0.012244   0.035985  40.761191  -0.340   0.7354    
Count_RSSI_GE_Minus50          0.006633   0.003261 496.218844   2.034   0.0425 *  
Period:Count_RSSI_GE_Minus50  -0.003454   0.001698 495.668407  -2.034   0.0425 *  
---
Signif. codes:  0***0.001**0.01*0.05.0.1 ‘ ’ 1

Correlation of Fixed Effects:
            (Intr) Period C_RSSI
Period      -0.951              
C_RSSI_GE_M -0.476  0.412       
P:C_RSSI_GE  0.489 -0.448 -0.987

My equations from my R output

Period 1 
1.274369 (intercept) + 0.006682 (slope)* (RSSI)
Period 2
1.274369 - 0.012289 = 1.26208 (intercept) + 0.006682 -0.003470 = 0.003212 (slope)* (RSSI)
1.26208 (intercept) +0.003212 (slope)* (RSSI) – To me my slope in period 2 should be negative 

This is the graph that I get from plot_model
Screenshot 2024-04-15 120151

Can someone please give me a hand with that?

@sj-perry
Copy link

I think the confusion here just stems from what the estimates in the regression table mean. This can get sort of tricky, which is why I always tell students to visualize model predictions with something like sjPlot ;)

This is what the four lines of the fixed effects are telling you:

(Intercept) - The predicted value of log_Grazing_Percentage when Period and Count_RSSI_GE_Minus50 are at zero. I'm assuming that Period is only really ever 1 or 2 in your data, so this value isn't easily interpretable.

Period - The difference in log_Grazing_Percentage for a change of 1 in Period when Count_RSSI_GE_Minus50 is at zero

Count_RSSI_GE_Minus50 - The estimated slope of this variable when Period is at zero. This isn't directly interpretable, because it is extrapolating beyond the observed values of Period to guess what the effect would be if the linear trends fit by the data continued.

Period:Count_RSSI_GE_Minus50 - Interactions can always be thought of in two reversible ways, so this is both how the effect of Period changes based on the values of Count_RSSI_GE_Minus50, and how the slope of the effect for Count_RSSI_GE_Minus50 changes based on Period.

So, to get the slope of Count_RSSI_GE_Minus50 when Period=1, you need to add multiply the interaction term by one (for Period = 1) and add it to the slope of your continuous variable estimated for when Period = 0.

# Slope of Count RSSI for Period = 1
0.006633 + (-0.003454*1)
[1] 0.003179

Then to get the slope when Period = 2, you multiply the interaction term by 2

# Slope of Count RSSI for Period = 2
> 0.006633 + (-0.003454*2)
[1] -0.000275

The following code tests some values to make predictions manually, and it looks like everything checks out:

Intercept <- 1.274435
Period <- -0.012244
Count_RSSI <- 0.006633
Interaction <- -0.003454

# Slope of Count RSSI for Period = 1
Count_RSSI + (Interaction*1)

# Slope of Count RSSI for Period = 2
Count_RSSI + (Interaction*2)

# Predicted log grazing when Period = 1 and RSSI = 0
Intercept + (Period*1) + (Count_RSSI*0) + (Interaction*0*1)

# Predicted log grazing when Period = 1 and RSSI = 80
Intercept + (Period*1) + (Count_RSSI*80) + (Interaction*80*1)

# Predicted log grazing when Period = 2 and RSSI = 0
Intercept + (Period*2) + (Count_RSSI*0) + (Interaction*0*2)

# Predicted log grazing when Period = 2 and RSSI = 80
Intercept + (Period*2) + (Count_RSSI*80) + (Interaction*80*2)

Also, I think you may want to fit Period as a random slope instead of the structure you currently have. Something like this may be more appropriate, although I would definitely check with someone in your field who has domain knowledge who is comfortable with using mixed models

lmer(log_Grazing_Percentage ~ Period + Count_RSSI_GE_Minus50 + Period:Count_RSSI_GE_Minus50 +
(1 + Period| Animal_id)

@strengejacke
Copy link
Owner

By default, predict() is used, however, you could also switch to a different package. The main points are:

  1. plot_model(type = "pred") plots adjusted predictions (sometimes "estimated marginal means") for the response value, i.e. it shows you the predicted value of your outcome, depending on meaningful (or representative) values of certain predictors (so called focal terms).

    The formula (by default) is just: y = b0 + b1*x1 + b2*x2 + ... bn*xn, where you vary the value of your focal terms, and hold the other values constant.

  2. Results may differ, depending on how you hold constant the non-focal terms. You can answer slightly different question based on this decision. I suggest reading this vignette for details.

Furthermore, if you're interested in predictions, I suggest using the ggeffects package directly. sjPlot::plot_model(type = "pred") is just a wrapper around that package, but less flexible.

The situation with mixed models is slightly more complicated, because you can condition on random effects or not. There's also a vignette here.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

3 participants