RICLPM_isolation_adhd.Rmd

---
title: "Separate mother and teacher RI-CLPM for social isolation and ADHD symptoms in childhood"
output:  
  html_document:
    df_print: paged
    toc: yes
    toc_depth: 2
    toc_float:
      collapsed: false
    number_sections: false
    highlight: monochrome
    theme: flatly
    code_folding: hide
    includes:
      after_body: footer.html
---

```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = TRUE,
                      comment = NA,
                      prompt = FALSE,
                      cache = FALSE,
                      message = FALSE,
                      warning = FALSE,
                      results = 'asis')

options(bitmapType = 'quartz') # to render fonts better
```

```{r Clear global environment, include=FALSE}
remove(list = ls())
```

```{r Load packages, include=FALSE}
library(knitr)
library(haven)
library(lavaan)
library(ChiBarSq.DiffTest)
library(tidyverse)
library(dplyr) #conflicts with tidyverse for e.g. rename and row_number
```

```{r source the data file path, include=FALSE}
#source raw data directory: data_raw and data included
source("../isolation_ADHD_data_path.R")
```

```{r read in dta data file}
dat.raw <- read_dta(paste0(data.raw_path, "Katie_19Jan22.dta"))

dat <- dat.raw %>%
  dplyr::select(id = atwinid,
         sampsex,
         seswq35,
         sisoem5,  # social isolation mother report
         sisoem7,
         sisoem10,
         sisoem12,
         sisoet5,  # social isolation teacher report
         sisoet7, 
         sisoet10,
         sisoet12,
         tadhdem5,  # total ADHD mother report
         tadhdem7,
         tadhdem10,
         tadhdem12,
         tadhdet5,  # total ADHD teacher report
         tadhdet7,
         tadhdet10,
         tadhdet12,
         hyem5,     # hyperactivity ADHD mother report
         hyem7,
         hyem10,
         hyem12, 
         hyet5,     # hyperactivity ADHD teacher report
         hyet7,
         hyet10,
         hyet12,
         inem5,    # inattention ADHD mother report
         inem7,
         inem10,
         inem12,
         inet5,    # inattention ADHD teacher report
         inet7,
         inet10,
         inet12,
         sisoe5,
         sisoe7,
         sisoe10,
         sisoe12
  )

colnames(dat)
```

***

Functions

```{r functions}
# Table of model fit 
table.model.fit <- function(model){
  model.fit <- as.data.frame(t(as.data.frame(model$FIT))) %>%
    dplyr::select(chisq, df, chisq.scaled, cfi.robust, tli.robust, aic, bic, bic2, rmsea.robust, rmsea.ci.lower.robust, rmsea.ci.upper.robust, srmr) #can only be used with "MLR" estimator
  return(model.fit)
}

# Table of regression and correlation (standardised covariance) coefficients
table.model.coef <- function(model, type, constraints){
  if (type == "RICLPM" & constraints == "No"){
    model.coef <- as.tibble(model$PE[c(17:32),]) %>% dplyr::select(-exo, -std.lv, -std.nox)
    return(model.coef)
  } else if(type == "RICLPM" & constraints == "Yes"){
    model.coef <- as.tibble(model$PE[c(17:32),]) %>% dplyr::select(-exo, -label, -std.lv, -std.nox)
    return(model.coef)
  } else if(type == "CLPM" & constraints == "No"){
    model.coef <- as.tibble(model$PE[c(1:16),]) %>% dplyr::select(-exo, -std.lv, -std.nox)
    return(model.coef)
  } else if(type == "CLPM" & constraints == "Yes"){
    model.coef <- as.tibble(model$PE[c(1:16),]) %>% dplyr::select(-exo, -std.lv, -std.nox)
    return(model.coef)
  } else {model.coef <- NULL}
}
```

***

# Normality checks

If the data is non-normal - we should use robust test statistics for all models going forward. Robust standard errors will also be calculated to account for the use of twins in our sample. 

```{r histogram social isolation}
# mother
hist(dat$sisoem5)
hist(dat$sisoem7)
hist(dat$sisoem10)
hist(dat$sisoem12)

# teacher
hist(dat$sisoet5)
hist(dat$sisoet7)
hist(dat$sisoet10)
hist(dat$sisoet12)
```

```{r histogram hyperactivity}
# mother
hist(dat$hyem5)
hist(dat$hyem7)
hist(dat$hyem10)
hist(dat$hyem12)

# mother
hist(dat$hyet5)
hist(dat$hyet7)
hist(dat$hyet10)
hist(dat$hyet12)
```

```{r histogram inattention}
# mother
hist(dat$inem5)
hist(dat$inem7)
hist(dat$inem10)
hist(dat$inem12)

# mother
hist(dat$inet5)
hist(dat$inet7)
hist(dat$inet10)
hist(dat$inet12)
```

***

The below introductory notes have been adapted from [Hamaker, Kuiper, and Grasman, 2015](https://www.researchgate.net/profile/Ellen-Hamaker/publication/274262847_A_Critique_of_the_Cross-Lagged_Panel_Model/links/5b80154c92851c1e122f351f/A-Critique-of-the-Cross-Lagged-Panel-Model.pdf)

The CLPM only accounts for temporal stability through the inclusion of autoregressive parameters. Thus, it is implicitly assumed that every person varies over time around the same means μt and πt, and that there are no trait-like individual differences that endure. this is a rather problematic assumption, as it is difficult to imagine a psychological construct – whether behavioral, cognitive, emotional or psychophysiological – that is not to some extent characterized by stable individual differences (if not for the entire lifespan, then at least for the duration of the study; Hamaker, ). Longitudinal data can actually be thought of as multilevel data, in which occasions are nested within individuals (or other systems, like dyads). When considering this perspective, it becomes clear that we need to separate the *within-person* level from the *between-person* level. 

The random intercepts cross lagged panel model (RI-CLPM) accounts not only for temporal stability, but also for time-invariant, trait-like stability through the inclusion of a random intercept. The cross-lagged parameters indicate the extent to which the two variables influence each other. The cross-lagged relationships pertain to a process that takes place at the within-person level and they are therefore of key interest when the interest is in *reciprocal influences over time within individuals or dyads*. The cross-lagged parameter indicates the extent to which the change in y can be predicted from the individual’s prior deviation from their expected score on the other variable, while controlling for the structural change in y and the prior deviation from one’s expected score on y. 

The CLPM requires only two waves of data, but the RI-CLPM requires at least *three waves* of data, in which case there is 1 degree of freedom (df). If the intervals are of the same size, and if we assume that the effects the variables have on each other remain stable over time, we could decide to constrain the lagged parameters over time, giving us an additional 4 df (i.e., 5 df in total).  If we are not willing to make these assumptions, and we are not sure whether the effect of the time-invariant stability components κi and ωi are equal over time, we may wish to remove the constraint on the factor loadings. This relaxation may especially be of interest when the observations are made further apart in time, and we expect that we are also measuring some structural changes. However, this would imply that κi and ωi no longer represents random intercepts (as in multilevel modeling), but rather represent latent variables or traits (as common in SEM). Even more so, it would imply we need more waves of data to estimate this model. **The gaps in our time waves may be a problem here.**


# All RI-CLPM models 

We conducted RI-CLPM to assess bidirectional associations between social isolation and ADHD symptoms. All models in this script have used sum scores for total ADHD symptoms, hyperactivity symptoms, inattention symptoms, and social isolation. Separate models have been conducted for mother and teacher scores. See the table below for a list of all models and how they are labeled throughout the code. 


| Model               | Description                                                                                                                       | 
| ------------------- | --------------------------------------------------------------------------------------------------------------------------------  |
|   CLPM              |  Cross-lagged panel model without random intercepts using **mother report** and **total ADHD** scores                             |
|   RICLPM            |  Basic RI-CLPM model using **mother report** ratings for AD and SI, and **total ADHD** scores                                     |
|   RICLPM2           |  RICLPM but with fixed autoregressive and cross-lagged relations over time                                                        |
|   RICLPM2a          |  RICLPM but with fixed ADHD autoregressive over time                                                                              |
|   RICLPM2b          |  RICLPM but with fixed social isolation autoregressive over time                                                                  |
|   RICLPM2c          |  RICLPM but with fixed ADHD to social isolation cross-lagged relations over time                                                  |
|   RICLPM2d          |  RICLPM but with fixed all cross-lagged relations over time                                                                       |
|   RICLPM3           |  RICLPM but with constrained grand means for AD and SI                                                                            |
|   RICLPM_hyp        |  Basic RI-CLPM model using **mother report** ratings for AD and SI, and **hyperactivity** scores                                  |
|   RICLPM_hyp2       |  RICLPM_hyp but with fixed autoregressive and cross-lagged relations over time                                                    |
|   RICLPM_hyp2a      |  RICLPM_hyp but with fixed ADHD autoregressive over time                                                                          |
|   RICLPM_hyp2b      |  RICLPM_hyp but with fixed social isolation autoregressive over time                                                              |
|   RICLPM_hyp2c      |  RICLPM_hyp but with fixed ADHD to social isolation cross-lagged relations over time                                              |
|   RICLPM_hyp2d      |  RICLPM_hyp but with fixed all cross-lagged relations over time                                                                   |
|   RICLPM_hyp3       |  RICLPM_hyp but with constrained grand means for AD and SI                                                                        |
|   RICLPM_inat       |  Basic RI-CLPM model using **mother report** ratings for AD and SI, and **inattention** scores                                    |
|   RICLPM_inat2      |  RICLPM_inat but with fixed autoregressive and cross-lagged relations over time                                                   |
|   RICLPM_inat2a     |  RICLPM_inat but with fixed ADHD autoregressive over time                                                                         |
|   RICLPM_inat2b     |  RICLPM_inat but with fixed social isolation autoregressive over time                                                             |
|   RICLPM_inat2c     |  RICLPM_inat but with fixed ADHD to social isolation cross-lagged relations over time                                             |
|   RICLPM_inat2d     |  RICLPM_inat but with fixed all cross-lagged relations over time                                                                  |
|   RICLPM_inat3      |  RICLPM_inat but with constrained grand means for AD and SI                                                                       |
|   CLPMt             |  Cross-lagged panel model without random intercepts using **teacher report** and **total ADHD** scores                            |
|   RICLPMt           |  Basic RI-CLPM model using **teacher report** ratings for AD and SI, and **total ADHD** scores                                    |
|   RICLPMt2          |  RICLPMt but with fixed autoregressive and cross-lagged relations over time                                                       |
|   RICLPMt3          |  RICLPMt but with constrained grand means for AD and SI                                                                           |
|   RICLPMt_hyp       |  Basic RI-CLPM model using **teacher report** ratings for AD and SI, and **hyperactivity** scores                                 |
|   RICLPMt_hyp2      |  RICLPMt_hyp but with fixed autoregressive and cross-lagged relations over time                                                   |
|   RICLPMt_hyp3      |  RICLPMt_hyp but with constrained grand means for AD and SI                                                                       |
|   RICLPMt_inat      |  Basic RI-CLPM model using **teacher report** ratings for AD and SI, and **inattention** scores                                   |
|   RICLPMt_inat2     |  RICLPMt_inat but with fixed autoregressive and cross-lagged relations over time                                                  |
|   RICLPMt_inat3     |  RICLPMt_inat but with constrained grand means for AD and SI                                                                      |

***

# RI-CLPM: Mother report only, total ADHD symptoms

For SEM in [lavaan](https://lavaan.ugent.be/index.html), we use three different formula types: latent variabele definitions (using the `=~` operator, which can be read as is "measured by"), regression formulas (using the `~` operator), and (co)variance formulas (using the `~~` operator). The regression formulas are similar to ordinary formulas in R. The expression `y1 ~~ y2` allows the residual variances of the two observed variables to be correlated. This is sometimes done if it is believed that the two variables have something in common that is not captured by the latent variables. Note that the two expressions `y2 ~~ y4` and `y2 ~~ y6`, can be combined into the expression `y2 ~~ y4 + y6`, because the variable on the left of the `~~` operator (y2) is the same. This is just a shorthand notation. In general, to fix a parameter in a lavaan formula, you need to pre-multiply the corresponding variable in the formula by a numerical value. This is called the pre-multiplication mechanism and will be used for many purposes e.g. `1*x2`. If you wish to fix the correlation (or covariance) between a pair of latent variables to zero, you need to explicity add a covariance-formula for this pair, and fix the parameter to zero.

To specify the RI-CLPM we need four parts: notes from [Mulder's webpage](https://jeroendmulder.github.io/RI-CLPM/lavaan.html).

* A between part, consisting of the random intercepts. It is specified using the =~ command, `RIx =~ 1*x1 1*x2 ...`, where 1* fixes the factor loading to one.

* A within part, consisting of within-unit fluctuations. It is also specified using the =~ command, `wx1 =~ 1*x1; wx2 =~ 1*x2;`

* The lagged regressions between the within-unit components, using `wx2 ~ wx1 wy1; wx3 ~ wx2 wy2; ...`. Here, this is written with two variables on the left hand side of the `~`, this is doing the same thing as calculating just regression paths - but combining DVs that have the same paths leading to them, e.g. `wx2 ~ wx1 + wy1` and `wy2 ~ wx1 + wy1` becomes `wx2+ wy2 ~ wx1 + wy1`. 

* Relevant covariances in both the between and within part. In the within part the components at wave 1, and their residuals at waves 2 and further are correlated within each wave, using `wx1 ~~ wy1; wx2 ~~ wy2;...`. We also need to specify their (residual) variances here using `wx1 ~~ wx1; wx2 ~~ wx2; ...`. For the between part we have to specify the variances and covariance of the random intercepts using `RIx ~~ RIy;`.

Tips for interpreting the output:
- Model Test User Model: which provides a test statistic, degrees of freedom, and a p-value for the model that was specified by the user.
- Model Test Baseline Model: The baseline is a null model, typically in which all of your observed variables are constrained to covary with no other variables (put another way, the covariances are fixed to 0)--just individual variances are estimated.

## Cross-lagged panel model (CLPM)

```{r specify CLPM}
CLPM <- '
  # Estimate the lagged effects between the observed variables.
  tadhdem7 + sisoem7 ~ tadhdem5 + sisoem5
  tadhdem10 + sisoem10 ~ tadhdem7 + sisoem7
  tadhdem12 + sisoem12 ~ tadhdem10 + sisoem10

  # Estimate the covariance between the observed variables at the first wave. 
  tadhdem5 ~~ sisoem5 # Covariance
  
  # Estimate the covariances between the residuals of the observed variables.
  tadhdem7 ~~ sisoem7
  tadhdem10 ~~ sisoem10
  tadhdem12 ~~ sisoem12
  
  # Estimate the (residual) variance of the observed variables.
  tadhdem5 ~~ tadhdem5 # Variances
  sisoem5 ~~ sisoem5 
  tadhdem7 ~~ tadhdem7 # Residual variances
  sisoem7 ~~ sisoem7 
  tadhdem10 ~~ tadhdem10 
  sisoem10 ~~ sisoem10 
  tadhdem12 ~~ tadhdem12 
  sisoem12 ~~ sisoem12 
'
```

```{r fit CLPM}
CLPM.fit <- lavaan(CLPM, 
                   data = dat, 
                   missing = 'ML',
                   meanstructure = TRUE, 
                   int.ov.free = TRUE,
                   se = "robust",
                   estimator = "MLR" #maximum likelihood with robust (Huber-White) standard errors and a scaled (Yuan-Bentler) and robust test statistic 
                   ) 

CLPM.fit.summary <- summary(CLPM.fit, 
                            fit.measures = TRUE,
                            standardized = TRUE)

#Table of model fit 
CLPM.fit.summary.fit <- table.model.fit(CLPM.fit.summary)
CLPM.fit.summary.fit
#Table of regression coefficients and covariances (concurrent associations)
CLPM.fit.summary.reg <- table.model.coef(model = CLPM.fit.summary, type = "CLPM", constraints = "No")
CLPM.fit.summary.reg
```

## The basic RI-CLPM model (RICLPM)
The code for specifying the basic RI-CLPM is given below.

```{r specify RICLPM, class.source = 'fold-show'}
RICLPM <- '
  # Create between components (random intercepts treated as factors here)
  RIad =~ 1*tadhdem5 + 1*tadhdem7 + 1*tadhdem10 + 1*tadhdem12 #x
  RIsi =~ 1*sisoem5 + 1*sisoem7 + 1*sisoem10 + 1*sisoem12 #y

  # Create within-person centered variables
  wad5 =~ 1*tadhdem5
  wad7 =~ 1*tadhdem7
  wad10 =~ 1*tadhdem10 
  wad12 =~ 1*tadhdem12
  wsi5 =~ 1*sisoem5
  wsi7 =~ 1*sisoem7
  wsi10 =~ 1*sisoem10
  wsi12 =~ 1*sisoem12
  
  # Estimate the lagged effects between the within-person centered variables
  wad7 + wsi7 ~ wad5 + wsi5
  wad10 + wsi10 ~ wad7 + wsi7
  wad12 + wsi12 ~ wad10 + wsi10
  
  # Estimate the covariance between the within-person centered variables at the first wave
  wad5 ~~ wsi5 # Covariance
  
  # Estimate the covariances between the residuals of the within-person centered variables (the innovations)
  wad7 ~~ wsi7
  wad10 ~~ wsi10
  wad12 ~~ wsi12
  
  # Estimate the variance and covariance of the random intercepts
  RIad ~~ RIad
  RIsi ~~ RIsi
  RIad ~~ RIsi
  
  # Estimate the (residual) variance of the within-person centered variables.
  wad5 ~~ wad5 # Variances
  wsi5 ~~ wsi5 
  wad7 ~~ wad7 # Residual variances
  wsi7 ~~ wsi7 
  wad10 ~~ wad10 
  wsi10 ~~ wsi10 
  wad12 ~~ wad12 
  wsi12 ~~ wsi12
'
```

```{r fit RICLPM}
RICLPM.fit <- lavaan(RICLPM,               # model
                     data = dat,           # data
                     missing = 'ML',       # how to handle missing data 
                     meanstructure = TRUE, # adds intercepts/means to the model for both observed and latent variables
                     se = "robust",        # robust standard errors
                     int.ov.free = TRUE,   # if FALSE, the intercepts of the observed variables are fixed to zero
                     estimator = "MLR" #maximum likelihood with robust (Huber-White) standard errors and a scaled (Yuan-Bentler) and robust test statistic
)

RICLPM.fit.summary <- summary(RICLPM.fit, 
                              fit.measures = TRUE,
                              standardized = TRUE)

#Table of model fit 
RICLPM.fit.summary.fit <- table.model.fit(RICLPM.fit.summary)
#Table of regression coefficients and covariances (concurrent associations)
RICLPM.fit.summary.reg <- table.model.coef(model = RICLPM.fit.summary, type = "RICLPM", constraints = "No")
```

## Comparison between CLPM and RI-CLPM

The use of the chi-square difference test is wide-spread in the SEM community to test constraints on parameters. However, when constraints are placed on the bound of the parameter space, we should use the chi-bar-square test (χ¯2-test) (Stoel et al. 2006). For example, if we constrain the variances of all random intercepts (and their covariance) in the RI-CLPM to zero, we obtain a model that is nested under the RI-CLPM, and that is statistically equivalent to the traditional cross-lagged panel model (CLPM). Below you can find R code for performing the chi-bar-square test (code by Rebecca M. Kuiper) for comparing these two models. It involves

  1. fitting both the RI-CLPM (RICLPM.fit) and CLPM (CLPM.fit); (already done - this is the non-constrained RI-CLPM - RICLPM.fit)
  2. extracting the covariance matrix of the random intercepts;
  3. extracting the χ2 and degrees of freedom of both models; and
  4. performing the χ¯2-test using the ChiBarSq.DiffTest package (Kuiper 2020).
  
```{r chi square bar test CLPM RICLPM}
# # Step 2: check which indices you need to get the covariance matrix of the random intercepts. 
# vcov(RICLPM.fit) # Full covariance matrix
# indices <- c(17, 18) # From the matrix above, you need the variance of each intercept. Here, the 17th and the 18th indices regard the random intercepts. 
# number <- length(indices) # Number of random intercepts
# cov.matrix <- vcov(RICLPM.fit)[indices, indices] # Extract full covariance matrix of the random intercepts
# 
# # Step 3: call Chi square and degrees of freedom from each model
# Chi2_CLPM <- summary(CLPM.fit, fit.measures = TRUE)[1]$FIT[c("chisq")] # Extract chi-square value of CLPM
# Chi2_RICLPM <- summary(RICLPM.fit, fit.measures = TRUE)[1]$FIT[c("chisq")] # Extract chi-square value of RI-CLPM
# 
# df_CLPM <- summary(CLPM.fit, fit.measures = TRUE)[1]$FIT[c("df")] # Extract df of CLPM
# df_RICLPM <- summary(RICLPM.fit, fit.measures = TRUE)[1]$FIT[c("df")] # Extract df of RI-CLPM
# 
# # Step 4: run function to do chi-bar-square test (and also obtain Chi-bar-square weigths)
# ChiBar2DiffTest <- ChiBarSq.DiffTest(number, cov.matrix, Chi2_CLPM, Chi2_RICLPM, df_CLPM, df_RICLPM)
# ChiBar2DiffTest
# ChiBar2DiffTest$p_value
```

## RI-CLPM Constraints over time {.tabset .tabset-fade}

Imposing constraints to the model can be achieved through **pre-multiplication**. It means that we have to prepend the number that we want to fix the parameter to, and an asterisk, to the parameter in the model specification. For example, `F =~ 0*x1` fixes the factor loading of item `x1` to factor `F` to 0. Using pre-multiplication we can also constrain parameters to be the same by giving them the same label. Below we specify an RI-CLPM with the following constraints: 

1) fixed auto-regressive and cross-lagged relations over time, `wx2 ~ a*wx1 + b*wy1; ...`
2) time-invariant (residual) (co-)variances in the within-person part `wx2 ~~ cov*wy2; ...`, `wx2 ~~ vx*wx2; ...`, and `wy2 ~~ vy*wy2; ...`
3) constrained grand means over time, `x1 + ... ~ mx*1` and `y1 + ... ~ my*1`

### Fixed autoregressive and cross-lagged relations over time (RICLPM2)

a = lag in ad
b = lag in si
c = cross lag ad->si
d = cross lag si->ad

From [Mulder and Hamaker (2021)](https://www.tandfonline.com/doi/full/10.1080/10705511.2020.1784738): We may consider testing if the lagged regression coefficients are time-invariant. This can be done by comparing the fit of a model with constrained regression coefficients (over time), with the fit of a model where these parameters are freely estimated (i.e., the unconstrained model). If this chi-square difference test is non-significant, this implies the constraints are tenable and the dynamics of the process are time-invariant. If the constraints are not tenable, this could be indicative of some kind of developmental process taking place during the time span covered by the study.
In this context, it is important to realize that the lagged regression coefficients depend critically on the time interval between the repeated measures. Hence, constraining the lagged parameters to be invariant across consecutive waves only makes sense when the time interval between the occasions is *(approximately) equal*. 

#### Constraining all lag and crosslag parameters at once (RICLPM2)

```{r specify RICLPM2}
RICLPM2 <- '
  # Create between components (random intercepts treated as factors here)
  RIad =~ 1*tadhdem5 + 1*tadhdem7 + 1*tadhdem10 + 1*tadhdem12 #x
  RIsi =~ 1*sisoem5 + 1*sisoem7 + 1*sisoem10 + 1*sisoem12 #y

  # Create within-person centered variables
  wad5 =~ 1*tadhdem5
  wad7 =~ 1*tadhdem7
  wad10 =~ 1*tadhdem10 
  wad12 =~ 1*tadhdem12
  wsi5 =~ 1*sisoem5
  wsi7 =~ 1*sisoem7
  wsi10 =~ 1*sisoem10
  wsi12 =~ 1*sisoem12
  
  
  # Constrained lagged effects between the within-person centered variables. 
  wad7 ~ a*wad5 + d*wsi5 
  wsi7 ~ c*wad5 + b*wsi5
  
  wad10 ~ a*wad7 + d*wsi7
  wsi10 ~ c*wad7 + b*wsi7
  
  wad12 ~ a*wad10 + d*wsi10
  wsi12 ~ c*wad10 + b*wsi10
  
  # Estimate the covariance between the within-person centered variables at the first wave
  wad5 ~~ wsi5 # Covariance
  
  # Estimate the covariances between the residuals of the within-person centered variables (the innovations)
  wad7 ~~ wsi7
  wad10 ~~ wsi10
  wad12 ~~ wsi12
  
  # Estimate the variance and covariance of the random intercepts
  RIad ~~ RIad
  RIsi ~~ RIsi
  RIad ~~ RIsi
  
  # Estimate the (residual) variance of the within-person centered variables.
  wad5 ~~ wad5 # Variances
  wsi5 ~~ wsi5 
  wad7 ~~ wad7 # Residual variances
  wsi7 ~~ wsi7 
  wad10 ~~ wad10 
  wsi10 ~~ wsi10 
  wad12 ~~ wad12 
  wsi12 ~~ wsi12
'
```

```{r fit RICLPM2}
RICLPM2.fit <- lavaan(RICLPM2, 
                      data = dat, 
                      missing = 'ML', 
                      meanstructure = TRUE, 
                      int.ov.free = TRUE,
                      se = "robust",
                      estimator = "MLR" #maximum likelihood with robust (Huber-White) standard errors and a scaled (Yuan-Bentler) and robust test statistic
                      ) 
```

```{r LRT for RICLPM and RICLPM2}
lavTestLRT(RICLPM.fit, RICLPM2.fit, method = "satorra.bentler.2010")
```

RICLPM2 is significantly *worse* fit compared to RICLPM according to the Chi square test.

#### Constraining lag in ADHD parameter (RICLPM2a)

Now we will test the constraints based on the following steps recommended by Curran and Bauer:
1. Estimate the baseline model where all parameters are freely estimated (The basic RI-CLPM)
2. Impose equlaity constraints on one set of stabilities (within variable lags). Use LRT tests to see if he fit becomes significantly worse. 
3. Impose equality constraints on the next set of stabilities. Compare this to wither model 1 (baseline/basic) if step 2 was significant. Compare to model 2 if step 2 was non-significant. *Non-sig LRT = not significantly worse fit*
4. Repeat for first set of cross-lags
5. Repeat for second set of cross-lags

```{r specify RICLPM2a}
RICLPM2a <- '
  # Create between components (random intercepts treated as factors here)
  RIad =~ 1*tadhdem5 + 1*tadhdem7 + 1*tadhdem10 + 1*tadhdem12 #x
  RIsi =~ 1*sisoem5 + 1*sisoem7 + 1*sisoem10 + 1*sisoem12 #y

  # Create within-person centered variables
  wad5 =~ 1*tadhdem5
  wad7 =~ 1*tadhdem7
  wad10 =~ 1*tadhdem10 
  wad12 =~ 1*tadhdem12
  wsi5 =~ 1*sisoem5
  wsi7 =~ 1*sisoem7
  wsi10 =~ 1*sisoem10
  wsi12 =~ 1*sisoem12
  
  
  # Constrained lagged effects between the within-person centered variables. 
  wad7 ~ a*wad5 + wsi5 
  wsi7 ~ wad5 + wsi5
  
  wad10 ~ a*wad7 + wsi7
  wsi10 ~ wad7 + wsi7
  
  wad12 ~ a*wad10 + wsi10
  wsi12 ~ wad10 + wsi10
  
  # Estimate the covariance between the within-person centered variables at the first wave
  wad5 ~~ wsi5 # Covariance
  
  # Estimate the covariances between the residuals of the within-person centered variables (the innovations)
  wad7 ~~ wsi7
  wad10 ~~ wsi10
  wad12 ~~ wsi12
  
  # Estimate the variance and covariance of the random intercepts
  RIad ~~ RIad
  RIsi ~~ RIsi
  RIad ~~ RIsi
  
  # Estimate the (residual) variance of the within-person centered variables.
  wad5 ~~ wad5 # Variances
  wsi5 ~~ wsi5 
  wad7 ~~ wad7 # Residual variances
  wsi7 ~~ wsi7 
  wad10 ~~ wad10 
  wsi10 ~~ wsi10 
  wad12 ~~ wad12 
  wsi12 ~~ wsi12
'
```

```{r fit RICLPM2a}
RICLPM2a.fit <- lavaan(RICLPM2a, 
                      data = dat, 
                      missing = 'ML', 
                      meanstructure = TRUE, 
                      int.ov.free = TRUE,
                      se = "robust",
                      estimator = "MLR" #maximum likelihood with robust (Huber-White) standard errors and a scaled (Yuan-Bentler) and robust test statistic
                      ) 
```

```{r LRT for RICLPM and RICLPM2a}
lavTestLRT(RICLPM.fit, RICLPM2a.fit, method = "satorra.bentler.2010")
```

RICLPM2a is significantly *worse* fit compared to RICLPM. Now we will constrain the other set of stability coefficients: lag in si (b)

#### Constraining lag in social isolation parameter (RICLPM2b)

```{r specify RICLPM2b}
RICLPM2b <- '
  # Create between components (random intercepts treated as factors here)
  RIad =~ 1*tadhdem5 + 1*tadhdem7 + 1*tadhdem10 + 1*tadhdem12 #x
  RIsi =~ 1*sisoem5 + 1*sisoem7 + 1*sisoem10 + 1*sisoem12 #y

  # Create within-person centered variables
  wad5 =~ 1*tadhdem5
  wad7 =~ 1*tadhdem7
  wad10 =~ 1*tadhdem10 
  wad12 =~ 1*tadhdem12
  wsi5 =~ 1*sisoem5
  wsi7 =~ 1*sisoem7
  wsi10 =~ 1*sisoem10
  wsi12 =~ 1*sisoem12
  
  
  # Constrained lagged effects between the within-person centered variables. 
  wad7 ~ wad5 + wsi5 
  wsi7 ~ wad5 + b*wsi5
  
  wad10 ~ wad7 + wsi7
  wsi10 ~ wad7 + b*wsi7
  
  wad12 ~ wad10 + wsi10
  wsi12 ~ wad10 + b*wsi10
  
  # Estimate the covariance between the within-person centered variables at the first wave
  wad5 ~~ wsi5 # Covariance
  
  # Estimate the covariances between the residuals of the within-person centered variables (the innovations)
  wad7 ~~ wsi7
  wad10 ~~ wsi10
  wad12 ~~ wsi12
  
  # Estimate the variance and covariance of the random intercepts
  RIad ~~ RIad
  RIsi ~~ RIsi
  RIad ~~ RIsi
  
  # Estimate the (residual) variance of the within-person centered variables.
  wad5 ~~ wad5 # Variances
  wsi5 ~~ wsi5 
  wad7 ~~ wad7 # Residual variances
  wsi7 ~~ wsi7 
  wad10 ~~ wad10 
  wsi10 ~~ wsi10 
  wad12 ~~ wad12 
  wsi12 ~~ wsi12
'
```

```{r fit RICLPM2b}
RICLPM2b.fit <- lavaan(RICLPM2b, 
                      data = dat, 
                      missing = 'ML', 
                      meanstructure = TRUE, 
                      int.ov.free = TRUE,
                      se = "robust",
                      estimator = "MLR" #maximum likelihood with robust (Huber-White) standard errors and a scaled (Yuan-Bentler) and robust test statistic
                      ) 
```

```{r LRT for RICLPM and RICLPM2b}
lavTestLRT(RICLPM.fit, RICLPM2b.fit, method = "satorra.bentler.2010")
```

RICLPM2b is significantly *worse* fit compared to RICLPM. Now we will constrain the other set of stability coefficients: cross lag ad->si (c)

#### Constraining cross lag in ADHD to social isolation parameter (RICLPM2c)

```{r specify RICLPM2c}
RICLPM2c <- '
  # Create between components (random intercepts treated as factors here)
  RIad =~ 1*tadhdem5 + 1*tadhdem7 + 1*tadhdem10 + 1*tadhdem12 #x
  RIsi =~ 1*sisoem5 + 1*sisoem7 + 1*sisoem10 + 1*sisoem12 #y

  # Create within-person centered variables
  wad5 =~ 1*tadhdem5
  wad7 =~ 1*tadhdem7
  wad10 =~ 1*tadhdem10 
  wad12 =~ 1*tadhdem12
  wsi5 =~ 1*sisoem5
  wsi7 =~ 1*sisoem7
  wsi10 =~ 1*sisoem10
  wsi12 =~ 1*sisoem12
  
  # Constrained lagged effects between the within-person centered variables. 
  wad7 ~ wad5 + wsi5 
  wsi7 ~ c*wad5 + wsi5
  
  wad10 ~ wad7 + wsi7
  wsi10 ~ c*wad7 + wsi7
  
  wad12 ~ wad10 + wsi10
  wsi12 ~ c*wad10 + wsi10
  
  # Estimate the covariance between the within-person centered variables at the first wave
  wad5 ~~ wsi5 # Covariance
  
  # Estimate the covariances between the residuals of the within-person centered variables (the innovations)
  wad7 ~~ wsi7
  wad10 ~~ wsi10
  wad12 ~~ wsi12
  
  # Estimate the variance and covariance of the random intercepts
  RIad ~~ RIad
  RIsi ~~ RIsi
  RIad ~~ RIsi
  
  # Estimate the (residual) variance of the within-person centered variables.
  wad5 ~~ wad5 # Variances
  wsi5 ~~ wsi5 
  wad7 ~~ wad7 # Residual variances
  wsi7 ~~ wsi7 
  wad10 ~~ wad10 
  wsi10 ~~ wsi10 
  wad12 ~~ wad12 
  wsi12 ~~ wsi12
'
```

```{r fit RICLPM2c}
RICLPM2c.fit <- lavaan(RICLPM2c, 
                      data = dat, 
                      missing = 'ML', 
                      meanstructure = TRUE, 
                      int.ov.free = TRUE,
                      se = "robust",
                      estimator = "MLR" #maximum likelihood with robust (Huber-White) standard errors and a scaled (Yuan-Bentler) and robust test statistic
                      ) 
```

```{r LRT for RICLPM and RICLPM2c}
lavTestLRT(RICLPM.fit, RICLPM2c.fit, method = "satorra.bentler.2010")
```

RICLPM2c is *NOT* significantly worse fit compared to RICLPM (p = 0.9616). Now we will constrain the other set of stability coefficients: cross lag si->ad (d) *but comparing that model to this one*. 

#### Constraining cross lag in social isolation to ADHD parameter (RICLPM2d) - **compared to RICLPM2c** 

Thus, constraining both sets of lags are tested in this model. 

```{r specify RICLPM2d}
RICLPM2d <- '
  # Create between components (random intercepts treated as factors here)
  RIad =~ 1*tadhdem5 + 1*tadhdem7 + 1*tadhdem10 + 1*tadhdem12 #x
  RIsi =~ 1*sisoem5 + 1*sisoem7 + 1*sisoem10 + 1*sisoem12 #y

  # Create within-person centered variables
  wad5 =~ 1*tadhdem5
  wad7 =~ 1*tadhdem7
  wad10 =~ 1*tadhdem10 
  wad12 =~ 1*tadhdem12
  wsi5 =~ 1*sisoem5
  wsi7 =~ 1*sisoem7
  wsi10 =~ 1*sisoem10
  wsi12 =~ 1*sisoem12
  
  
  # Constrained lagged effects between the within-person centered variables. 
  wad7 ~ wad5 + d*wsi5 
  wsi7 ~ c*wad5 + wsi5
  
  wad10 ~ wad7 + d*wsi7
  wsi10 ~ c*wad7 + wsi7
  
  wad12 ~ wad10 + d*wsi10
  wsi12 ~ c*wad10 + wsi10
  
  # Estimate the covariance between the within-person centered variables at the first wave
  wad5 ~~ wsi5 # Covariance
  
  # Estimate the covariances between the residuals of the within-person centered variables (the innovations)
  wad7 ~~ wsi7
  wad10 ~~ wsi10
  wad12 ~~ wsi12
  
  # Estimate the variance and covariance of the random intercepts
  RIad ~~ RIad
  RIsi ~~ RIsi
  RIad ~~ RIsi
  
  # Estimate the (residual) variance of the within-person centered variables.
  wad5 ~~ wad5 # Variances
  wsi5 ~~ wsi5 
  wad7 ~~ wad7 # Residual variances
  wsi7 ~~ wsi7 
  wad10 ~~ wad10 
  wsi10 ~~ wsi10 
  wad12 ~~ wad12 
  wsi12 ~~ wsi12
'
```

```{r fit RICLPM2d}
RICLPM2d.fit <- lavaan(RICLPM2d, 
                      data = dat, 
                      missing = 'ML', 
                      meanstructure = TRUE, 
                      int.ov.free = TRUE,
                      se = "robust",
                      estimator = "MLR" #maximum likelihood with robust (Huber-White) standard errors and a scaled (Yuan-Bentler) and robust test statistic
                      ) 

RICLPM2d.fit.summary <- summary(RICLPM2d.fit, 
                                fit.measures = TRUE,
                                standardized = TRUE)

#Table of model fit 
RICLPM2d.fit.summary.fit <- table.model.fit(RICLPM2d.fit.summary)
RICLPM2d.fit.summary.fit
#Table of regression coefficients and covariances (concurrent associations)
RICLPM2d.fit.summary.reg <- table.model.coef(model = RICLPM2d.fit.summary, type = "RICLPM", constraints = "Yes")
```

```{r LRT for RICLPM2c and RICLPM2d}
lavTestLRT(RICLPM2c.fit, RICLPM2d.fit, method = "satorra.bentler.2010") #comparing to other best fitting model
```

RICLPM2d is *NOT* significantly worse fit compared to RICLPM2c (p = 0.4225). 

**The best fitting model (mother report and total ADHD symptoms) is currently RICLPM2d where cross-lags are constrained to be equal across time.**

***

### Constrained grand means (RICLPM3)

The grand means are the means over all units per occasion. These grand means may be time-varying, or may be fixed to be invariant over time. 

```{r specify RICLPM3}
RICLPM3 <- '
  # Create between components (random intercepts treated as factors here)
  RIad =~ 1*tadhdem5 + 1*tadhdem7 + 1*tadhdem10 + 1*tadhdem12 #x
  RIsi =~ 1*sisoem5 + 1*sisoem7 + 1*sisoem10 + 1*sisoem12 #y

  # Create within-person centered variables
  wad5 =~ 1*tadhdem5
  wad7 =~ 1*tadhdem7
  wad10 =~ 1*tadhdem10 
  wad12 =~ 1*tadhdem12
  wsi5 =~ 1*sisoem5
  wsi7 =~ 1*sisoem7
  wsi10 =~ 1*sisoem10
  wsi12 =~ 1*sisoem12
  
  # Constrained lagged effects between the within-person centered variables. 
  wad7 ~ wad5 + wsi5 
  wsi7 ~ wad5 + wsi5
  wad10 ~ wad7 + wsi7
  wsi10 ~ wad7 + wsi7
  wad12 ~ wad10 + wsi10
  wsi12 ~ wad10 + wsi10
  
  # Estimate the covariance between the within-person centered variables at the first wave
  wad5 ~~ wsi5 # Covariance
  
  # Estimate the covariances between the residuals of the within-person centered variables (the innovations)
  wad7 ~~ wsi7
  wad10 ~~ wsi10
  wad12 ~~ wsi12
  
  # Estimate the variance and covariance of the random intercepts
  RIad ~~ RIad
  RIsi ~~ RIsi
  RIad ~~ RIsi
  
  # Estimate the (residual) variance of the within-person centered variables
  wad5 ~~ wad5 # Variances
  wsi5 ~~ wsi5 
  wad7 ~~ wad7 # Residual variances
  wsi7 ~~ wsi7 
  wad10 ~~ wad10 
  wsi10 ~~ wsi10 
  wad12 ~~ wad12 
  wsi12 ~~ wsi12
  
  # Constrain the grand means over time
  tadhdem5 + tadhdem7 + tadhdem10 + tadhdem12 ~ mad*1
  sisoem5 + sisoem7 + sisoem10 + sisoem12 ~ msi*1
'
```

```{r fit RICLPM3}
RICLPM3.fit <- lavaan(RICLPM3, 
                      data = dat, 
                      missing = 'ML', 
                      meanstructure = TRUE, 
                      int.ov.free = TRUE,
                      se = "robust",
                      estimator = "MLR" #maximum likelihood with robust (Huber-White) standard errors and a scaled (Yuan-Bentler) and robust test statistic
                      ) 
```

```{r LRT for RICLPM and RICLPM3}
lavTestLRT(RICLPM.fit, RICLPM3.fit, method = "satorra.bentler.2010") 
```

If the grand means cannot be constrained to be invariant over time, this implies that on average there is some change in this variable over time, which may reflect some occasion-specific effect, or a developmental trend. By allowing the means to freely vary over time, we account for such average changes over time.

*This makes sense as we have shown variation in social isolation over time and other literature has shown variation in ADHD over time.*

***

# RI-CLPM: Mother report only, Hyperactive/Impulsive ADHD symptoms

## The basic RI-CLPM model (RICLPM_hyp)
The code for specifying the basic RI-CLPM is given below.

```{r specify RICLPM_hyp}
RICLPM_hyp <- '
  # Create between components (random intercepts treated as factors here)
  RIad =~ 1*hyem5 + 1*hyem7 + 1*hyem10 + 1*hyem12 #x
  RIsi =~ 1*sisoem5 + 1*sisoem7 + 1*sisoem10 + 1*sisoem12 #y

  # Create within-person centered variables
  wad5 =~ 1*hyem5
  wad7 =~ 1*hyem7
  wad10 =~ 1*hyem10 
  wad12 =~ 1*hyem12
  wsi5 =~ 1*sisoem5
  wsi7 =~ 1*sisoem7
  wsi10 =~ 1*sisoem10
  wsi12 =~ 1*sisoem12
  
  # Estimate the lagged effects between the within-person centered variables
  wad7 + wsi7 ~ wad5 + wsi5
  wad10 + wsi10 ~ wad7 + wsi7
  wad12 + wsi12 ~ wad10 + wsi10
  
  # Estimate the covariance between the within-person centered variables at the first wave
  wad5 ~~ wsi5 # Covariance
  
  # Estimate the covariances between the residuals of the within-person centered variables (the innovations)
  wad7 ~~ wsi7
  wad10 ~~ wsi10
  wad12 ~~ wsi12
  
  # Estimate the variance and covariance of the random intercepts
  RIad ~~ RIad
  RIsi ~~ RIsi
  RIad ~~ RIsi
  
  # Estimate the (residual) variance of the within-person centered variables.
  wad5 ~~ wad5 # Variances
  wsi5 ~~ wsi5 
  wad7 ~~ wad7 # Residual variances
  wsi7 ~~ wsi7 
  wad10 ~~ wad10 
  wsi10 ~~ wsi10 
  wad12 ~~ wad12 
  wsi12 ~~ wsi12
'
```

```{r fit RICLPM_hyp}
RICLPM_hyp.fit <- lavaan(RICLPM_hyp,       # model
                     data = dat,           # data
                     missing = 'ML',       # how to handle missing data 
                     meanstructure = TRUE, # adds intercepts/means to the model for both observed and latent variables
                     se = "robust",        # robust standard errors
                     int.ov.free = TRUE,   # if FALSE, the intercepts of the observed variables are fixed to zero
                     estimator = "MLR" #maximum likelihood with robust (Huber-White) standard errors and a scaled (Yuan-Bentler) and robust test statistic
)

RICLPM_hyp.fit.summary <- summary(RICLPM_hyp.fit, 
                                  fit.measures = TRUE,
                                  standardized = TRUE)

#Table of model fit 
RICLPM_hyp.fit.summary.fit <- table.model.fit(RICLPM_hyp.fit.summary)
#Table of regression coefficients and covariances (concurrent associations)
RICLPM_hyp.fit.summary.reg <- table.model.coef(model = RICLPM_hyp.fit.summary, type = "RICLPM", constraints = "No")
```

## RI-CLPM Constraints over time {.tabset .tabset-fade}

Imposing constraints to the model can be achieved through **pre-multiplication**. It means that we have to prepend the number that we want to fix the parameter to, and an asterisk, to the parameter in the model specification. For example, `F =~ 0*x1` fixes the factor loading of item `x1` to factor `F` to 0. Using pre-multiplication we can also constrain parameters to be the same by giving them the same label. Below we specify an RI-CLPM with the following constraints: 

1) fixed auto-regressive and cross-lagged relations over time, `wx2 ~ a*wx1 + b*wy1; ...`
2) time-invariant (residual) (co-)variances in the within-person part `wx2 ~~ cov*wy2; ...`, `wx2 ~~ vx*wx2; ...`, and `wy2 ~~ vy*wy2; ...`
3) constrained grand means over time, `x1 + ... ~ mx*1` and `y1 + ... ~ my*1`

Here we will use the Chi square significance value to indicate a significantly worse fitting model. 

You can also use fit indices for large sample sizes as below, we will use this method when looking at invariance testing in the measirement models:
 - RMSEA difference = increase smaller than 0.01 is acceptable
 - CFI difference = decrease smaller than 0.01 is acceptable
 - TLI difference = decrease smaller than 0.01 is acceptable

### Fixed autoregressive and cross-lagged relations over time (RICLPM_hyp2)

a = lag in ad
b = lag in si
c = cross lag ad->si
d = cross lag si->ad

#### Constraining all lag and crosslag parameters at once (RICLPM2)

```{r specify RICLPM_hyp2}
RICLPM_hyp2 <- '
  # Create between components (random intercepts treated as factors here)
  RIad =~ 1*hyem5 + 1*hyem7 + 1*hyem10 + 1*hyem12 #x
  RIsi =~ 1*sisoem5 + 1*sisoem7 + 1*sisoem10 + 1*sisoem12 #y

  # Create within-person centered variables
  wad5 =~ 1*hyem5
  wad7 =~ 1*hyem7
  wad10 =~ 1*hyem10 
  wad12 =~ 1*hyem12
  wsi5 =~ 1*sisoem5
  wsi7 =~ 1*sisoem7
  wsi10 =~ 1*sisoem10
  wsi12 =~ 1*sisoem12
  
  
  # Constrained lagged effects between the within-person centered variables. 
  wad7 ~ a*wad5 + d*wsi5 
  wsi7 ~ c*wad5 + b*wsi5
  
  wad10 ~ a*wad7 + d*wsi7
  wsi10 ~ c*wad7 + b*wsi7
  
  wad12 ~ a*wad10 + d*wsi10
  wsi12 ~ c*wad10 + b*wsi10
  
  # Estimate the covariance between the within-person centered variables at the first wave
  wad5 ~~ wsi5 # Covariance
  
  # Estimate the covariances between the residuals of the within-person centered variables (the innovations)
  wad7 ~~ wsi7
  wad10 ~~ wsi10
  wad12 ~~ wsi12
  
  # Estimate the variance and covariance of the random intercepts
  RIad ~~ RIad
  RIsi ~~ RIsi
  RIad ~~ RIsi
  
  # Estimate the (residual) variance of the within-person centered variables.
  wad5 ~~ wad5 # Variances
  wsi5 ~~ wsi5 
  wad7 ~~ wad7 # Residual variances
  wsi7 ~~ wsi7 
  wad10 ~~ wad10 
  wsi10 ~~ wsi10 
  wad12 ~~ wad12 
  wsi12 ~~ wsi12
'
```

```{r fit RICLPM_hyp2}
RICLPM_hyp2.fit <- lavaan(RICLPM_hyp2, 
                      data = dat, 
                      missing = 'ML', 
                      meanstructure = TRUE, 
                      int.ov.free = TRUE,
                      se = "robust",
                      estimator = "MLR" #maximum likelihood with robust (Huber-White) standard errors and a scaled (Yuan-Bentler) and robust test statistic
                      ) 
```

```{r LRT for RICLPM_hyp and RICLPM_hyp2}
lavTestLRT(RICLPM_hyp.fit, RICLPM_hyp2.fit, method = "satorra.bentler.2010")
```

RICLPM_hyp2 is significantly *worse* fit compared to RICLPM_hyp. 

#### Constraining lag in ADHD parameter (RICLPM_hyp2a)

Now we will test the constraints based on the following steps recommended by Curran and Bauer:
1. Estimate the baseline model where all parameters are freely estimated (The basic RI-CLPM)
2. Impose equlaity constraints on one set of stabilities (within variable lags). Use LRT tests to see if he fit becomes significantly worse. 
3. Impose equality constraints on the next set of stabilities. Compare this to wither model 1 (baseline/basic) if step 2 was significant. Compare to model 2 if step 2 was non-significant. *Non-sig LRT = not significantly worse fit*
4. Repeat for first set of cross-lags
5. Repeat for second set of cross-lags

```{r specify RICLPM_hyp2a}
RICLPM_hyp2a <- '
  # Create between components (random intercepts treated as factors here)
  RIad =~ 1*hyem5 + 1*hyem7 + 1*hyem10 + 1*hyem12 #x
  RIsi =~ 1*sisoem5 + 1*sisoem7 + 1*sisoem10 + 1*sisoem12 #y

  # Create within-person centered variables
  wad5 =~ 1*hyem5
  wad7 =~ 1*hyem7
  wad10 =~ 1*hyem10 
  wad12 =~ 1*hyem12
  wsi5 =~ 1*sisoem5
  wsi7 =~ 1*sisoem7
  wsi10 =~ 1*sisoem10
  wsi12 =~ 1*sisoem12
  
  
  # Constrained lagged effects between the within-person centered variables. 
  wad7 ~ a*wad5 + wsi5 
  wsi7 ~ wad5 + wsi5
  
  wad10 ~ a*wad7 + wsi7
  wsi10 ~ wad7 + wsi7
  
  wad12 ~ a*wad10 + wsi10
  wsi12 ~ wad10 + wsi10
  
  # Estimate the covariance between the within-person centered variables at the first wave
  wad5 ~~ wsi5 # Covariance
  
  # Estimate the covariances between the residuals of the within-person centered variables (the innovations)
  wad7 ~~ wsi7
  wad10 ~~ wsi10
  wad12 ~~ wsi12
  
  # Estimate the variance and covariance of the random intercepts
  RIad ~~ RIad
  RIsi ~~ RIsi
  RIad ~~ RIsi
  
  # Estimate the (residual) variance of the within-person centered variables.
  wad5 ~~ wad5 # Variances
  wsi5 ~~ wsi5 
  wad7 ~~ wad7 # Residual variances
  wsi7 ~~ wsi7 
  wad10 ~~ wad10 
  wsi10 ~~ wsi10 
  wad12 ~~ wad12 
  wsi12 ~~ wsi12
'
```

```{r fit RICLPM_hyp2a}
RICLPM_hyp2a.fit <- lavaan(RICLPM_hyp2a, 
                      data = dat, 
                      missing = 'ML', 
                      meanstructure = TRUE, 
                      int.ov.free = TRUE,
                      se = "robust",
                      estimator = "MLR" #maximum likelihood with robust (Huber-White) standard errors and a scaled (Yuan-Bentler) and robust test statistic
                      ) 
```

```{r LRT for RICLPM_hyp and RICLPM_hyp2a}
lavTestLRT(RICLPM_hyp.fit, RICLPM_hyp2a.fit, method = "satorra.bentler.2010")
```

RICLPM_hyp2a is significantly *worse* fit compared to RICLPM_hyp. Now we will constrain the other set of stability coefficients: lag in si (b)

#### Constraining lag in social isolation parameter (RICLPM_hyp2b)

```{r specify RICLPM_hyp2b}
RICLPM_hyp2b <- '
  # Create between components (random intercepts treated as factors here)
  RIad =~ 1*hyem5 + 1*hyem7 + 1*hyem10 + 1*hyem12 #x
  RIsi =~ 1*sisoem5 + 1*sisoem7 + 1*sisoem10 + 1*sisoem12 #y

  # Create within-person centered variables
  wad5 =~ 1*hyem5
  wad7 =~ 1*hyem7
  wad10 =~ 1*hyem10 
  wad12 =~ 1*hyem12
  wsi5 =~ 1*sisoem5
  wsi7 =~ 1*sisoem7
  wsi10 =~ 1*sisoem10
  wsi12 =~ 1*sisoem12
  
  # Constrained lagged effects between the within-person centered variables. 
  wad7 ~ wad5 + wsi5 
  wsi7 ~ wad5 + b*wsi5
  
  wad10 ~ wad7 + wsi7
  wsi10 ~ wad7 + b*wsi7
  
  wad12 ~ wad10 + wsi10
  wsi12 ~ wad10 + b*wsi10
  
  # Estimate the covariance between the within-person centered variables at the first wave
  wad5 ~~ wsi5 # Covariance
  
  # Estimate the covariances between the residuals of the within-person centered variables (the innovations)
  wad7 ~~ wsi7
  wad10 ~~ wsi10
  wad12 ~~ wsi12
  
  # Estimate the variance and covariance of the random intercepts
  RIad ~~ RIad
  RIsi ~~ RIsi
  RIad ~~ RIsi
  
  # Estimate the (residual) variance of the within-person centered variables.
  wad5 ~~ wad5 # Variances
  wsi5 ~~ wsi5 
  wad7 ~~ wad7 # Residual variances
  wsi7 ~~ wsi7 
  wad10 ~~ wad10 
  wsi10 ~~ wsi10 
  wad12 ~~ wad12 
  wsi12 ~~ wsi12
'
```

```{r fit RICLPM_hyp2b}
RICLPM_hyp2b.fit <- lavaan(RICLPM_hyp2b, 
                      data = dat, 
                      missing = 'ML', 
                      meanstructure = TRUE, 
                      int.ov.free = TRUE,
                      se = "robust",
                      estimator = "MLR" #maximum likelihood with robust (Huber-White) standard errors and a scaled (Yuan-Bentler) and robust test statistic
                      ) 
```

```{r LRT for RICLPM_hyp and RICLPM_hyp2b}
lavTestLRT(RICLPM_hyp.fit, RICLPM_hyp2b.fit, method = "satorra.bentler.2010")
```

RICLPM_hyp2b is significantly *worse* fit compared to RICLPM_hyp. Now we will constrain the other set of stability coefficients: cross lag ad->si (c)

#### Constraining cross lag in ADHD to social isolation parameter (RICLPM_hyp2c)

```{r specify RICLPM_hyp2c}
RICLPM_hyp2c <- '
  # Create between components (random intercepts treated as factors here)
  RIad =~ 1*hyem5 + 1*hyem7 + 1*hyem10 + 1*hyem12 #x
  RIsi =~ 1*sisoem5 + 1*sisoem7 + 1*sisoem10 + 1*sisoem12 #y

  # Create within-person centered variables
  wad5 =~ 1*hyem5
  wad7 =~ 1*hyem7
  wad10 =~ 1*hyem10 
  wad12 =~ 1*hyem12
  wsi5 =~ 1*sisoem5
  wsi7 =~ 1*sisoem7
  wsi10 =~ 1*sisoem10
  wsi12 =~ 1*sisoem12
  
  
  # Constrained lagged effects between the within-person centered variables. 
  wad7 ~ wad5 + wsi5 
  wsi7 ~ c*wad5 + wsi5
  
  wad10 ~ wad7 + wsi7
  wsi10 ~ c*wad7 + wsi7
  
  wad12 ~ wad10 + wsi10
  wsi12 ~ c*wad10 + wsi10
  
  # Estimate the covariance between the within-person centered variables at the first wave
  wad5 ~~ wsi5 # Covariance
  
  # Estimate the covariances between the residuals of the within-person centered variables (the innovations)
  wad7 ~~ wsi7
  wad10 ~~ wsi10
  wad12 ~~ wsi12
  
  # Estimate the variance and covariance of the random intercepts
  RIad ~~ RIad
  RIsi ~~ RIsi
  RIad ~~ RIsi
  
  # Estimate the (residual) variance of the within-person centered variables.
  wad5 ~~ wad5 # Variances
  wsi5 ~~ wsi5 
  wad7 ~~ wad7 # Residual variances
  wsi7 ~~ wsi7 
  wad10 ~~ wad10 
  wsi10 ~~ wsi10 
  wad12 ~~ wad12 
  wsi12 ~~ wsi12
'
```

```{r fit RICLPM_hyp2c}
RICLPM_hyp2c.fit <- lavaan(RICLPM_hyp2c, 
                      data = dat, 
                      missing = 'ML', 
                      meanstructure = TRUE, 
                      int.ov.free = TRUE,
                      se = "robust",
                      estimator = "MLR" #maximum likelihood with robust (Huber-White) standard errors and a scaled (Yuan-Bentler) and robust test statistic
                      ) 
```

```{r LRT for RICLPM_hyp and RICLPM_hyp2c}
lavTestLRT(RICLPM_hyp.fit, RICLPM_hyp2c.fit, method = "satorra.bentler.2010")
```

RICLPM_hyp2c is *NOT* significantly worse fit compared to RICLPM_hyp (p = 0.9412). Now we will constrain the other set of stability coefficients: cross lag si->ad (d) *but comparing that model to this one*. 

#### Constraining cross lag in social isolation to ADHD parameter (RICLPM_hyp2d) - **compared to RICLPM_hyp2c** 

Thus, constraining both sets of lags are tested in this model. 

```{r specify RICLPM_hyp2d}
RICLPM_hyp2d <- '
  # Create between components (random intercepts treated as factors here)
  RIad =~ 1*hyem5 + 1*hyem7 + 1*hyem10 + 1*hyem12 #x
  RIsi =~ 1*sisoem5 + 1*sisoem7 + 1*sisoem10 + 1*sisoem12 #y

  # Create within-person centered variables
  wad5 =~ 1*hyem5
  wad7 =~ 1*hyem7
  wad10 =~ 1*hyem10 
  wad12 =~ 1*hyem12
  wsi5 =~ 1*sisoem5
  wsi7 =~ 1*sisoem7
  wsi10 =~ 1*sisoem10
  wsi12 =~ 1*sisoem12
  
  # Constrained lagged effects between the within-person centered variables. 
  wad7 ~ wad5 + d*wsi5 
  wsi7 ~ c*wad5 + wsi5
  
  wad10 ~ wad7 + d*wsi7
  wsi10 ~ c*wad7 + wsi7
  
  wad12 ~ wad10 + d*wsi10
  wsi12 ~ c*wad10 + wsi10
  
  # Estimate the covariance between the within-person centered variables at the first wave
  wad5 ~~ wsi5 # Covariance
  
  # Estimate the covariances between the residuals of the within-person centered variables (the innovations)
  wad7 ~~ wsi7
  wad10 ~~ wsi10
  wad12 ~~ wsi12
  
  # Estimate the variance and covariance of the random intercepts
  RIad ~~ RIad
  RIsi ~~ RIsi
  RIad ~~ RIsi
  
  # Estimate the (residual) variance of the within-person centered variables.
  wad5 ~~ wad5 # Variances
  wsi5 ~~ wsi5 
  wad7 ~~ wad7 # Residual variances
  wsi7 ~~ wsi7 
  wad10 ~~ wad10 
  wsi10 ~~ wsi10 
  wad12 ~~ wad12 
  wsi12 ~~ wsi12
'
```

```{r fit RICLPM_hyp2d}
RICLPM_hyp2d.fit <- lavaan(RICLPM_hyp2d, 
                      data = dat, 
                      missing = 'ML', 
                      meanstructure = TRUE, 
                      int.ov.free = TRUE,
                      se = "robust",
                      estimator = "MLR" #maximum likelihood with robust (Huber-White) standard errors and a scaled (Yuan-Bentler) and robust test statistic
                      ) 

RICLPM_hyp2d.fit.summary <- summary(RICLPM_hyp2d.fit, 
                                    fit.measures = TRUE,
                                    standardized = TRUE)

#Table of model fit 
RICLPM_hyp2d.fit.summary.fit <- table.model.fit(RICLPM_hyp2d.fit.summary)
RICLPM_hyp2d.fit.summary.fit
#Table of regression coefficients and covariances (concurrent associations)
RICLPM_hyp2d.fit.summary.reg <- table.model.coef(model = RICLPM_hyp2d.fit.summary, type = "RICLPM", constraints = "Yes")
RICLPM_hyp2d.fit.summary.reg
```

```{r LRT for RICLPM_hyp2c and RICLPM_hyp2d}
lavTestLRT(RICLPM_hyp2c.fit, RICLPM_hyp2d.fit, method = "satorra.bentler.2010") #comparing to other best fitting model
```

RICLPM_hyp2d is *NOT* significantly worse fit compared to RICLPM_hyp2c.fit (p = 0.3911). 

**The best fitting model (mother report and Hyperactivity/impulsivity ADHD symptoms) is currently RICLPM2d where cross-lags are constrained to be equal across time.**

***

### Constrained grand means (RICLPM_hyp3)

The grand means are the means over all units per occasion. These grand means may be time-varying, or may be fixed to be invariant over time. 

```{r specify RICLPM_hyp3}
RICLPM_hyp3 <- '
  # Create between components (random intercepts treated as factors here)
  RIad =~ 1*hyem5 + 1*hyem7 + 1*hyem10 + 1*hyem12 #x
  RIsi =~ 1*sisoem5 + 1*sisoem7 + 1*sisoem10 + 1*sisoem12 #y

  # Create within-person centered variables
  wad5 =~ 1*hyem5
  wad7 =~ 1*hyem7
  wad10 =~ 1*hyem10 
  wad12 =~ 1*hyem12
  wsi5 =~ 1*sisoem5
  wsi7 =~ 1*sisoem7
  wsi10 =~ 1*sisoem10
  wsi12 =~ 1*sisoem12
  
  # Constrained lagged effects between the within-person centered variables. 
  wad7 ~ wad5 + wsi5 
  wsi7 ~ wad5 + wsi5
  wad10 ~ wad7 + wsi7
  wsi10 ~ wad7 + wsi7
  wad12 ~ wad10 + wsi10
  wsi12 ~ wad10 + wsi10
  
  # Estimate the covariance between the within-person centered variables at the first wave
  wad5 ~~ wsi5 # Covariance
  
  # Estimate the covariances between the residuals of the within-person centered variables (the innovations)
  wad7 ~~ wsi7
  wad10 ~~ wsi10
  wad12 ~~ wsi12
  
  # Estimate the variance and covariance of the random intercepts
  RIad ~~ RIad
  RIsi ~~ RIsi
  RIad ~~ RIsi
  
  # Estimate the (residual) variance of the within-person centered variables
  wad5 ~~ wad5 # Variances
  wsi5 ~~ wsi5 
  wad7 ~~ wad7 # Residual variances
  wsi7 ~~ wsi7 
  wad10 ~~ wad10 
  wsi10 ~~ wsi10 
  wad12 ~~ wad12 
  wsi12 ~~ wsi12
  
  # Constrain the grand means over time
  hyem5 + hyem7 + hyem10 + hyem12 ~ mad*1
  sisoem5 + sisoem7 + sisoem10 + sisoem12 ~ msi*1
'
```

```{r fit RICLPM_hyp3}
RICLPM_hyp3.fit <- lavaan(RICLPM_hyp3, 
                      data = dat, 
                      missing = 'ML', 
                      meanstructure = TRUE, 
                      int.ov.free = TRUE,
                      se = "robust",
                      estimator = "MLR" #maximum likelihood with robust (Huber-White) standard errors and a scaled (Yuan-Bentler) and robust test statistic
                      ) 
```

```{r LRT for RICLPM_hyp and RICLPM_hyp3}
lavTestLRT(RICLPM_hyp.fit, RICLPM_hyp3.fit, method = "satorra.bentler.2010") 
```

If the grand means cannot be constrained to be invariant over time, this implies that on average there is some change in this variable over time, which may reflect some occasion-specific effect, or a developmental trend. By allowing the means to freely vary over time, we account for such average changes over time.

*This makes sense as we have shown variation in social isolation over time and other literature has shown variation in ADHD over time.*

***

# RI-CLPM: Mother report only, Inattention ADHD symptoms

## The basic RI-CLPM model (RICLPM_inat)
The code for specifying the basic RI-CLPM is given below.

```{r specify RICLPM_inat}
RICLPM_inat <- '
  # Create between components (random intercepts treated as factors here)
  RIad =~ 1*inem5 + 1*inem7 + 1*inem10 + 1*inem12 #x
  RIsi =~ 1*sisoem5 + 1*sisoem7 + 1*sisoem10 + 1*sisoem12 #y

  # Create within-person centered variables
  wad5 =~ 1*inem5
  wad7 =~ 1*inem7
  wad10 =~ 1*inem10 
  wad12 =~ 1*inem12
  wsi5 =~ 1*sisoem5
  wsi7 =~ 1*sisoem7
  wsi10 =~ 1*sisoem10
  wsi12 =~ 1*sisoem12
  
  # Estimate the lagged effects between the within-person centered variables
  wad7 + wsi7 ~ wad5 + wsi5
  wad10 + wsi10 ~ wad7 + wsi7
  wad12 + wsi12 ~ wad10 + wsi10
  
  # Estimate the covariance between the within-person centered variables at the first wave
  wad5 ~~ wsi5 # Covariance
  
  # Estimate the covariances between the residuals of the within-person centered variables (the innovations)
  wad7 ~~ wsi7
  wad10 ~~ wsi10
  wad12 ~~ wsi12
  
  # Estimate the variance and covariance of the random intercepts
  RIad ~~ RIad
  RIsi ~~ RIsi
  RIad ~~ RIsi
  
  # Estimate the (residual) variance of the within-person centered variables.
  wad5 ~~ wad5 # Variances
  wsi5 ~~ wsi5 
  wad7 ~~ wad7 # Residual variances
  wsi7 ~~ wsi7 
  wad10 ~~ wad10 
  wsi10 ~~ wsi10 
  wad12 ~~ wad12 
  wsi12 ~~ wsi12
'
```

```{r fit RICLPM_inat}
RICLPM_inat.fit <- lavaan(RICLPM_inat,     # model
                     data = dat,           # data
                     missing = 'ML',       # how to handle missing data 
                     meanstructure = TRUE, # adds intercepts/means to the model for both observed and latent variables
                     se = "robust",        # robust standard errors
                     int.ov.free = TRUE,   # if FALSE, the intercepts of the observed variables are fixed to zero
                     estimator = "MLR" #maximum likelihood with robust (Huber-White) standard errors and a scaled (Yuan-Bentler) and robust test statistic
)

RICLPM_inat.fit.summary <- summary(RICLPM_inat.fit, 
                                    fit.measures = TRUE,
                                    standardized = TRUE)

#Table of model fit 
RICLPM_inat.fit.summary.fit <- table.model.fit(RICLPM_inat.fit.summary)
#Table of regression coefficients and covariances (concurrent associations)
RICLPM_inat.fit.summary.reg <- table.model.coef(model = RICLPM_inat.fit.summary, type = "RICLPM", constraints = "No")
```

## RI-CLPM Constraints over time {.tabset .tabset-fade}

Imposing constraints to the model can be achieved through **pre-multiplication**. It means that we have to prepend the number that we want to fix the parameter to, and an asterisk, to the parameter in the model specification. For example, `F =~ 0*x1` fixes the factor loading of item `x1` to factor `F` to 0. Using pre-multiplication we can also constrain parameters to be the same by giving them the same label. Below we specify an RI-CLPM with the following constraints: 

1) fixed auto-regressive and cross-lagged relations over time, `wx2 ~ a*wx1 + b*wy1; ...`
2) time-invariant (residual) (co-)variances in the within-person part `wx2 ~~ cov*wy2; ...`, `wx2 ~~ vx*wx2; ...`, and `wy2 ~~ vy*wy2; ...`
3) constrained grand means over time, `x1 + ... ~ mx*1` and `y1 + ... ~ my*1`


### Fixed autoregressive and cross-lagged relations over time (RICLPM_inat2)

a = lag in ad
b = lag in si
c = cross lag ad->si
d = cross lag si->ad

#### Constraining all lag and crosslag parameters at once (RICLPM2)

```{r specify RICLPM_inat2}
RICLPM_inat2 <- '
  # Create between components (random intercepts treated as factors here)
  RIad =~ 1*inem5 + 1*inem7 + 1*inem10 + 1*inem12 #x
  RIsi =~ 1*sisoem5 + 1*sisoem7 + 1*sisoem10 + 1*sisoem12 #y

  # Create within-person centered variables
  wad5 =~ 1*inem5
  wad7 =~ 1*inem7
  wad10 =~ 1*inem10 
  wad12 =~ 1*inem12
  wsi5 =~ 1*sisoem5
  wsi7 =~ 1*sisoem7
  wsi10 =~ 1*sisoem10
  wsi12 =~ 1*sisoem12
  
  
  # Constrained lagged effects between the within-person centered variables. 
  wad7 ~ a*wad5 + d*wsi5 
  wsi7 ~ c*wad5 + b*wsi5
  
  wad10 ~ a*wad7 + d*wsi7
  wsi10 ~ c*wad7 + b*wsi7
  
  wad12 ~ a*wad10 + d*wsi10
  wsi12 ~ c*wad10 + b*wsi10
  
  # Estimate the covariance between the within-person centered variables at the first wave
  wad5 ~~ wsi5 # Covariance
  
  # Estimate the covariances between the residuals of the within-person centered variables (the innovations)
  wad7 ~~ wsi7
  wad10 ~~ wsi10
  wad12 ~~ wsi12
  
  # Estimate the variance and covariance of the random intercepts
  RIad ~~ RIad
  RIsi ~~ RIsi
  RIad ~~ RIsi
  
  # Estimate the (residual) variance of the within-person centered variables.
  wad5 ~~ wad5 # Variances
  wsi5 ~~ wsi5 
  wad7 ~~ wad7 # Residual variances
  wsi7 ~~ wsi7 
  wad10 ~~ wad10 
  wsi10 ~~ wsi10 
  wad12 ~~ wad12 
  wsi12 ~~ wsi12
'
```

```{r fit RICLPM_inat2}
RICLPM_inat2.fit <- lavaan(RICLPM_inat2, 
                      data = dat, 
                      missing = 'ML', 
                      meanstructure = TRUE, 
                      int.ov.free = TRUE,
                      se = "robust",
                      estimator = "MLR" #maximum likelihood with robust (Huber-White) standard errors and a scaled (Yuan-Bentler) and robust test statistic
                      ) 
```

```{r LRT for RICLPM_inat and RICLPM_inat2}
lavTestLRT(RICLPM_inat.fit, RICLPM_inat2.fit, method = "satorra.bentler.2010")
```

RICLPM_inat2 is significantly *worse* fit compared to RICLPM_inat. 

#### Constraining lag in ADHD parameter (RICLPM_inat2a)

Now we will test the constraints based on the following steps recommended by Curran and Bauer:
1. Estimate the baseline model where all parameters are freely estimated (The basic RI-CLPM)
2. Impose equlaity constraints on one set of stabilities (within variable lags). Use LRT tests to see if he fit becomes significantly worse. 
3. Impose equality constraints on the next set of stabilities. Compare this to wither model 1 (baseline/basic) if step 2 was significant. Compare to model 2 if step 2 was non-significant. *Non-sig LRT = not significantly worse fit*
4. Repeat for first set of cross-lags
5. Repeat for second set of cross-lags

```{r specify RICLPM_inat2a}
RICLPM_inat2a <- '
  # Create between components (random intercepts treated as factors here)
  RIad =~ 1*inem5 + 1*inem7 + 1*inem10 + 1*inem12 #x
  RIsi =~ 1*sisoem5 + 1*sisoem7 + 1*sisoem10 + 1*sisoem12 #y

  # Create within-person centered variables
  wad5 =~ 1*inem5
  wad7 =~ 1*inem7
  wad10 =~ 1*inem10 
  wad12 =~ 1*inem12
  wsi5 =~ 1*sisoem5
  wsi7 =~ 1*sisoem7
  wsi10 =~ 1*sisoem10
  wsi12 =~ 1*sisoem12
  
  # Constrained lagged effects between the within-person centered variables. 
  wad7 ~ a*wad5 + wsi5 
  wsi7 ~ wad5 + wsi5
  
  wad10 ~ a*wad7 + wsi7
  wsi10 ~ wad7 + wsi7
  
  wad12 ~ a*wad10 + wsi10
  wsi12 ~ wad10 + wsi10
  
  # Estimate the covariance between the within-person centered variables at the first wave
  wad5 ~~ wsi5 # Covariance
  
  # Estimate the covariances between the residuals of the within-person centered variables (the innovations)
  wad7 ~~ wsi7
  wad10 ~~ wsi10
  wad12 ~~ wsi12
  
  # Estimate the variance and covariance of the random intercepts
  RIad ~~ RIad
  RIsi ~~ RIsi
  RIad ~~ RIsi
  
  # Estimate the (residual) variance of the within-person centered variables.
  wad5 ~~ wad5 # Variances
  wsi5 ~~ wsi5 
  wad7 ~~ wad7 # Residual variances
  wsi7 ~~ wsi7 
  wad10 ~~ wad10 
  wsi10 ~~ wsi10 
  wad12 ~~ wad12 
  wsi12 ~~ wsi12
'
```

```{r fit RICLPM_inat2a}
RICLPM_inat2a.fit <- lavaan(RICLPM_inat2a, 
                      data = dat, 
                      missing = 'ML', 
                      meanstructure = TRUE, 
                      int.ov.free = TRUE,
                      se = "robust",
                      estimator = "MLR" #maximum likelihood with robust (Huber-White) standard errors and a scaled (Yuan-Bentler) and robust test statistic
                      ) 
```

```{r LRT for RICLPM_inat and RICLPM_inat2a}
lavTestLRT(RICLPM_inat.fit, RICLPM_inat2a.fit, method = "satorra.bentler.2010")
```

RICLPM_inat2a is significantly *worse* fit compared to RICLPM_inat. Now we will constrain the other set of stability coefficients: lag in si (b)

#### Constraining lag in social isolation parameter (RICLPM_inat2b)

```{r specify RICLPM_inat2b}
RICLPM_inat2b <- '
  # Create between components (random intercepts treated as factors here)
  RIad =~ 1*inem5 + 1*inem7 + 1*inem10 + 1*inem12 #x
  RIsi =~ 1*sisoem5 + 1*sisoem7 + 1*sisoem10 + 1*sisoem12 #y

  # Create within-person centered variables
  wad5 =~ 1*inem5
  wad7 =~ 1*inem7
  wad10 =~ 1*inem10 
  wad12 =~ 1*inem12
  wsi5 =~ 1*sisoem5
  wsi7 =~ 1*sisoem7
  wsi10 =~ 1*sisoem10
  wsi12 =~ 1*sisoem12
  
  # Constrained lagged effects between the within-person centered variables. 
  wad7 ~ wad5 + wsi5 
  wsi7 ~ wad5 + b*wsi5
  
  wad10 ~ wad7 + wsi7
  wsi10 ~ wad7 + b*wsi7
  
  wad12 ~ wad10 + wsi10
  wsi12 ~ wad10 + b*wsi10
  
  # Estimate the covariance between the within-person centered variables at the first wave
  wad5 ~~ wsi5 # Covariance
  
  # Estimate the covariances between the residuals of the within-person centered variables (the innovations)
  wad7 ~~ wsi7
  wad10 ~~ wsi10
  wad12 ~~ wsi12
  
  # Estimate the variance and covariance of the random intercepts
  RIad ~~ RIad
  RIsi ~~ RIsi
  RIad ~~ RIsi
  
  # Estimate the (residual) variance of the within-person centered variables.
  wad5 ~~ wad5 # Variances
  wsi5 ~~ wsi5 
  wad7 ~~ wad7 # Residual variances
  wsi7 ~~ wsi7 
  wad10 ~~ wad10 
  wsi10 ~~ wsi10 
  wad12 ~~ wad12 
  wsi12 ~~ wsi12
'
```

```{r fit RICLPM_inat2b}
RICLPM_inat2b.fit <- lavaan(RICLPM_inat2b, 
                      data = dat, 
                      missing = 'ML', 
                      meanstructure = TRUE, 
                      int.ov.free = TRUE,
                      se = "robust",
                      estimator = "MLR" #maximum likelihood with robust (Huber-White) standard errors and a scaled (Yuan-Bentler) and robust test statistic
                      ) 
```

```{r LRT for RICLPM_inat and RICLPM_inat2b}
lavTestLRT(RICLPM_inat.fit, RICLPM_inat2b.fit, method = "satorra.bentler.2010")
```

RICLPM_inat2b is significantly *worse* fit compared to RICLPM_inat. Now we will constrain the other set of stability coefficients: cross lag ad->si (c)

#### Constraining cross lag in ADHD to social isolation parameter (RICLPM_inat2c)

```{r specify RICLPM_inat2c}
RICLPM_inat2c <- '
  # Create between components (random intercepts treated as factors here)
  RIad =~ 1*inem5 + 1*inem7 + 1*inem10 + 1*inem12 #x
  RIsi =~ 1*sisoem5 + 1*sisoem7 + 1*sisoem10 + 1*sisoem12 #y

  # Create within-person centered variables
  wad5 =~ 1*inem5
  wad7 =~ 1*inem7
  wad10 =~ 1*inem10 
  wad12 =~ 1*inem12
  wsi5 =~ 1*sisoem5
  wsi7 =~ 1*sisoem7
  wsi10 =~ 1*sisoem10
  wsi12 =~ 1*sisoem12
  
  # Constrained lagged effects between the within-person centered variables. 
  wad7 ~ wad5 + wsi5 
  wsi7 ~ c*wad5 + wsi5
  
  wad10 ~ wad7 + wsi7
  wsi10 ~ c*wad7 + wsi7
  
  wad12 ~ wad10 + wsi10
  wsi12 ~ c*wad10 + wsi10
  
  # Estimate the covariance between the within-person centered variables at the first wave
  wad5 ~~ wsi5 # Covariance
  
  # Estimate the covariances between the residuals of the within-person centered variables (the innovations)
  wad7 ~~ wsi7
  wad10 ~~ wsi10
  wad12 ~~ wsi12
  
  # Estimate the variance and covariance of the random intercepts
  RIad ~~ RIad
  RIsi ~~ RIsi
  RIad ~~ RIsi
  
  # Estimate the (residual) variance of the within-person centered variables.
  wad5 ~~ wad5 # Variances
  wsi5 ~~ wsi5 
  wad7 ~~ wad7 # Residual variances
  wsi7 ~~ wsi7 
  wad10 ~~ wad10 
  wsi10 ~~ wsi10 
  wad12 ~~ wad12 
  wsi12 ~~ wsi12
'
```

```{r fit RICLPM_inat2c}
RICLPM_inat2c.fit <- lavaan(RICLPM_inat2c, 
                      data = dat, 
                      missing = 'ML', 
                      meanstructure = TRUE, 
                      int.ov.free = TRUE,
                      se = "robust",
                      estimator = "MLR" #maximum likelihood with robust (Huber-White) standard errors and a scaled (Yuan-Bentler) and robust test statistic
                      ) 
```

```{r LRT for RICLPM_inat and RICLPM_inat2c}
lavTestLRT(RICLPM_inat.fit, RICLPM_inat2c.fit, method = "satorra.bentler.2010")
```

RICLPM_inat2c is *NOT* significantly worse fit compared to RICLPM_inat (p = 0.8485). Now we will constrain the other set of stability coefficients: cross lag si->ad (d) *but comparing that model to this one*. 

#### Constraining cross lag in social isolation to ADHD parameter (RICLPM_inat2d) - **compared to RICLPM_inat2c** 

Thus, constraining both sets of lags are tested in this model. 

```{r specify RICLPM_inat2d}
RICLPM_inat2d <- '
  # Create between components (random intercepts treated as factors here)
  RIad =~ 1*inem5 + 1*inem7 + 1*inem10 + 1*inem12 #x
  RIsi =~ 1*sisoem5 + 1*sisoem7 + 1*sisoem10 + 1*sisoem12 #y

  # Create within-person centered variables
  wad5 =~ 1*inem5
  wad7 =~ 1*inem7
  wad10 =~ 1*inem10 
  wad12 =~ 1*inem12
  wsi5 =~ 1*sisoem5
  wsi7 =~ 1*sisoem7
  wsi10 =~ 1*sisoem10
  wsi12 =~ 1*sisoem12
  
  # Constrained lagged effects between the within-person centered variables. 
  wad7 ~ wad5 + d*wsi5 
  wsi7 ~ c*wad5 + wsi5
  
  wad10 ~ wad7 + d*wsi7
  wsi10 ~ c*wad7 + wsi7
  
  wad12 ~ wad10 + d*wsi10
  wsi12 ~ c*wad10 + wsi10
  
  # Estimate the covariance between the within-person centered variables at the first wave
  wad5 ~~ wsi5 # Covariance
  
  # Estimate the covariances between the residuals of the within-person centered variables (the innovations)
  wad7 ~~ wsi7
  wad10 ~~ wsi10
  wad12 ~~ wsi12
  
  # Estimate the variance and covariance of the random intercepts
  RIad ~~ RIad
  RIsi ~~ RIsi
  RIad ~~ RIsi
  
  # Estimate the (residual) variance of the within-person centered variables.
  wad5 ~~ wad5 # Variances
  wsi5 ~~ wsi5 
  wad7 ~~ wad7 # Residual variances
  wsi7 ~~ wsi7 
  wad10 ~~ wad10 
  wsi10 ~~ wsi10 
  wad12 ~~ wad12 
  wsi12 ~~ wsi12
'
```

```{r fit RICLPM_inat2d}
RICLPM_inat2d.fit <- lavaan(RICLPM_inat2d, 
                      data = dat, 
                      missing = 'ML', 
                      meanstructure = TRUE, 
                      int.ov.free = TRUE,
                      se = "robust",
                      estimator = "MLR" #maximum likelihood with robust (Huber-White) standard errors and a scaled (Yuan-Bentler) and robust test statistic
                      ) 

RICLPM_inat2d.fit.summary <- summary(RICLPM_inat2d.fit, 
                                      fit.measures = TRUE,
                                      standardized = TRUE)

#Table of model fit 
RICLPM_inat2d.fit.summary.fit <- table.model.fit(RICLPM_inat2d.fit.summary)
RICLPM_inat2d.fit.summary.fit
#Table of regression coefficients and covariances (concurrent associations)
RICLPM_inat2d.fit.summary.reg <- table.model.coef(model = RICLPM_inat2d.fit.summary, type = "RICLPM", constraints = "Yes")
RICLPM_inat2d.fit.summary.reg
```

```{r LRT for RICLPM_inat2c and RICLPM_inat2d}
lavTestLRT(RICLPM_inat2c.fit, RICLPM_inat2d.fit, method = "satorra.bentler.2010") #comparing to other best fitting model
```

RICLPM_inat2d is *NOT* significantly worse fit compared to RICLPM_inat2c.fit (p = 0.8625). 

**The best fitting model (mother report and inattention ADHD symptoms) is currently RICLPM2d where cross-lags are constrained to be equal across time.**

***

### Constrained grand means (RICLPM_inat3)

The grand means are the means over all units per occasion. These grand means may be time-varying, or may be fixed to be invariant over time. 

```{r specify RICLPM_inat3}
RICLPM_inat3 <- '
  # Create between components (random intercepts treated as factors here)
  RIad =~ 1*inem5 + 1*inem7 + 1*inem10 + 1*inem12 #x
  RIsi =~ 1*sisoem5 + 1*sisoem7 + 1*sisoem10 + 1*sisoem12 #y

  # Create within-person centered variables
  wad5 =~ 1*inem5
  wad7 =~ 1*inem7
  wad10 =~ 1*inem10 
  wad12 =~ 1*inem12
  wsi5 =~ 1*sisoem5
  wsi7 =~ 1*sisoem7
  wsi10 =~ 1*sisoem10
  wsi12 =~ 1*sisoem12
  
  # Constrained lagged effects between the within-person centered variables. 
  wad7 ~ wad5 + wsi5 
  wsi7 ~ wad5 + wsi5
  wad10 ~ wad7 + wsi7
  wsi10 ~ wad7 + wsi7
  wad12 ~ wad10 + wsi10
  wsi12 ~ wad10 + wsi10
  
  # Estimate the covariance between the within-person centered variables at the first wave
  wad5 ~~ wsi5 # Covariance
  
  # Estimate the covariances between the residuals of the within-person centered variables (the innovations)
  wad7 ~~ wsi7
  wad10 ~~ wsi10
  wad12 ~~ wsi12
  
  # Estimate the variance and covariance of the random intercepts
  RIad ~~ RIad
  RIsi ~~ RIsi
  RIad ~~ RIsi
  
  # Estimate the (residual) variance of the within-person centered variables
  wad5 ~~ wad5 # Variances
  wsi5 ~~ wsi5 
  wad7 ~~ wad7 # Residual variances
  wsi7 ~~ wsi7 
  wad10 ~~ wad10 
  wsi10 ~~ wsi10 
  wad12 ~~ wad12 
  wsi12 ~~ wsi12
  
  # Constrain the grand means over time
  inem5 + inem7 + inem10 + inem12 ~ mad*1
  sisoem5 + sisoem7 + sisoem10 + sisoem12 ~ msi*1
'
```

```{r fit RICLPM_inat3}
RICLPM_inat3.fit <- lavaan(RICLPM_inat3, 
                      data = dat, 
                      missing = 'ML', 
                      meanstructure = TRUE, 
                      int.ov.free = TRUE,
                      se = "robust",
                      estimator = "MLR" #maximum likelihood with robust (Huber-White) standard errors and a scaled (Yuan-Bentler) and robust test statistic
                      ) 
```

```{r LRT for RICLPM_inat and RICLPM_inat3}
lavTestLRT(RICLPM_inat.fit, RICLPM_inat3.fit, method = "satorra.bentler.2010") 
```

If the grand means cannot be constrained to be invariant over time, this implies that on average there is some change in this variable over time, which may reflect some occasion-specific effect, or a developmental trend. By allowing the means to freely vary over time, we account for such average changes over time.

*This makes sense as we have shown variation in social isolation over time and other literature has shown variation in ADHD over time.*

***

# RI-CLPMt: Teacher report only, total ADHD symptoms

## The basic cross-lagged panel model (CLPM) without including random intercepts 

```{r specify CLPMt}
CLPMt <- '
  # Estimate the lagged effects between the observed variables.
  tadhdet7 + sisoet7 ~ tadhdet5 + sisoet5
  tadhdet10 + sisoet10 ~ tadhdet7 + sisoet7
  tadhdet12 + sisoet12 ~ tadhdet10 + sisoet10

  # Estimate the covariance between the observed variables at the first wave. 
  tadhdet5 ~~ sisoet5 # Covariance
  
  # Estimate the covariances between the residuals of the observed variables.
  tadhdet7 ~~ sisoet7
  tadhdet10 ~~ sisoet10
  tadhdet12 ~~ sisoet12
  
  # Estimate the (residual) variance of the observed variables.
  tadhdet5 ~~ tadhdet5 # Variances
  sisoet5 ~~ sisoet5 
  tadhdet7 ~~ tadhdet7 # Residual variances
  sisoet7 ~~ sisoet7 
  tadhdet10 ~~ tadhdet10 
  sisoet10 ~~ sisoet10 
  tadhdet12 ~~ tadhdet12 
  sisoet12 ~~ sisoet12 
'
```

```{r fit CLPMt}
CLPMt.fit <- lavaan(CLPMt, 
                   data = dat, 
                   missing = 'ML',
                   meanstructure = TRUE, 
                   int.ov.free = TRUE,
                   se = "robust",
                   estimator = "MLR" #maximum likelihood with robust (Huber-White) standard errors and a scaled (Yuan-Bentler) and robust test statistic
                   ) 

CLPMt.fit.summary <- summary(CLPMt.fit, 
                              fit.measures = TRUE,
                              standardized = TRUE)

#Table of model fit 
CLPMt.fit.summary.fit <- table.model.fit(CLPMt.fit.summary)
CLPMt.fit.summary.fit
#Table of regression coefficients and covariances (concurrent associations)
CLPMt.fit.summary.reg <- table.model.coef(model = CLPMt.fit.summary, type = "CLPM", constraints = "No")
CLPMt.fit.summary.reg
```

## The basic RI-CLPM model teacher report (RICLPMt)
The code for specifying the basic RI-CLPM teacher report is given below.

```{r specify RICLPMt, class.source = 'fold-show'}
RICLPMt <- '
  # Create between components (random intercepts treated as factors here)
  RIad =~ 1*tadhdet5 + 1*tadhdet7 + 1*tadhdet10 + 1*tadhdet12 #x
  RIsi =~ 1*sisoet5 + 1*sisoet7 + 1*sisoet10 + 1*sisoet12 #y

  # Create within-person centered variables
  wad5 =~ 1*tadhdet5
  wad7 =~ 1*tadhdet7
  wad10 =~ 1*tadhdet10 
  wad12 =~ 1*tadhdet12
  wsi5 =~ 1*sisoet5
  wsi7 =~ 1*sisoet7
  wsi10 =~ 1*sisoet10
  wsi12 =~ 1*sisoet12
  
  # Estimate the lagged effects between the within-person centered variables
  wad7 + wsi7 ~ wad5 + wsi5
  wad10 + wsi10 ~ wad7 + wsi7
  wad12 + wsi12 ~ wad10 + wsi10
  
  # Estimate the covariance between the within-person centered variables at the first wave
  wad5 ~~ wsi5 # Covariance
  
  # Estimate the covariances between the residuals of the within-person centered variables (the innovations)
  wad7 ~~ wsi7
  wad10 ~~ wsi10
  wad12 ~~ wsi12
  
  # Estimate the variance and covariance of the random intercepts
  RIad ~~ RIad
  RIsi ~~ RIsi
  RIad ~~ RIsi
  
  # Estimate the (residual) variance of the within-person centered variables.
  wad5 ~~ wad5 # Variances
  wsi5 ~~ wsi5 
  wad7 ~~ wad7 # Residual variances
  wsi7 ~~ wsi7 
  wad10 ~~ wad10 
  wsi10 ~~ wsi10 
  wad12 ~~ wad12 
  wsi12 ~~ wsi12
'
```

```{r fit RICLPMt}
RICLPMt.fit <- lavaan(RICLPMt,             # model
                     data = dat,           # data
                     missing = 'ML',       # how to handle missing data 
                     meanstructure = TRUE, # adds intercepts/means to the model for both observed and latent variables
                     se = "robust",        # robust standard errors
                     int.ov.free = TRUE,   # if FALSE, the intercepts of the observed variables are fixed to zero
                     estimator = "MLR" #maximum likelihood with robust (Huber-White) standard errors and a scaled (Yuan-Bentler) and robust test statistic
)

RICLPMt.fit.summary <- summary(RICLPMt.fit, 
                              fit.measures = TRUE,
                              standardized = TRUE)

#Table of model fit 
RICLPMt.fit.summary.fit <- table.model.fit(RICLPMt.fit.summary)
#Table of regression coefficients and covariances (concurrent associations)
RICLPMt.fit.summary.reg <- table.model.coef(model = RICLPMt.fit.summary, type = "RICLPM", constraints = "No")
```

## RI-CLPMt Constraints over time {.tabset .tabset-fade}

1) fixed auto-regressive and cross-lagged relations over time, `wx2 ~ a*wx1 + b*wy1; ...`
2) time-invariant (residual) (co-)variances in the within-person part `wx2 ~~ cov*wy2; ...`, `wx2 ~~ vx*wx2; ...`, and `wy2 ~~ vy*wy2; ...`
3) constrained grand means over time, `x1 + ... ~ mx*1` and `y1 + ... ~ my*1`

### Fixed autoregressive and cross-lagged relations over time (RICLPMt2)

a = lag in ad
b = lag in si
c = cross lag ad->si
d = cross lag si->ad

#### Constraining all lag and crosslag parameters at once (RICLPMt2)

```{r specify RICLPMt2}
RICLPMt2 <- '
  # Create between components (random intercepts treated as factors here)
  RIad =~ 1*tadhdet5 + 1*tadhdet7 + 1*tadhdet10 + 1*tadhdet12 #x
  RIsi =~ 1*sisoet5 + 1*sisoet7 + 1*sisoet10 + 1*sisoet12 #y

  # Create within-person centered variables
  wad5 =~ 1*tadhdet5
  wad7 =~ 1*tadhdet7
  wad10 =~ 1*tadhdet10 
  wad12 =~ 1*tadhdet12
  wsi5 =~ 1*sisoet5
  wsi7 =~ 1*sisoet7
  wsi10 =~ 1*sisoet10
  wsi12 =~ 1*sisoet12
  
  # Constrained lagged effects between the within-person centered variables. 
  wad7 ~ a*wad5 + d*wsi5 
  wsi7 ~ c*wad5 + b*wsi5
  
  wad10 ~ a*wad7 + d*wsi7
  wsi10 ~ c*wad7 + b*wsi7
  
  wad12 ~ a*wad10 + d*wsi10
  wsi12 ~ c*wad10 + b*wsi10
  
  # Estimate the covariance between the within-person centered variables at the first wave
  wad5 ~~ wsi5 # Covariance
  
  # Estimate the covariances between the residuals of the within-person centered variables (the innovations)
  wad7 ~~ wsi7
  wad10 ~~ wsi10
  wad12 ~~ wsi12
  
  # Estimate the variance and covariance of the random intercepts
  RIad ~~ RIad
  RIsi ~~ RIsi
  RIad ~~ RIsi
  
  # Estimate the (residual) variance of the within-person centered variables.
  wad5 ~~ wad5 # Variances
  wsi5 ~~ wsi5 
  wad7 ~~ wad7 # Residual variances
  wsi7 ~~ wsi7 
  wad10 ~~ wad10 
  wsi10 ~~ wsi10 
  wad12 ~~ wad12 
  wsi12 ~~ wsi12
'
```

```{r fit RICLPMt2}
RICLPMt2.fit <- lavaan(RICLPMt2, 
                      data = dat, 
                      missing = 'ML', 
                      meanstructure = TRUE, 
                      int.ov.free = TRUE,
                      se = "robust",
                      estimator = "MLR" #maximum likelihood with robust (Huber-White) standard errors and a scaled (Yuan-Bentler) and robust test statistic
                      ) 

RICLPMt2.fit.summary <- summary(RICLPMt2.fit, 
                              fit.measures = TRUE,
                              standardized = TRUE)

#Table of model fit 
RICLPMt2.fit.summary.fit <- table.model.fit(RICLPMt2.fit.summary)
RICLPMt2.fit.summary.fit
#Table of regression coefficients and covariances (concurrent associations)
RICLPMt2.fit.summary.reg <- table.model.coef(model = RICLPMt2.fit.summary, type = "RICLPM", constraints = "Yes")
RICLPMt2.fit.summary.reg
```

```{r LRT for RICLPMt and RICLPMt2}
lavTestLRT(RICLPMt.fit, RICLPMt2.fit, method = "satorra.bentler.2010")
```

RICLPMt2 is *NOT* significantly worse fit compared to RICLPMt.fit (p = 0.7661). 

**The best fitting model (teacher report and total ADHD symptoms) is currently RICLPMt2 where autoregressive AND cross-lags are constrained to be equal across time.**

***

### Constrained grand means teacher report (RICLPMt3)

The grand means are the means over all units per occasion. These grand means may be time-varying, or may be fixed to be invariant over time. 

```{r specify RICLPMt3}
RICLPMt3 <- '
  # Create between components (random intercepts treated as factors here)
  RIad =~ 1*tadhdet5 + 1*tadhdet7 + 1*tadhdet10 + 1*tadhdet12 #x
  RIsi =~ 1*sisoet5 + 1*sisoet7 + 1*sisoet10 + 1*sisoet12 #y

  # Create within-person centered variables
  wad5 =~ 1*tadhdet5
  wad7 =~ 1*tadhdet7
  wad10 =~ 1*tadhdet10 
  wad12 =~ 1*tadhdet12
  wsi5 =~ 1*sisoet5
  wsi7 =~ 1*sisoet7
  wsi10 =~ 1*sisoet10
  wsi12 =~ 1*sisoet12
  
  # Constrained lagged effects between the within-person centered variables. 
  wad7 ~ wad5 + wsi5 
  wsi7 ~ wad5 + wsi5
  wad10 ~ wad7 + wsi7
  wsi10 ~ wad7 + wsi7
  wad12 ~ wad10 + wsi10
  wsi12 ~ wad10 + wsi10
  
  # Estimate the covariance between the within-person centered variables at the first wave
  wad5 ~~ wsi5 # Covariance
  
  # Estimate the covariances between the residuals of the within-person centered variables (the innovations)
  wad7 ~~ wsi7
  wad10 ~~ wsi10
  wad12 ~~ wsi12
  
  # Estimate the variance and covariance of the random intercepts
  RIad ~~ RIad
  RIsi ~~ RIsi
  RIad ~~ RIsi
  
  # Estimate the (residual) variance of the within-person centered variables
  wad5 ~~ wad5 # Variances
  wsi5 ~~ wsi5 
  wad7 ~~ wad7 # Residual variances
  wsi7 ~~ wsi7 
  wad10 ~~ wad10 
  wsi10 ~~ wsi10 
  wad12 ~~ wad12 
  wsi12 ~~ wsi12
  
  # Constrain the grand means over time
  tadhdet5 + tadhdet7 + tadhdet10 + tadhdet12 ~ mad*1
  sisoet5 + sisoet7 + sisoet10 + sisoet12 ~ msi*1
'
```

```{r fit RICLPMt3}
RICLPMt3.fit <- lavaan(RICLPMt3, 
                      data = dat, 
                      missing = 'ML', 
                      meanstructure = TRUE, 
                      int.ov.free = TRUE,
                      se = "robust",
                      estimator = "MLR" #maximum likelihood with robust (Huber-White) standard errors and a scaled (Yuan-Bentler) and robust test statistic
                      ) 
```

```{r LRT for RICLPMt and RICLPMt3}
lavTestLRT(RICLPMt.fit, RICLPMt3.fit, method = "satorra.bentler.2010") 
```

RICLPMt3 is significantly *worse* fit than RICLPMt. 

If the grand means cannot be constrained to be invariant over time, this implies that on average there is some change in this variable over time, which may reflect some occasion-specific effect, or a developmental trend. By allowing the means to freely vary over time, we account for such average changes over time.

*This makes sense as we have shown variation in social isolation over time and other literature has shown variation in ADHD over time.*

***

# RI-CLPM: Teacher report only, Hyperactive/Impulsive ADHD symptoms

## The basic RI-CLPM model (RICLPMt_hyp)
The code for specifying the basic RI-CLPM is given below.

```{r specify RICLPMt_hyp}
RICLPMt_hyp <- '
  # Create between components (random intercepts treated as factors here)
  RIad =~ 1*hyet5 + 1*hyet7 + 1*hyet10 + 1*hyet12 #x
  RIsi =~ 1*sisoet5 + 1*sisoet7 + 1*sisoet10 + 1*sisoet12 #y

  # Create within-person centered variables
  wad5 =~ 1*hyet5
  wad7 =~ 1*hyet7
  wad10 =~ 1*hyet10 
  wad12 =~ 1*hyet12
  wsi5 =~ 1*sisoet5
  wsi7 =~ 1*sisoet7
  wsi10 =~ 1*sisoet10
  wsi12 =~ 1*sisoet12
  
  # Estimate the lagged effects between the within-person centered variables
  wad7 + wsi7 ~ wad5 + wsi5
  wad10 + wsi10 ~ wad7 + wsi7
  wad12 + wsi12 ~ wad10 + wsi10
  
  # Estimate the covariance between the within-person centered variables at the first wave
  wad5 ~~ wsi5 # Covariance
  
  # Estimate the covariances between the residuals of the within-person centered variables (the innovations)
  wad7 ~~ wsi7
  wad10 ~~ wsi10
  wad12 ~~ wsi12
  
  # Estimate the variance and covariance of the random intercepts
  RIad ~~ RIad
  RIsi ~~ RIsi
  RIad ~~ RIsi
  
  # Estimate the (residual) variance of the within-person centered variables.
  wad5 ~~ wad5 # Variances
  wsi5 ~~ wsi5 
  wad7 ~~ wad7 # Residual variances
  wsi7 ~~ wsi7 
  wad10 ~~ wad10 
  wsi10 ~~ wsi10 
  wad12 ~~ wad12 
  wsi12 ~~ wsi12
'
```

```{r fit RICLPMt_hyp}
RICLPMt_hyp.fit <- lavaan(RICLPMt_hyp,     # model
                     data = dat,           # data
                     missing = 'ML',       # how to handle missing data 
                     meanstructure = TRUE, # adds intercepts/means to the model for both observed and latent variables
                     se = "robust",        # robust standard errors
                     int.ov.free = TRUE,   # if FALSE, the intercepts of the observed variables are fixed to zero
                     estimator = "MLR" #maximum likelihood with robust (Huber-White) standard errors and a scaled (Yuan-Bentler) and robust test statistic
)

RICLPMt_hyp.fit.summary <- summary(RICLPMt_hyp.fit, 
                                  fit.measures = TRUE,
                                  standardized = TRUE)

#Table of model fit 
RICLPMt_hyp.fit.summary.fit <- table.model.fit(RICLPMt_hyp.fit.summary)
#Table of regression coefficients and covariances (concurrent associations)
RICLPMt_hyp.fit.summary.reg <- table.model.coef(model = RICLPMt_hyp.fit.summary, type = "RICLPM", constraints = "No")
```

***

## RI-CLPM Constraints over time teacher report {.tabset .tabset-fade}

Imposing constraints to the model can be achieved through **pre-multiplication**. It means that we have to prepend the number that we want to fix the parameter to, and an asterisk, to the parameter in the model specification. For example, `F =~ 0*x1` fixes the factor loading of item `x1` to factor `F` to 0. Using pre-multiplication we can also constrain parameters to be the same by giving them the same label. Below we specify an RI-CLPM with the following constraints: 

1) fixed auto-regressive and cross-lagged relations over time, `wx2 ~ a*wx1 + b*wy1; ...`
2) time-invariant (residual) (co-)variances in the within-person part `wx2 ~~ cov*wy2; ...`, `wx2 ~~ vx*wx2; ...`, and `wy2 ~~ vy*wy2; ...`
3) constrained grand means over time, `x1 + ... ~ mx*1` and `y1 + ... ~ my*1`

### Fixed autoregressive and cross-lagged relations over time (RICLPMt_hyp2)

a = lag in ad
b = lag in si
c = cross lag ad->si
d = cross lag si->ad

#### Constraining all lag and crosslag parameters at once (RICLPMt2)

```{r specify RICLPMt_hyp2}
RICLPMt_hyp2 <- '
  # Create between components (random intercepts treated as factors here)
  RIad =~ 1*hyet5 + 1*hyet7 + 1*hyet10 + 1*hyet12 #x
  RIsi =~ 1*sisoet5 + 1*sisoet7 + 1*sisoet10 + 1*sisoet12 #y

  # Create within-person centered variables
  wad5 =~ 1*hyet5
  wad7 =~ 1*hyet7
  wad10 =~ 1*hyet10 
  wad12 =~ 1*hyet12
  wsi5 =~ 1*sisoet5
  wsi7 =~ 1*sisoet7
  wsi10 =~ 1*sisoet10
  wsi12 =~ 1*sisoet12
  
  # Constrained lagged effects between the within-person centered variables. 
  wad7 ~ a*wad5 + d*wsi5 
  wsi7 ~ c*wad5 + b*wsi5
  
  wad10 ~ a*wad7 + d*wsi7
  wsi10 ~ c*wad7 + b*wsi7
  
  wad12 ~ a*wad10 + d*wsi10
  wsi12 ~ c*wad10 + b*wsi10
  
  # Estimate the covariance between the within-person centered variables at the first wave
  wad5 ~~ wsi5 # Covariance
  
  # Estimate the covariances between the residuals of the within-person centered variables (the innovations)
  wad7 ~~ wsi7
  wad10 ~~ wsi10
  wad12 ~~ wsi12
  
  # Estimate the variance and covariance of the random intercepts
  RIad ~~ RIad
  RIsi ~~ RIsi
  RIad ~~ RIsi
  
  # Estimate the (residual) variance of the within-person centered variables.
  wad5 ~~ wad5 # Variances
  wsi5 ~~ wsi5 
  wad7 ~~ wad7 # Residual variances
  wsi7 ~~ wsi7 
  wad10 ~~ wad10 
  wsi10 ~~ wsi10 
  wad12 ~~ wad12 
  wsi12 ~~ wsi12
'
```

```{r fit RICLPMt_hyp2}
RICLPMt_hyp2.fit <- lavaan(RICLPMt_hyp2, 
                      data = dat, 
                      missing = 'ML', 
                      meanstructure = TRUE, 
                      int.ov.free = TRUE,
                      se = "robust",
                      estimator = "MLR" #maximum likelihood with robust (Huber-White) standard errors and a scaled (Yuan-Bentler) and robust test statistic
                      ) 

RICLPMt_hyp2.fit.summary <- summary(RICLPMt_hyp2.fit, 
                                    fit.measures = TRUE,
                                    standardized = TRUE)

#Table of model fit 
RICLPMt_hyp2.fit.summary.fit <- table.model.fit(RICLPMt_hyp2.fit.summary)
RICLPMt_hyp2.fit.summary.fit
#Table of regression coefficients and covariances (concurrent associations)
RICLPMt_hyp2.fit.summary.reg <- table.model.coef(model = RICLPMt_hyp2.fit.summary, type = "RICLPM", constraints = "Yes")
RICLPMt_hyp2.fit.summary.reg
```

```{r LRT for RICLPMt_hyp and RICLPMt_hyp2}
lavTestLRT(RICLPMt_hyp.fit, RICLPMt_hyp2.fit, method = "satorra.bentler.2010")
```

RICLPMt_hyp2 is *not* significantly worse fit compared to RICLPMt_hyp (p = 0.3735) 

***

### Constrained grand means (RICLPMt_hyp3)

The grand means are the means over all units per occasion. These grand means may be time-varying, or may be fixed to be invariant over time. 

```{r specify RICLPMt_hyp3}
RICLPMt_hyp3 <- '
  # Create between components (random intercepts treated as factors here)
  RIad =~ 1*hyet5 + 1*hyet7 + 1*hyet10 + 1*hyet12 #x
  RIsi =~ 1*sisoet5 + 1*sisoet7 + 1*sisoet10 + 1*sisoet12 #y

  # Create within-person centered variables
  wad5 =~ 1*hyet5
  wad7 =~ 1*hyet7
  wad10 =~ 1*hyet10 
  wad12 =~ 1*hyet12
  wsi5 =~ 1*sisoet5
  wsi7 =~ 1*sisoet7
  wsi10 =~ 1*sisoet10
  wsi12 =~ 1*sisoet12
  
  # Constrained lagged effects between the within-person centered variables. 
  wad7 ~ wad5 + wsi5 
  wsi7 ~ wad5 + wsi5
  wad10 ~ wad7 + wsi7
  wsi10 ~ wad7 + wsi7
  wad12 ~ wad10 + wsi10
  wsi12 ~ wad10 + wsi10
  
  # Estimate the covariance between the within-person centered variables at the first wave
  wad5 ~~ wsi5 # Covariance
  
  # Estimate the covariances between the residuals of the within-person centered variables (the innovations)
  wad7 ~~ wsi7
  wad10 ~~ wsi10
  wad12 ~~ wsi12
  
  # Estimate the variance and covariance of the random intercepts
  RIad ~~ RIad
  RIsi ~~ RIsi
  RIad ~~ RIsi
  
  # Estimate the (residual) variance of the within-person centered variables
  wad5 ~~ wad5 # Variances
  wsi5 ~~ wsi5 
  wad7 ~~ wad7 # Residual variances
  wsi7 ~~ wsi7 
  wad10 ~~ wad10 
  wsi10 ~~ wsi10 
  wad12 ~~ wad12 
  wsi12 ~~ wsi12
  
  # Constrain the grand means over time
  hyet5 + hyet7 + hyet10 + hyet12 ~ mad*1
  sisoet5 + sisoet7 + sisoet10 + sisoet12 ~ msi*1
'
```

```{r fit RICLPMt_hyp3}
RICLPMt_hyp3.fit <- lavaan(RICLPMt_hyp3, 
                      data = dat, 
                      missing = 'ML', 
                      meanstructure = TRUE, 
                      int.ov.free = TRUE,
                      se = "robust",
                      estimator = "MLR" #maximum likelihood with robust (Huber-White) standard errors and a scaled (Yuan-Bentler) and robust test statistic
                      ) 
```

```{r LRT for RICLPMt_hyp and RICLPMt_hyp3}
lavTestLRT(RICLPMt_hyp.fit, RICLPMt_hyp3.fit, method = "satorra.bentler.2010") 
```

RICLPMt_hyp3 is significantly *worse* fit than RICLPMt_hyp. 

If the grand means cannot be constrained to be invariant over time, this implies that on average there is some change in this variable over time, which may reflect some occasion-specific effect, or a developmental trend. By allowing the means to freely vary over time, we account for such average changes over time.

*This makes sense as we have shown variation in social isolation over time and other literature has shown variation in ADHD over time.*

***

# RI-CLPM: Teacher report only, Inattention ADHD symptoms

## The basic RI-CLPM model (RICLPMt_inat)
The code for specifying the basic RI-CLPM is given below.

```{r specify RICLPMt_inat}
RICLPMt_inat <- '
  # Create between components (random intercepts treated as factors here)
  RIad =~ 1*inet5 + 1*inet7 + 1*inet10 + 1*inet12 #x
  RIsi =~ 1*sisoet5 + 1*sisoet7 + 1*sisoet10 + 1*sisoet12 #y

  # Create within-person centered variables
  wad5 =~ 1*inet5
  wad7 =~ 1*inet7
  wad10 =~ 1*inet10 
  wad12 =~ 1*inet12
  wsi5 =~ 1*sisoet5
  wsi7 =~ 1*sisoet7
  wsi10 =~ 1*sisoet10
  wsi12 =~ 1*sisoet12
  
  # Estimate the lagged effects between the within-person centered variables
  wad7 + wsi7 ~ wad5 + wsi5
  wad10 + wsi10 ~ wad7 + wsi7
  wad12 + wsi12 ~ wad10 + wsi10
  
  # Estimate the covariance between the within-person centered variables at the first wave
  wad5 ~~ wsi5 # Covariance
  
  # Estimate the covariances between the residuals of the within-person centered variables (the innovations)
  wad7 ~~ wsi7
  wad10 ~~ wsi10
  wad12 ~~ wsi12
  
  # Estimate the variance and covariance of the random intercepts
  RIad ~~ RIad
  RIsi ~~ RIsi
  RIad ~~ RIsi
  
  # Estimate the (residual) variance of the within-person centered variables.
  wad5 ~~ wad5 # Variances
  wsi5 ~~ wsi5 
  wad7 ~~ wad7 # Residual variances
  wsi7 ~~ wsi7 
  wad10 ~~ wad10 
  wsi10 ~~ wsi10 
  wad12 ~~ wad12 
  wsi12 ~~ wsi12
'
```

```{r fit RICLPMt_inat}
RICLPMt_inat.fit <- lavaan(RICLPMt_inat,   # model
                     data = dat,           # data
                     missing = 'ML',       # how to handle missing data 
                     meanstructure = TRUE, # adds intercepts/means to the model for both observed and latent variables
                     se = "robust",        # robust standard errors
                     int.ov.free = TRUE,   # if FALSE, the intercepts of the observed variables are fixed to zero
                     estimator = "MLR" #maximum likelihood with robust (Huber-White) standard errors and a scaled (Yuan-Bentler) and robust test statistic
)

RICLPMt_inat.fit.summary <- summary(RICLPMt_inat.fit, 
                                    fit.measures = TRUE,
                                    standardized = TRUE)

#Table of model fit 
RICLPMt_inat.fit.summary.fit <- table.model.fit(RICLPMt_inat.fit.summary)
#Table of regression coefficients and covariances (concurrent associations)
RICLPMt_inat.fit.summary.reg <- table.model.coef(model = RICLPMt_inat.fit.summary, type = "RICLPM", constraints = "No")
```

## RI-CLPM Constraints over time {.tabset .tabset-fade}

Imposing constraints to the model can be achieved through **pre-multiplication**. It means that we have to prepend the number that we want to fix the parameter to, and an asterisk, to the parameter in the model specification. For example, `F =~ 0*x1` fixes the factor loading of item `x1` to factor `F` to 0. Using pre-multiplication we can also constrain parameters to be the same by giving them the same label. Below we specify an RI-CLPM with the following constraints: 

1) fixed auto-regressive and cross-lagged relations over time, `wx2 ~ a*wx1 + b*wy1; ...`
2) time-invariant (residual) (co-)variances in the within-person part `wx2 ~~ cov*wy2; ...`, `wx2 ~~ vx*wx2; ...`, and `wy2 ~~ vy*wy2; ...`
3) constrained grand means over time, `x1 + ... ~ mx*1` and `y1 + ... ~ my*1`

### Fixed autoregressive and cross-lagged relations over time (RICLPMt_inat2)

a = lag in ad
b = lag in si
c = cross lag ad->si
d = cross lag si->ad

#### Constraining all lag and crosslag parameters at once (RICLPMt2)

```{r specify RICLPMt_inat2}
RICLPMt_inat2 <- '
  # Create between components (random intercepts treated as factors here)
  RIad =~ 1*inet5 + 1*inet7 + 1*inet10 + 1*inet12 #x
  RIsi =~ 1*sisoet5 + 1*sisoet7 + 1*sisoet10 + 1*sisoet12 #y

  # Create within-person centered variables
  wad5 =~ 1*inet5
  wad7 =~ 1*inet7
  wad10 =~ 1*inet10 
  wad12 =~ 1*inet12
  wsi5 =~ 1*sisoet5
  wsi7 =~ 1*sisoet7
  wsi10 =~ 1*sisoet10
  wsi12 =~ 1*sisoet12
  
  # Constrained lagged effects between the within-person centered variables. 
  wad7 ~ a*wad5 + d*wsi5 
  wsi7 ~ c*wad5 + b*wsi5
  
  wad10 ~ a*wad7 + d*wsi7
  wsi10 ~ c*wad7 + b*wsi7
  
  wad12 ~ a*wad10 + d*wsi10
  wsi12 ~ c*wad10 + b*wsi10
  
  # Estimate the covariance between the within-person centered variables at the first wave
  wad5 ~~ wsi5 # Covariance
  
  # Estimate the covariances between the residuals of the within-person centered variables (the innovations)
  wad7 ~~ wsi7
  wad10 ~~ wsi10
  wad12 ~~ wsi12
  
  # Estimate the variance and covariance of the random intercepts
  RIad ~~ RIad
  RIsi ~~ RIsi
  RIad ~~ RIsi
  
  # Estimate the (residual) variance of the within-person centered variables.
  wad5 ~~ wad5 # Variances
  wsi5 ~~ wsi5 
  wad7 ~~ wad7 # Residual variances
  wsi7 ~~ wsi7 
  wad10 ~~ wad10 
  wsi10 ~~ wsi10 
  wad12 ~~ wad12 
  wsi12 ~~ wsi12
'
```

```{r fit RICLPMt_inat2}
RICLPMt_inat2.fit <- lavaan(RICLPMt_inat2, 
                      data = dat, 
                      missing = 'ML', 
                      meanstructure = TRUE, 
                      int.ov.free = TRUE,
                      se = "robust",
                      estimator = "MLR" #maximum likelihood with robust (Huber-White) standard errors and a scaled (Yuan-Bentler) and robust test statistic
                      ) 

RICLPMt_inat2.fit.summary <- summary(RICLPMt_inat2.fit, 
                                    fit.measures = TRUE,
                                    standardized = TRUE)

#Table of model fit 
RICLPMt_inat2.fit.summary.fit <- table.model.fit(RICLPMt_inat2.fit.summary)
RICLPMt_inat2.fit.summary.fit
#Table of regression coefficients and covariances (concurrent associations)
RICLPMt_inat2.fit.summary.reg <- table.model.coef(model = RICLPMt_inat2.fit.summary, type = "RICLPM", constraints = "Yes")
RICLPMt_inat2.fit.summary.reg
```

```{r LRT for RICLPMt_inat and RICLPMt_inat2}
lavTestLRT(RICLPMt_inat.fit, RICLPMt_inat2.fit, method = "satorra.bentler.2010")
```

RICLPMt_inat2 is *not* significantly worse fit compared to RICLPMt_inat (0.6158)

***

### Constrained grand means (RICLPMt_inat3)

The grand means are the means over all units per occasion. These grand means may be time-varying, or may be fixed to be invariant over time. 

```{r specify RICLPMt_inat3}
RICLPMt_inat3 <- '
  # Create between components (random intercepts treated as factors here)
  RIad =~ 1*inem5 + 1*inem7 + 1*inem10 + 1*inem12 #x
  RIsi =~ 1*sisoem5 + 1*sisoem7 + 1*sisoem10 + 1*sisoem12 #y

  # Create within-person centered variables
  wad5 =~ 1*inem5
  wad7 =~ 1*inem7
  wad10 =~ 1*inem10 
  wad12 =~ 1*inem12
  wsi5 =~ 1*sisoem5
  wsi7 =~ 1*sisoem7
  wsi10 =~ 1*sisoem10
  wsi12 =~ 1*sisoem12
  
  # Constrained lagged effects between the within-person centered variables. 
  wad7 ~ wad5 + wsi5 
  wsi7 ~ wad5 + wsi5
  wad10 ~ wad7 + wsi7
  wsi10 ~ wad7 + wsi7
  wad12 ~ wad10 + wsi10
  wsi12 ~ wad10 + wsi10
  
  # Estimate the covariance between the within-person centered variables at the first wave
  wad5 ~~ wsi5 # Covariance
  
  # Estimate the covariances between the residuals of the within-person centered variables (the innovations)
  wad7 ~~ wsi7
  wad10 ~~ wsi10
  wad12 ~~ wsi12
  
  # Estimate the variance and covariance of the random intercepts
  RIad ~~ RIad
  RIsi ~~ RIsi
  RIad ~~ RIsi
  
  # Estimate the (residual) variance of the within-person centered variables
  wad5 ~~ wad5 # Variances
  wsi5 ~~ wsi5 
  wad7 ~~ wad7 # Residual variances
  wsi7 ~~ wsi7 
  wad10 ~~ wad10 
  wsi10 ~~ wsi10 
  wad12 ~~ wad12 
  wsi12 ~~ wsi12
  
  # Constrain the grand means over time
  inem5 + inem7 + inem10 + inem12 ~ mad*1
  sisoem5 + sisoem7 + sisoem10 + sisoem12 ~ msi*1
'
```

```{r fit RICLPMt_inat3}
RICLPMt_inat3.fit <- lavaan(RICLPMt_inat3, 
                      data = dat, 
                      missing = 'ML', 
                      meanstructure = TRUE, 
                      int.ov.free = TRUE,
                      se = "robust",
                      estimator = "MLR" #maximum likelihood with robust (Huber-White) standard errors and a scaled (Yuan-Bentler) and robust test statistic
                      ) 
```

```{r LRT for RICLPMt_inat and RICLPMt_inat3}
lavTestLRT(RICLPMt_inat.fit, RICLPMt_inat3.fit, method = "satorra.bentler.2010") 
```

If the grand means cannot be constrained to be invariant over time, this implies that on average there is some change in this variable over time, which may reflect some occasion-specific effect, or a developmental trend. By allowing the means to freely vary over time, we account for such average changes over time.

*This makes sense as we have shown variation in social isolation over time and other literature has shown variation in ADHD over time.*

***