Skip to contents
library(modsem)
#> This is modsem (1.0.22). Please report any bugs!

Here we show an example of a random intercept cross-lagged panel model with a within ×\times within interaction effect. The example below is based on Testing for Within ×\times Within and Between ×\times Within Moderation Using Random Intercept Cross-Lagged Panel Models. In particular, it is based on Section 8.4. Within ×\times Within Interaction Using the Within-Person Component of a Time-Varying Moderator, and the corresponding Mplus output files on OSF. In the original article, the Bayes estimator was used, but here we use LMS. The LMS approach in Mplus integrates along four dimensions, making estimation quite impractical for this model. The LMS estimator in modsem, however, only has to integrate along two dimensions, making estimation viable.

Simulate Data

library(mvtnorm)
set.seed(235910)

N <- 10000

# Simulate between-person random intercepts
phi <- matrix(c(
  0.806, 0.481, 0.438,
  0.481, 0.758, 0.469,
  0.438, 0.469, 1.290
), nrow = 3, byrow = TRUE)

Xi <- rmvnorm(
  N,
  mean = c(0, 0, 0),
  sigma = phi
)

RIemo <- Xi[, 1]
RIper <- Xi[, 2]
RIcon <- Xi[, 3]


# Time 3
psi3 <- matrix(c(
  1.234, 0.328, 0.478,
  0.328, 1.614, 0.427,
  0.478, 0.427, 2.743
), nrow = 3, byrow = TRUE)

zeta3 <- rmvnorm(
  N,
  mean = c(0, 0, 0),
  sigma = psi3
)

wemo_3 <- zeta3[, 1]
wper_3 <- zeta3[, 2]
wcon_3 <- zeta3[, 3]


# Time 5
psi5 <- matrix(c(
  1.489, 0.397, 0.301,
  0.397, 1.143, 0.202,
  0.301, 0.202, 0.764
), nrow = 3, byrow = TRUE)

zeta5 <- rmvnorm(
  N,
  mean = c(0, 0, 0),
  sigma = psi5
)

wemo_5 <- 0.142 * wemo_3 + 0.071 * wper_3 + 0.086 * wcon_3 - 0.010 * (wcon_3 * wper_3) + zeta5[,1]
wper_5 <- 0.018 * wemo_3 + 0.100 * wper_3 + 0.050 * wcon_3 + zeta5[,2]
wcon_5 <- -0.008 * wemo_3 + 0.006 * wper_3 + 0.105 * wcon_3 + zeta5[,3]

# Time 7
psi7 <- matrix(c(
  1.808, 0.521, 0.475,
  0.521, 1.287, 0.311,
  0.475, 0.311, 0.881
), nrow = 3, byrow = TRUE)

zeta7 <- rmvnorm(
  N,
  mean = c(0, 0, 0),
  sigma = psi7
)

wemo_7 <- 0.361 * wemo_5 + 0.097 * wper_5 + 0.147 * wcon_5 + 0.105 * (wcon_5 * wper_5) + zeta7[,1]
wper_7 <- 0.030 * wemo_5 + 0.348 * wper_5 + 0.123 * wcon_5 + zeta7[,2]
wcon_7 <- 0.074 * wemo_5 + 0.056 * wper_5 + 0.086 * wcon_5 + zeta7[,3]

# Indicators/Observed Variables
data <- data.frame(
  emo_3 = 1.385 + RIemo + wemo_3 + rnorm(N, 0, sqrt(.200)),
  emo_5 = 1.407 + RIemo + wemo_5 + rnorm(N, 0, sqrt(.200)),
  emo_7 = 1.528 + RIemo + wemo_7 + rnorm(N, 0, sqrt(.200)),

  per_3 = 1.563 + RIper + wper_3 + rnorm(N, 0, sqrt(.200)),
  per_5 = 1.176 + RIper + wper_5 + rnorm(N, 0, sqrt(.200)),
  per_7 = 1.242 + RIper + wper_7 + rnorm(N, 0, sqrt(.200)),

  con_3 = 2.827 + RIcon + wcon_3 + rnorm(N, 0, sqrt(.200)),
  con_5 = 1.517 + RIcon + wcon_5 + rnorm(N, 0, sqrt(.200)),
  con_7 = 1.402 + RIcon + wcon_7 + rnorm(N, 0, sqrt(.200))
)

Fit The Model

model.inp <- '
  # Create between components (random intercepts)
  RIemo =~ 1 * emo_3 + 1 * emo_5 + 1 * emo_7;
  RIcon =~ 1 * con_3 + 1 * con_5 + 1 * con_7;

  # Create within-person centered variables
  wemo_3 =~ 1 * emo_3;
  wemo_5 =~ 1 * emo_5;
  wemo_7 =~ 1 * emo_7;

  wcon_3 =~ 1 * con_3;
  wcon_5 =~ 1 * con_5;
  wcon_7 =~ 1 * con_7;

  # Moderator also needs to be decomposed into within and between parts
  RIper =~ 1 * per_3 + 1 * per_5 + 1 * per_7;

  wper_3 =~ 1 * per_3;
  wper_5 =~ 1 * per_5;
  wper_7 =~ 1 * per_7;

  # Constrain the measurement error variances close to zero
  # to allow for reasonable imputation times
  con_3 ~~ 0.2 * con_3
  con_5 ~~ 0.2 * con_5
  con_7 ~~ 0.2 * con_7
  emo_3 ~~ 0.2 * emo_3
  emo_5 ~~ 0.2 * emo_5
  emo_7 ~~ 0.2 * emo_7
  per_3 ~~ 0.2 * per_3
  per_5 ~~ 0.2 * per_5
  per_7 ~~ 0.2 * per_7

  # Estimate the covariance between the random intercepts
  RIemo ~~ RIper
  RIemo ~~ RIcon
  RIper ~~ RIcon

  # Estimate the lagged effects between
  # the within-person centered variables
  wemo_7 ~ wemo_5 + wper_5 + wcon_5
  wper_7 ~ wemo_5 + wper_5 + wcon_5
  wcon_7 ~ wemo_5 + wper_5 + wcon_5

  wemo_5 ~ wemo_3 + wper_3 + wcon_3
  wper_5 ~ wemo_3 + wper_3 + wcon_3
  wcon_5 ~ wemo_3 + wper_3 + wcon_3

  # Specify interaction terms between within-person centred
  # conduct problems and the within-person centered peer problems
  # predict within-person centred emotional problems with interaction
  wemo_7 ~ wcon_5:wper_5
  wemo_5 ~ wcon_3:wper_3

  # Estimate the covariance between the within-person
  # components at the first wave
  wemo_3 ~~ wper_3
  wemo_3 ~~ wcon_3
  wper_3 ~~ wcon_3

  # Estimate the covariances between the residuals of
  # the within-person components (the innovations)
  wemo_5 ~~ wper_5
  wemo_5 ~~ wcon_5
  wper_5 ~~ wcon_5

  wemo_7 ~~ wper_7
  wemo_7 ~~ wcon_7
  wper_7 ~~ wcon_7
'

fit.lms <- modsem(
  model.syntax = model.inp,
  data              = data,
  method            = "lms",
  nodes             = 32,
  optimize          = FALSE, # we're currently unable to optimize starting parameters here
  orthogonal.x      = TRUE,  # make sure the model is identifiable
  orthogonal.y      = TRUE,  # not strictly necessary for this model in particular
  auto.split.syntax = TRUE   # allow eta x eta interactions
)

summary(fit.lms)
#> Warning: modsem->fit_modsem_da():  
#>    Unable to retrieve fit measures for lavaan model! Message: no applicable 
#>    method for 'fitMeasures' applied to an object of class "NULL"
#> 
#> modsem (1.0.22) ended normally after 72 iterations
#> 
#>   Estimator                                        LMS
#>   Optimization method                       EMA-NLMINB
#>   Number of model parameters                        53
#> 
#>   Number of observations                         10000
#> 
#> Loglikelihood and Information Criteria:
#>   Loglikelihood                             -157580.72
#>   Akaike (AIC)                               315267.44
#>   Bayesian (BIC)                             315649.59
#>  
#> Numerical Integration:
#>   Points of integration (per dim)                   32
#>   Dimensions                                         2
#>   Total points of integration                     1024
#> 
#> Fit Measures for Baseline Model (H0):
#>                                               Standard
#>   Chi-square                                      2.36
#>   Degrees of Freedom (Chi-square)                    3
#>   P-value (Chi-square)                           0.501
#>   RMSEA                                          0.000
#>                                                       
#>   Loglikelihood                             -157592.49
#>   Akaike (AIC)                               315286.98
#>   Bayesian (BIC)                             315654.71
#>  
#> Comparative Fit to H0 (LRT test):
#>   Loglikelihood change                           11.77
#>   Difference test (D)                            23.53
#>   Degrees of freedom (D)                             2
#>   P-value (D)                                    0.000
#>  
#> R-Squared Interaction Model (H1):
#>   wper_5                                         0.029
#>   wcon_5                                         0.044
#>   wemo_5                                         0.056
#>   wcon_7                                         0.040
#>   wper_7                                         0.136
#>   wemo_7                                         0.155
#> R-Squared Baseline Model (H0):
#>   wper_5                                         0.029
#>   wcon_5                                         0.046
#>   wemo_5                                         0.058
#>   wcon_7                                         0.043
#>   wper_7                                         0.136
#>   wemo_7                                         0.149
#> R-Squared Change (H1 - H0):
#>   wper_5                                         0.000
#>   wcon_5                                        -0.002
#>   wemo_5                                        -0.002
#>   wcon_7                                        -0.004
#>   wper_7                                         0.000
#>   wemo_7                                         0.006
#> 
#> Parameter Estimates:
#>   Coefficients                          unstandardized
#>   Information                                 observed
#>   Standard errors                             standard
#>  
#> Latent Variables:
#>                  Estimate  Std.Error  z.value  P(>|z|)
#>   RIemo =~      
#>     emo_3           1.000                             
#>     emo_5           1.000                             
#>     emo_7           1.000                             
#>   RIcon =~      
#>     con_5           1.000                             
#>     con_3           1.000                             
#>     con_7           1.000                             
#>   wemo_3 =~     
#>     emo_3           1.000                             
#>   wcon_3 =~     
#>     con_3           1.000                             
#>   RIper =~      
#>     per_3           1.000                             
#>     per_5           1.000                             
#>     per_7           1.000                             
#>   wper_3 =~     
#>     per_3           1.000                             
#>   wcon_5 =~     
#>     con_5           1.000                             
#>   wper_5 =~     
#>     per_5           1.000                             
#>   wemo_7 =~     
#>     emo_7           1.000                             
#>   wper_7 =~     
#>     per_7           1.000                             
#>   wcon_7 =~     
#>     con_7           1.000                             
#>   wemo_5 =~     
#>     emo_5           1.000                             
#> 
#> Regressions:
#>                  Estimate  Std.Error  z.value  P(>|z|)
#>   wcon_5 ~      
#>     wemo_3         -0.009      0.020   -0.439    0.661
#>     wcon_3          0.120      0.013    8.962    0.000
#>     wper_3         -0.018      0.015   -1.169    0.242
#>   wper_5 ~      
#>     wemo_3          0.007      0.025    0.281    0.778
#>     wcon_3          0.050      0.013    4.017    0.000
#>     wper_3          0.115      0.023    5.032    0.000
#>   wemo_7 ~      
#>     wcon_5          0.173      0.031    5.649    0.000
#>     wper_5          0.060      0.022    2.723    0.006
#>     wemo_5          0.370      0.021   17.840    0.000
#>     wcn_5:wp_5      0.112      0.024    4.696    0.000
#>   wper_7 ~      
#>     wcon_5          0.146      0.029    5.087    0.000
#>     wper_5          0.345      0.023   14.705    0.000
#>     wemo_5          0.025      0.017    1.463    0.143
#>   wcon_7 ~      
#>     wcon_5          0.087      0.047    1.862    0.063
#>     wper_5          0.076      0.024    3.144    0.002
#>     wemo_5          0.076      0.020    3.834    0.000
#>   wemo_5 ~      
#>     wemo_3          0.158      0.033    4.826    0.000
#>     wcon_3          0.097      0.013    7.481    0.000
#>     wper_3          0.058      0.019    2.960    0.003
#>     wcn_3:wp_3     -0.003      0.008   -0.327    0.744
#> 
#> Intercepts:
#>                  Estimate  Std.Error  z.value  P(>|z|)
#>    .con_5           1.525      0.015  100.223    0.000
#>    .con_3           2.815      0.021  136.997    0.000
#>    .emo_3           1.371      0.015   91.978    0.000
#>    .emo_5           1.393      0.016   84.913    0.000
#>    .emo_7           1.518      0.019   80.404    0.000
#>    .con_7           1.383      0.015   89.842    0.000
#>    .per_3           1.562      0.016   98.151    0.000
#>    .per_5           1.176      0.015   80.542    0.000
#>    .per_7           1.226      0.015   79.216    0.000
#> 
#> Covariances:
#>                  Estimate  Std.Error  z.value  P(>|z|)
#>   RIemo ~~      
#>     RIcon           0.414      0.028   14.720    0.000
#>     RIper           0.499      0.032   15.769    0.000
#>   RIcon ~~      
#>     RIper           0.434      0.028   15.669    0.000
#>   wemo_3 ~~     
#>     wcon_3          0.499      0.033   14.961    0.000
#>     wper_3          0.316      0.034    9.273    0.000
#>   wper_3 ~~     
#>     wcon_3          0.438      0.034   12.803    0.000
#>  .wper_5 ~~     
#>    .wcon_5          0.240      0.028    8.616    0.000
#>  .wper_7 ~~     
#>    .wemo_7          0.483      0.021   23.308    0.000
#>  .wcon_7 ~~     
#>    .wemo_7          0.488      0.022   22.206    0.000
#>    .wper_7          0.318      0.020   15.851    0.000
#>  .wemo_5 ~~     
#>    .wcon_5          0.334      0.028   12.060    0.000
#>    .wper_5          0.419      0.029   14.284    0.000
#> 
#> Variances:
#>                  Estimate  Std.Error  z.value  P(>|z|)
#>    .con_5           0.200                             
#>    .con_3           0.200                             
#>    .emo_3           0.200                             
#>    .emo_5           0.200                             
#>    .emo_7           0.200                             
#>    .con_7           0.200                             
#>    .per_3           0.200                             
#>    .per_5           0.200                             
#>    .per_7           0.200                             
#>     RIemo           0.806      0.045   17.980    0.000
#>     RIcon           1.267      0.043   29.176    0.000
#>     wemo_3          1.217      0.047   25.795    0.000
#>     wcon_3          2.756      0.054   50.871    0.000
#>     RIper           0.746      0.039   19.182    0.000
#>     wper_3          1.587      0.045   35.270    0.000
#>    .wcon_5          0.810      0.042   19.402    0.000
#>    .wper_5          1.152      0.039   29.686    0.000
#>    .wemo_7          1.797      0.032   55.711    0.000
#>    .wper_7          1.252      0.025   49.724    0.000
#>    .wcon_7          0.864      0.035   24.755    0.000
#>    .wemo_5          1.466      0.041   35.869    0.000