Intro

Miscellaneous functions are included to fill a need.

Examples

Find typical values

In many cases, predictions are made by holding the random effects to zero, which can be seen as the typical case. When the clusters carry substantive meaning, for example, school, country, or hospital, it might be of interest to note which are most ‘typical’, or at any other value relative to their peers. This will work for any model to which extract_random_effects applies.

library(lme4)

lmer_model <- lmer(Reaction ~ Days + (1 + Days | Subject), data = sleepstudy)

library(mixedup)

find_typical(lmer_model)  # closest RE to zero

[38;5;246m# A tibble: 2 × 7
[39m
  group_var effect    group  value    se lower_2.5 upper_97.5
  
[3m
[38;5;246m<chr>
[39m
[23m     
[3m
[38;5;246m<chr>
[39m
[23m     
[3m
[38;5;246m<fct>
[39m
[23m  
[3m
[38;5;246m<dbl>
[39m
[23m 
[3m
[38;5;246m<dbl>
[39m
[23m     
[3m
[38;5;246m<dbl>
[39m
[23m      
[3m
[38;5;246m<dbl>
[39m
[23m

[38;5;250m1
[39m Subject   Days      333   -
[31m0
[39m
[31m.
[39m
[31m224
[39m  2.31     -
[31m4
[39m
[31m.
[39m
[31m74
[39m       4.29

[38;5;250m2
[39m Subject   Intercept 335   -
[31m0
[39m
[31m.
[39m
[31m334
[39m 12.1     -
[31m24
[39m
[31m.
[39m
[31m0
[39m       23.3 


find_typical(lmer_model, probs = c(.1, .9))

[38;5;246m# A tibble: 4 × 8
[39m
  group_var effect    group  value    se lower_2.5 upper_97.5 probs
  
[3m
[38;5;246m<chr>
[39m
[23m     
[3m
[38;5;246m<chr>
[39m
[23m     
[3m
[38;5;246m<fct>
[39m
[23m  
[3m
[38;5;246m<dbl>
[39m
[23m 
[3m
[38;5;246m<dbl>
[39m
[23m     
[3m
[38;5;246m<dbl>
[39m
[23m      
[3m
[38;5;246m<dbl>
[39m
[23m 
[3m
[38;5;246m<chr>
[39m
[23m

[38;5;250m1
[39m Subject   Days      310    -
[31m5
[39m
[31m.
[39m
[31m45
[39m  2.31     -
[31m9
[39m
[31m.
[39m
[31m97
[39m     -
[31m0
[39m
[31m.
[39m
[31m931
[39m 10%  

[38;5;250m2
[39m Subject   Days      350     6.61  2.31      2.10     11.1   90%  

[38;5;250m3
[39m Subject   Intercept 370   -
[31m25
[39m
[31m.
[39m
[31m6
[39m  12.1     -
[31m49
[39m
[31m.
[39m
[31m3
[39m      -
[31m1
[39m
[31m.
[39m
[31m96
[39m  10%  

[38;5;250m4
[39m Subject   Intercept 331    22.3  12.1      -
[31m1
[39m
[31m.
[39m
[31m40
[39m     45.9   90%  

Extract heterogeneous variances

For nlme, only models with varIdent have been tested. For glmmTMB, this applies to models that use the diag function in the formula.

library(nlme)

lme_1 <- lme(
  distance ~ scale(age) + Sex,
  data = Orthodont,
  random = ~ 1 | Subject,
  weights = varIdent(form = ~ 1 | Sex)
)

extract_het_var(lme_1)

[38;5;246m# A tibble: 2 × 2
[39m
  group  variance
  
[3m
[38;5;246m<chr>
[39m
[23m     
[3m
[38;5;246m<dbl>
[39m
[23m

[38;5;250m1
[39m Male       3.10

[38;5;250m2
[39m Female     0.64

May include brms at some point, but for these, you can extract any parameters associated with sigma model for distributional models using the component argument of the other standard functions.

Extract correlation structure

With some packages one can estimate the residual correlation structure. Right now this works for nlme, glmmTMB, and brms objects.

base_model <-
  lme(Reaction ~ Days,
      random = ~ 1 + Days | Subject,
      data = sleepstudy)

lme_corSymm <-
  update(
    base_model,
    corr = corSymm(form = ~ 1 | Subject),
    data = dplyr::filter(sleepstudy, Days < 5)
  )

lme_corAR <-
  update(
    base_model,
    corr = corAR1(form = ~ Days),
    data = dplyr::filter(sleepstudy, Days < 5)
  )


extract_cor_structure(lme_corSymm)
       V1    V2     V3     V4     V5
V1  1.000 0.263 -0.460 -0.222 -0.092
V2  0.263 1.000  0.272  0.356  0.081
V3 -0.460 0.272  1.000  0.377 -0.417
V4 -0.222 0.356  0.377  1.000  0.572
V5 -0.092 0.081 -0.417  0.572  1.000


extract_cor_structure(lme_corAR)

[38;5;246m# A tibble: 1 × 1
[39m
    Phi
  
[3m
[38;5;246m<dbl>
[39m
[23m

[38;5;250m1
[39m 0.662

Extract model data

Just a wrapper for model.frame.

head(extract_model_data(lmer_model))

[38;5;246m# A tibble: 6 × 3
[39m
  Reaction  Days Subject
     
[3m
[38;5;246m<dbl>
[39m
[23m 
[3m
[38;5;246m<dbl>
[39m
[23m 
[3m
[38;5;246m<fct>
[39m
[23m  

[38;5;250m1
[39m     250.     0 308    

[38;5;250m2
[39m     259.     1 308    

[38;5;250m3
[39m     251.     2 308    

[38;5;250m4
[39m     321.     3 308    

[38;5;250m5
[39m     357.     4 308    

[38;5;250m6
[39m     415.     5 308