Spline Modeling

2018/05/05

A few spline models (also known as piecewise models). As in previous posts, ‘affect’ is the name given to values of \(y\) throughout.

1) Growth and Even More Growth

A model that captures a process that increases initially and then increases at an even greater rate once it reaches time point 5. The data generating process:

\[\begin{equation} y_{it} = \begin{cases} 4 + 0.3_{t} + error_{t}, & \text{if time < 5}\\ 8 + 0.9_{t} + error_{t}, & \text{otherwise} \end{cases} \end{equation}\]

The data generating code and plot

library(tidyverse)
library(lavaan)
library(ggplot2)
library(MASS)

N <- 400
time <- 10

intercept_1_mu <- 4
intercept_2_mu <- 8

growth1_mu <- 0.3
growth2_mu <- 0.9

sigma <- matrix(c(1.0, 0.4, 0.5, 0.6, 
                  0.4, 1.0, 0.6, 0.3,
                  0.5, 0.6, 1.0, 0.9,
                  0.6, 0.3, 0.9, 1.0), 4, 4, byrow = T)

df_matrix <- matrix(, ncol = 3, nrow = N*time)


count <- 0

for(i in 1:N){
  
  unob_het_y <- rnorm(1,0,1)
  
  parameter_space <- mvrnorm(1, c(intercept_1_mu, intercept_2_mu, growth1_mu, growth2_mu), sigma)
  
  intercept_1 <- parameter_space[1]
  intercept_2 <- parameter_space[2]
  growth1 <- parameter_space[3]
  growth2 <- parameter_space[4]
  
  for(j in 1:time){
    
    count <- count + 1
    
    if(j < 5){
    df_matrix[count, 1] <- i
    df_matrix[count, 2] <- j
    df_matrix[count, 3] <- intercept_1 + growth1*j + unob_het_y + rnorm(1,0,1)
    
    }else{
      
      df_matrix[count, 1] <- i
      df_matrix[count, 2] <- j
      df_matrix[count, 3] <- intercept_2 + growth2*j + unob_het_y + rnorm(1,0,1)
      
      
    }
  }
  
}

df <- data.frame(df_matrix)

names(df) <- c('id', 'time', 'affect')

df1 <- df %>%
  filter(time < 5)

df2 <- df %>%
  filter(time >= 5)

df_sum1 <- df1 %>%
  group_by(time) %>%
  summarise(
    affect = mean(affect)
  )

df_sum2 <- df2 %>%
  group_by(time) %>%
  summarise(
    affect = mean(affect)
  )

ggplot() + 
  geom_point(data = df1, aes(x = time, y = affect, group = id), color = 'gray85') + 
  geom_line(data = df1, aes(x = time, y = affect, group = id), color = 'gray85') + 
  geom_point(data = df2, aes(x = time, y = affect, group = id), color = 'gray85') + 
  geom_line(data = df2, aes(x = time, y = affect, group = id), color = 'gray85') + 
  geom_line(data = df_sum1, aes(x = time, y = affect)) + 
  geom_line(data = df_sum2, aes(x = time, y = affect))

Estimating the parameters using SEM:

df_wide <- reshape(df, idvar = 'id', timevar = 'time', direction = 'wide')


spline_string <- '

# latent intercept for first half

level1_affect =~ 1*affect.1 + 1*affect.2 + 1*affect.3 + 1*affect.4 + 0*affect.5 + 0*affect.6 + 0*affect.7 + 0*affect.8 + 0*affect.9 + 0*affect.10

# latent intercept for second half

level2_affect =~ 0*affect.1 + 0*affect.2 + 0*affect.3 + 0*affect.4 + 1*affect.5 + 1*affect.6 + 1*affect.7 + 1*affect.8 + 1*affect.9 + 1*affect.10

# latent slope for first half basis coefficients

slope1_affect =~ 1*affect.1 + 2*affect.2 + 3*affect.3 + 4*affect.4 + 0*affect.5 + 0*affect.6 + 0*affect.7 + 0*affect.8 + 0*affect.9 + 0*affect.10

# latent slope for second half basis coefficients

slope2_affect =~ 0*affect.1 + 0*affect.2 + 0*affect.3 + 0*affect.4 + 5*affect.5 + 6*affect.6 + 7*affect.7 + 8*affect.8 + 9*affect.9 + 10*affect.10

# means and variance of latent factors

level1_affect ~~ level1_affect
level2_affect ~~ level2_affect
slope1_affect ~~ slope1_affect
slope2_affect ~~ slope2_affect

# covariance between latent factors

level1_affect ~~ level2_affect
level1_affect ~~ slope1_affect
level1_affect ~~ slope2_affect

level2_affect ~~ slope1_affect
level2_affect ~~ slope2_affect

slope1_affect ~~ slope2_affect

# constrain means of indicators to zero across time

affect.1 ~ 0
affect.2 ~ 0
affect.3 ~ 0
affect.4 ~ 0
affect.5 ~ 0
affect.6 ~ 0
affect.7 ~ 0
affect.8 ~ 0
affect.9 ~ 0
affect.10 ~ 0

# constrain residual variance to equality across time

affect.1 ~~ res_var*affect.1
affect.2 ~~ res_var*affect.2
affect.3 ~~ res_var*affect.3
affect.4 ~~ res_var*affect.4
affect.5 ~~ res_var*affect.5
affect.6 ~~ res_var*affect.6
affect.7 ~~ res_var*affect.7
affect.8 ~~ res_var*affect.8
affect.9 ~~ res_var*affect.9
affect.10 ~~ res_var*affect.10

'

spline_model <- growth(spline_string, data = df_wide)
summary(spline_model, fit.measures = T)
## lavaan (0.5-23.1097) converged normally after 107 iterations
## 
##   Number of observations                           400
## 
##   Estimator                                         ML
##   Minimum Function Test Statistic               67.359
##   Degrees of freedom                                50
##   P-value (Chi-square)                           0.051
## 
## Model test baseline model:
## 
##   Minimum Function Test Statistic            10010.730
##   Degrees of freedom                                45
##   P-value                                        0.000
## 
## User model versus baseline model:
## 
##   Comparative Fit Index (CFI)                    0.998
##   Tucker-Lewis Index (TLI)                       0.998
## 
## Loglikelihood and Information Criteria:
## 
##   Loglikelihood user model (H0)              -7343.079
##   Loglikelihood unrestricted model (H1)      -7309.399
## 
##   Number of free parameters                         15
##   Akaike (AIC)                               14716.157
##   Bayesian (BIC)                             14776.029
##   Sample-size adjusted Bayesian (BIC)        14728.433
## 
## Root Mean Square Error of Approximation:
## 
##   RMSEA                                          0.029
##   90 Percent Confidence Interval          0.000  0.046
##   P-value RMSEA <= 0.05                          0.980
## 
## Standardized Root Mean Square Residual:
## 
##   SRMR                                           0.016
## 
## Parameter Estimates:
## 
##   Information                                 Expected
##   Standard Errors                             Standard
## 
## Latent Variables:
##                    Estimate  Std.Err  z-value  P(>|z|)
##   level1_affect =~                                    
##     affect.1          1.000                           
##     affect.2          1.000                           
##     affect.3          1.000                           
##     affect.4          1.000                           
##     affect.5          0.000                           
##     affect.6          0.000                           
##     affect.7          0.000                           
##     affect.8          0.000                           
##     affect.9          0.000                           
##     affect.10         0.000                           
##   level2_affect =~                                    
##     affect.1          0.000                           
##     affect.2          0.000                           
##     affect.3          0.000                           
##     affect.4          0.000                           
##     affect.5          1.000                           
##     affect.6          1.000                           
##     affect.7          1.000                           
##     affect.8          1.000                           
##     affect.9          1.000                           
##     affect.10         1.000                           
##   slope1_affect =~                                    
##     affect.1          1.000                           
##     affect.2          2.000                           
##     affect.3          3.000                           
##     affect.4          4.000                           
##     affect.5          0.000                           
##     affect.6          0.000                           
##     affect.7          0.000                           
##     affect.8          0.000                           
##     affect.9          0.000                           
##     affect.10         0.000                           
##   slope2_affect =~                                    
##     affect.1          0.000                           
##     affect.2          0.000                           
##     affect.3          0.000                           
##     affect.4          0.000                           
##     affect.5          5.000                           
##     affect.6          6.000                           
##     affect.7          7.000                           
##     affect.8          8.000                           
##     affect.9          9.000                           
##     affect.10        10.000                           
## 
## Covariances:
##                    Estimate  Std.Err  z-value  P(>|z|)
##   level1_affect ~~                                    
##     level2_affect     1.162    0.210    5.521    0.000
##     slope1_affect     0.456    0.094    4.838    0.000
##     slope2_affect     0.547    0.096    5.704    0.000
##   level2_affect ~~                                    
##     slope1_affect     0.206    0.115    1.790    0.073
##     slope2_affect     0.014    0.115    0.123    0.902
##   slope1_affect ~~                                    
##     slope2_affect     0.835    0.067   12.518    0.000
## 
## Intercepts:
##                    Estimate  Std.Err  z-value  P(>|z|)
##    .affect.1          0.000                           
##    .affect.2          0.000                           
##    .affect.3          0.000                           
##    .affect.4          0.000                           
##    .affect.5          0.000                           
##    .affect.6          0.000                           
##    .affect.7          0.000                           
##    .affect.8          0.000                           
##    .affect.9          0.000                           
##    .affect.10         0.000                           
##     level1_affect     3.935    0.091   43.348    0.000
##     level2_affect     8.002    0.111   71.813    0.000
##     slope1_affect     0.296    0.051    5.757    0.000
##     slope2_affect     0.891    0.051   17.581    0.000
## 
## Variances:
##                    Estimate  Std.Err  z-value  P(>|z|)
##     lvl1_ff           1.795    0.237    7.571    0.000
##     lvl2_ff           1.583    0.365    4.342    0.000
##     slp1_ff           0.854    0.075   11.424    0.000
##     slp2_ff           0.970    0.073   13.351    0.000
##    .affct.1 (rs_v)    1.001    0.029   34.641    0.000
##    .affct.2 (rs_v)    1.001    0.029   34.641    0.000
##    .affct.3 (rs_v)    1.001    0.029   34.641    0.000
##    .affct.4 (rs_v)    1.001    0.029   34.641    0.000
##    .affct.5 (rs_v)    1.001    0.029   34.641    0.000
##    .affct.6 (rs_v)    1.001    0.029   34.641    0.000
##    .affct.7 (rs_v)    1.001    0.029   34.641    0.000
##    .affct.8 (rs_v)    1.001    0.029   34.641    0.000
##    .affct.9 (rs_v)    1.001    0.029   34.641    0.000
##    .affc.10 (rs_v)    1.001    0.029   34.641    0.000

The structure of the basis coefficients is the important piece that allows us to capture the change in slope:

'
# latent slope for first half basis coefficients

slope1_affect =~ 1*affect.1 + 2*affect.2 + 3*affect.3 + 4*affect.4 + 0*affect.5 + 0*affect.6 + 0*affect.7 + 0*affect.8 + 0*affect.9 + 0*affect.10

# latent slope for second half basis coefficients

slope2_affect =~ 0*affect.1 + 0*affect.2 + 0*affect.3 + 0*affect.4 + 5*affect.5 + 6*affect.6 + 7*affect.7 + 8*affect.8 + 9*affect.9 + 10*affect.10

'

2) Growth and Negative Growth

A model that captures a process that goes up and then goes down. The data generating process:

\[\begin{equation} y_{it} = \begin{cases} 4 + 0.5_{t} + error_{t}, & \text{if time < 5}\\ 4 - 0.5_{t} + error_{t}, & \text{otherwise} \end{cases} \end{equation}\]

The data generating code and plot

library(tidyverse)
library(lavaan)
library(ggplot2)
library(MASS)

N <- 400
time <- 10

intercept_1_mu <- 4
intercept_2_mu <- 4

growth1_mu <- 0.8
growth2_mu <- -0.8

sigma <- matrix(c(1.0, 0.4, 0.5, 0.6, 
                  0.4, 1.0, 0.6, 0.3,
                  0.5, 0.6, 1.0, 0.9,
                  0.6, 0.3, 0.9, 1.0), 4, 4, byrow = T)

df_matrix_b <- matrix(, ncol = 3, nrow = N*time)


count <- 0

for(i in 1:N){
  
  unob_het_y <- rnorm(1,0,1)
  
  parameter_space <- mvrnorm(1, c(intercept_1_mu, intercept_2_mu, growth1_mu, growth2_mu), sigma)
  
  intercept_1 <- parameter_space[1]
  intercept_2 <- parameter_space[2]
  growth1 <- parameter_space[3]
  growth2 <- parameter_space[4]
  
  for(j in 1:time){
    
    count <- count + 1
    
    if(j < 5){
      df_matrix_b[count, 1] <- i
      df_matrix_b[count, 2] <- j
      df_matrix_b[count, 3] <- intercept_1 + growth1*j + unob_het_y + rnorm(1,0,1)
      
    }else{
      
      df_matrix_b[count, 1] <- i
      df_matrix_b[count, 2] <- j
      df_matrix_b[count, 3] <- intercept_2 + growth2*j + unob_het_y + rnorm(1,0,1)
      
      
    }
  }
  
}

df_b <- data.frame(df_matrix_b)

names(df_b) <- c('id', 'time', 'affect')

df1_b <- df_b %>%
  filter(time < 5)

df2_b <- df_b %>%
  filter(time >= 5)

df_sum1_b <- df1_b %>%
  group_by(time) %>%
  summarise(
    affect = mean(affect)
  )

df_sum2_b <- df2_b %>%
  group_by(time) %>%
  summarise(
    affect = mean(affect)
  )

ggplot() + 
  geom_point(data = df1_b, aes(x = time, y = affect, group = id), color = 'gray85') + 
  geom_line(data = df1_b, aes(x = time, y = affect, group = id), color = 'gray85') + 
  geom_point(data = df2_b, aes(x = time, y = affect, group = id), color = 'gray85') + 
  geom_line(data = df2_b, aes(x = time, y = affect, group = id), color = 'gray85') + 
  geom_line(data = df_sum1_b, aes(x = time, y = affect)) + 
  geom_line(data = df_sum2_b, aes(x = time, y = affect))

Estimating the parameters using SEM:

df_wide_b <- reshape(df_b, idvar = 'id', timevar = 'time', direction = 'wide')


spline_string_b <- '

# latent intercept for first half

level1_affect =~ 1*affect.1 + 1*affect.2 + 1*affect.3 + 1*affect.4 + 0*affect.5 + 0*affect.6 + 0*affect.7 + 0*affect.8 + 0*affect.9 + 0*affect.10

# latent intercept for second half

level2_affect =~ 0*affect.1 + 0*affect.2 + 0*affect.3 + 0*affect.4 + 1*affect.5 + 1*affect.6 + 1*affect.7 + 1*affect.8 + 1*affect.9 + 1*affect.10

# latent slope for first half basis coefficients

slope1_affect =~ 1*affect.1 + 2*affect.2 + 3*affect.3 + 4*affect.4 + 0*affect.5 + 0*affect.6 + 0*affect.7 + 0*affect.8 + 0*affect.9 + 0*affect.10

# latent slope for second half basis coefficients

slope2_affect =~ 0*affect.1 + 0*affect.2 + 0*affect.3 + 0*affect.4 + 5*affect.5 + 6*affect.6 + 7*affect.7 + 8*affect.8 + 9*affect.9 + 10*affect.10

# means and variance of latent factors

level1_affect ~~ level1_affect
level2_affect ~~ level2_affect
slope1_affect ~~ slope1_affect
slope2_affect ~~ slope2_affect

# covariance between latent factors

level1_affect ~~ level2_affect
level1_affect ~~ slope1_affect
level1_affect ~~ slope2_affect

level2_affect ~~ slope1_affect
level2_affect ~~ slope2_affect

slope1_affect ~~ slope2_affect

# constrain means of indicators to zero across time

affect.1 ~ 0
affect.2 ~ 0
affect.3 ~ 0
affect.4 ~ 0
affect.5 ~ 0
affect.6 ~ 0
affect.7 ~ 0
affect.8 ~ 0
affect.9 ~ 0
affect.10 ~ 0

# constrain residual variance to equality across time

affect.1 ~~ res_var*affect.1
affect.2 ~~ res_var*affect.2
affect.3 ~~ res_var*affect.3
affect.4 ~~ res_var*affect.4
affect.5 ~~ res_var*affect.5
affect.6 ~~ res_var*affect.6
affect.7 ~~ res_var*affect.7
affect.8 ~~ res_var*affect.8
affect.9 ~~ res_var*affect.9
affect.10 ~~ res_var*affect.10

'

spline_model_b <- growth(spline_string_b, data = df_wide_b)
summary(spline_model_b, fit.measures = T)
## lavaan (0.5-23.1097) converged normally after 100 iterations
## 
##   Number of observations                           400
## 
##   Estimator                                         ML
##   Minimum Function Test Statistic               47.866
##   Degrees of freedom                                50
##   P-value (Chi-square)                           0.559
## 
## Model test baseline model:
## 
##   Minimum Function Test Statistic            10656.159
##   Degrees of freedom                                45
##   P-value                                        0.000
## 
## User model versus baseline model:
## 
##   Comparative Fit Index (CFI)                    1.000
##   Tucker-Lewis Index (TLI)                       1.000
## 
## Loglikelihood and Information Criteria:
## 
##   Loglikelihood user model (H0)              -7364.913
##   Loglikelihood unrestricted model (H1)      -7340.980
## 
##   Number of free parameters                         15
##   Akaike (AIC)                               14759.827
##   Bayesian (BIC)                             14819.699
##   Sample-size adjusted Bayesian (BIC)        14772.103
## 
## Root Mean Square Error of Approximation:
## 
##   RMSEA                                          0.000
##   90 Percent Confidence Interval          0.000  0.030
##   P-value RMSEA <= 0.05                          1.000
## 
## Standardized Root Mean Square Residual:
## 
##   SRMR                                           0.012
## 
## Parameter Estimates:
## 
##   Information                                 Expected
##   Standard Errors                             Standard
## 
## Latent Variables:
##                    Estimate  Std.Err  z-value  P(>|z|)
##   level1_affect =~                                    
##     affect.1          1.000                           
##     affect.2          1.000                           
##     affect.3          1.000                           
##     affect.4          1.000                           
##     affect.5          0.000                           
##     affect.6          0.000                           
##     affect.7          0.000                           
##     affect.8          0.000                           
##     affect.9          0.000                           
##     affect.10         0.000                           
##   level2_affect =~                                    
##     affect.1          0.000                           
##     affect.2          0.000                           
##     affect.3          0.000                           
##     affect.4          0.000                           
##     affect.5          1.000                           
##     affect.6          1.000                           
##     affect.7          1.000                           
##     affect.8          1.000                           
##     affect.9          1.000                           
##     affect.10         1.000                           
##   slope1_affect =~                                    
##     affect.1          1.000                           
##     affect.2          2.000                           
##     affect.3          3.000                           
##     affect.4          4.000                           
##     affect.5          0.000                           
##     affect.6          0.000                           
##     affect.7          0.000                           
##     affect.8          0.000                           
##     affect.9          0.000                           
##     affect.10         0.000                           
##   slope2_affect =~                                    
##     affect.1          0.000                           
##     affect.2          0.000                           
##     affect.3          0.000                           
##     affect.4          0.000                           
##     affect.5          5.000                           
##     affect.6          6.000                           
##     affect.7          7.000                           
##     affect.8          8.000                           
##     affect.9          9.000                           
##     affect.10        10.000                           
## 
## Covariances:
##                    Estimate  Std.Err  z-value  P(>|z|)
##   level1_affect ~~                                    
##     level2_affect     1.635    0.240    6.820    0.000
##     slope1_affect     0.487    0.106    4.602    0.000
##     slope2_affect     0.547    0.103    5.318    0.000
##   level2_affect ~~                                    
##     slope1_affect     0.596    0.141    4.235    0.000
##     slope2_affect     0.292    0.131    2.233    0.026
##   slope1_affect ~~                                    
##     slope2_affect     1.002    0.079   12.741    0.000
## 
## Intercepts:
##                    Estimate  Std.Err  z-value  P(>|z|)
##    .affect.1          0.000                           
##    .affect.2          0.000                           
##    .affect.3          0.000                           
##    .affect.4          0.000                           
##    .affect.5          0.000                           
##    .affect.6          0.000                           
##    .affect.7          0.000                           
##    .affect.8          0.000                           
##    .affect.9          0.000                           
##    .affect.10         0.000                           
##     level1_affect     3.963    0.093   42.785    0.000
##     level2_affect     3.947    0.122   32.444    0.000
##     slope1_affect     0.749    0.057   13.248    0.000
##     slope2_affect    -0.832    0.054  -15.525    0.000
## 
## Variances:
##                    Estimate  Std.Err  z-value  P(>|z|)
##     lvl1_ff           1.990    0.246    8.081    0.000
##     lvl2_ff           2.669    0.429    6.222    0.000
##     slp1_ff           1.088    0.091   11.995    0.000
##     slp2_ff           1.093    0.081   13.463    0.000
##    .affct.1 (rs_v)    0.961    0.028   34.641    0.000
##    .affct.2 (rs_v)    0.961    0.028   34.641    0.000
##    .affct.3 (rs_v)    0.961    0.028   34.641    0.000
##    .affct.4 (rs_v)    0.961    0.028   34.641    0.000
##    .affct.5 (rs_v)    0.961    0.028   34.641    0.000
##    .affct.6 (rs_v)    0.961    0.028   34.641    0.000
##    .affct.7 (rs_v)    0.961    0.028   34.641    0.000
##    .affct.8 (rs_v)    0.961    0.028   34.641    0.000
##    .affct.9 (rs_v)    0.961    0.028   34.641    0.000
##    .affc.10 (rs_v)    0.961    0.028   34.641    0.000

Notice that the string syntax is the exact same because the process changes at the same point in time, it does not matter if the process changes to ‘more positive’ or ‘more negative.’

3) Negative Growth, Growth, and Negative Growth

Now a process that goes down, goes up, and then goes back down. The data generating process:

\[\begin{equation} y_{it} = \begin{cases} 4 - 0.5_{t} + error_{t}, & \text{if time < 5}\\ 4 + 0.5_{t} + error_{t}, & \text{if 5 < time < 10}\\ 4 - 0.5_{t} + error_{t}, & \text{otherwise} \end{cases} \end{equation}\]

The data generating code and plot

library(tidyverse)
library(lavaan)
library(ggplot2)
library(MASS)

N <- 400
time <- 15

intercept_1_mu <- 4
intercept_2_mu <- 4
intercept_3_mu <- 4

growth1_mu <- -0.5
growth2_mu <- 0.5
growth3_mu <- -0.5

sigma <- matrix(c(1.0, 0.4, -0.5, 0.2, 0.1, 0.6, 
                  0.4, 1.0, 0.1, 0.3, 0.2, 0.5,
                  -0.5, 0.1, 1.0, -0.1, 0.15, -0.1,
                  0.2, 0.3, -0.1, 1.0, 0.2, 0.3,
                  0.1, 0.2, 0.15, 0.2, 1.0, -0.4,
                  0.6, 0.5, -0.1, 0.3, -0.4, 1.0), 6, 6, byrow = T)

df_matrix_c <- matrix(, ncol = 3, nrow = N*time)


count <- 0

for(i in 1:N){
  
  unob_het_y <- rnorm(1,0,1)
  
  parameter_space <- mvrnorm(1, c(intercept_1_mu, intercept_2_mu, intercept_3_mu, growth1_mu, growth2_mu, growth3_mu), sigma)
  
  intercept_1 <- parameter_space[1]
  intercept_2 <- parameter_space[2]
  intercept_3 <- parameter_space[3]
  
  growth1 <- parameter_space[4]
  growth2 <- parameter_space[5]
  growth3 <- parameter_space[6]
  
  for(j in 1:time){
    
    count <- count + 1
    
    if(j < 5){
      df_matrix_c[count, 1] <- i
      df_matrix_c[count, 2] <- j
      df_matrix_c[count, 3] <- intercept_1 + growth1*j + unob_het_y + rnorm(1,0,1)
      
    }else if(j >= 5 && j < 10){
      
      df_matrix_c[count, 1] <- i
      df_matrix_c[count, 2] <- j
      df_matrix_c[count, 3] <- intercept_2 + growth2*j + unob_het_y + rnorm(1,0,1)
      
      
    }else{
      
      df_matrix_c[count, 1] <- i
      df_matrix_c[count, 2] <- j
      df_matrix_c[count, 3] <- intercept_3 + growth3*j + unob_het_y + rnorm(1,0,1)
      
    }
  }
  
}

df_c <- data.frame(df_matrix_c)

names(df_c) <- c('id', 'time', 'affect')

df1_c <- df_c %>%
  filter(time < 5)

df2_c <- df_c %>%
  filter(time >= 5 & time < 10)

df3_c <- df_c %>%
  filter(time >= 10)

df_sum1_c <- df1_c %>%
  group_by(time) %>%
  summarise(
    affect = mean(affect)
  )

df_sum2_c <- df2_c %>%
  group_by(time) %>%
  summarise(
    affect = mean(affect)
  )

df_sum3_c <- df3_c %>%
  group_by(time) %>%
  summarise(
    affect = mean(affect)
  )

ggplot() + 
  geom_point(data = df1_c, aes(x = time, y = affect, group = id), color = 'gray85') + 
  geom_line(data = df1_c, aes(x = time, y = affect, group = id), color = 'gray85') + 
  geom_point(data = df2_c, aes(x = time, y = affect, group = id), color = 'gray85') + 
  geom_line(data = df2_c, aes(x = time, y = affect, group = id), color = 'gray85') + 
  geom_line(data = df_sum1_c, aes(x = time, y = affect)) + 
  geom_line(data = df_sum2_c, aes(x = time, y = affect)) + 
  geom_point(data = df3_c, aes(x = time, y = affect, group = id), color = 'gray85') + 
  geom_line(data = df3_c, aes(x = time, y = affect, group = id), color = 'gray85') + 
  geom_line(data = df_sum3_c, aes(x = time, y = affect))

Now estimate the parameters using SEM:

df_wide_c <- reshape(df_c, idvar = 'id', timevar = 'time', direction = 'wide')


spline_string_c <- '

# latent intercept for first third

level1_affect =~ 1*affect.1 + 1*affect.2 + 1*affect.3 + 1*affect.4 + 0*affect.5 + 0*affect.6 + 0*affect.7 + 0*affect.8 + 0*affect.9 + 0*affect.10 + 0*affect.11 + 0*affect.12 + 0*affect.13 + 0*affect.14 + 0*affect.15

# latent intercept for second third

level2_affect =~ 0*affect.1 + 0*affect.2 + 0*affect.3 + 0*affect.4 + 1*affect.5 + 1*affect.6 + 1*affect.7 + 1*affect.8 + 1*affect.9 + 0*affect.10 + 0*affect.11 + 0*affect.12 + 0*affect.13 + 0*affect.14 + 0*affect.15

# latent intercept for final third

level3_affect =~ 0*affect.1 + 0*affect.2 + 0*affect.3 + 0*affect.4 + 0*affect.5 + 0*affect.6 + 0*affect.7 + 0*affect.8 + 0*affect.9 + 1*affect.10 + 1*affect.11 + 1*affect.12 + 1*affect.13 + 1*affect.14 + 1*affect.15


# latent slope for first third basis coefficients

slope1_affect =~ 1*affect.1 + 2*affect.2 + 3*affect.3 + 4*affect.4 + 0*affect.5 + 0*affect.6 + 0*affect.7 + 0*affect.8 + 0*affect.9 + 0*affect.10 + 0*affect.11 + 0*affect.12 + 0*affect.13 + 0*affect.14 + 0*affect.15

# latent slope for second third basis coefficients

slope2_affect =~ 0*affect.1 + 0*affect.2 + 0*affect.3 + 0*affect.4 + 5*affect.5 + 6*affect.6 + 7*affect.7 + 8*affect.8 + 9*affect.9 + 0*affect.10 + 0*affect.11 + 0*affect.12 + 0*affect.13 + 0*affect.14 + 0*affect.15

# latent slope for final third basis coefficients

slope3_affect =~ 0*affect.1 + 0*affect.2 + 0*affect.3 + 0*affect.4 + 0*affect.5 + 0*affect.6 + 0*affect.7 + 0*affect.8 + 0*affect.9 + 10*affect.10 + 11*affect.11 + 12*affect.12 + 13*affect.13 + 14*affect.14 + 15*affect.15



# means and variance of latent factors

level1_affect ~~ level1_affect
level2_affect ~~ level2_affect
level3_affect ~~ level3_affect
slope1_affect ~~ slope1_affect
slope2_affect ~~ slope2_affect
slope3_affect ~~ slope3_affect

# covariance between latent factors

level1_affect ~~ level2_affect
level1_affect ~~ level3_affect
level1_affect ~~ slope1_affect
level1_affect ~~ slope2_affect
level1_affect ~~ slope3_affect

level2_affect ~~ level3_affect
level2_affect ~~ slope1_affect
level2_affect ~~ slope2_affect
level2_affect ~~ slope3_affect

level3_affect ~~ slope1_affect
level3_affect ~~ slope2_affect
level3_affect ~~ slope3_affect

slope1_affect ~~ slope2_affect
slope1_affect ~~ slope3_affect

slope2_affect ~~ slope3_affect

# constrain means of indicators to zero across time

affect.1 ~ 0
affect.2 ~ 0
affect.3 ~ 0
affect.4 ~ 0
affect.5 ~ 0
affect.6 ~ 0
affect.7 ~ 0
affect.8 ~ 0
affect.9 ~ 0
affect.10 ~ 0

# constrain residual variance to equality across time

affect.1 ~~ res_var*affect.1
affect.2 ~~ res_var*affect.2
affect.3 ~~ res_var*affect.3
affect.4 ~~ res_var*affect.4
affect.5 ~~ res_var*affect.5
affect.6 ~~ res_var*affect.6
affect.7 ~~ res_var*affect.7
affect.8 ~~ res_var*affect.8
affect.9 ~~ res_var*affect.9
affect.10 ~~ res_var*affect.10

'

spline_model_c <- growth(spline_string_c, data = df_wide_c)
summary(spline_model_c, fit.measures = T)
## lavaan (0.5-23.1097) converged normally after 157 iterations
## 
##   Number of observations                           400
## 
##   Estimator                                         ML
##   Minimum Function Test Statistic              115.330
##   Degrees of freedom                               102
##   P-value (Chi-square)                           0.173
## 
## Model test baseline model:
## 
##   Minimum Function Test Statistic            17025.314
##   Degrees of freedom                               105
##   P-value                                        0.000
## 
## User model versus baseline model:
## 
##   Comparative Fit Index (CFI)                    0.999
##   Tucker-Lewis Index (TLI)                       0.999
## 
## Loglikelihood and Information Criteria:
## 
##   Loglikelihood user model (H0)             -11797.925
##   Loglikelihood unrestricted model (H1)     -11740.260
## 
##   Number of free parameters                         33
##   Akaike (AIC)                               23661.850
##   Bayesian (BIC)                             23793.569
##   Sample-size adjusted Bayesian (BIC)        23688.857
## 
## Root Mean Square Error of Approximation:
## 
##   RMSEA                                          0.018
##   90 Percent Confidence Interval          0.000  0.033
##   P-value RMSEA <= 0.05                          1.000
## 
## Standardized Root Mean Square Residual:
## 
##   SRMR                                           0.010
## 
## Parameter Estimates:
## 
##   Information                                 Expected
##   Standard Errors                             Standard
## 
## Latent Variables:
##                    Estimate  Std.Err  z-value  P(>|z|)
##   level1_affect =~                                    
##     affect.1          1.000                           
##     affect.2          1.000                           
##     affect.3          1.000                           
##     affect.4          1.000                           
##     affect.5          0.000                           
##     affect.6          0.000                           
##     affect.7          0.000                           
##     affect.8          0.000                           
##     affect.9          0.000                           
##     affect.10         0.000                           
##     affect.11         0.000                           
##     affect.12         0.000                           
##     affect.13         0.000                           
##     affect.14         0.000                           
##     affect.15         0.000                           
##   level2_affect =~                                    
##     affect.1          0.000                           
##     affect.2          0.000                           
##     affect.3          0.000                           
##     affect.4          0.000                           
##     affect.5          1.000                           
##     affect.6          1.000                           
##     affect.7          1.000                           
##     affect.8          1.000                           
##     affect.9          1.000                           
##     affect.10         0.000                           
##     affect.11         0.000                           
##     affect.12         0.000                           
##     affect.13         0.000                           
##     affect.14         0.000                           
##     affect.15         0.000                           
##   level3_affect =~                                    
##     affect.1          0.000                           
##     affect.2          0.000                           
##     affect.3          0.000                           
##     affect.4          0.000                           
##     affect.5          0.000                           
##     affect.6          0.000                           
##     affect.7          0.000                           
##     affect.8          0.000                           
##     affect.9          0.000                           
##     affect.10         1.000                           
##     affect.11         1.000                           
##     affect.12         1.000                           
##     affect.13         1.000                           
##     affect.14         1.000                           
##     affect.15         1.000                           
##   slope1_affect =~                                    
##     affect.1          1.000                           
##     affect.2          2.000                           
##     affect.3          3.000                           
##     affect.4          4.000                           
##     affect.5          0.000                           
##     affect.6          0.000                           
##     affect.7          0.000                           
##     affect.8          0.000                           
##     affect.9          0.000                           
##     affect.10         0.000                           
##     affect.11         0.000                           
##     affect.12         0.000                           
##     affect.13         0.000                           
##     affect.14         0.000                           
##     affect.15         0.000                           
##   slope2_affect =~                                    
##     affect.1          0.000                           
##     affect.2          0.000                           
##     affect.3          0.000                           
##     affect.4          0.000                           
##     affect.5          5.000                           
##     affect.6          6.000                           
##     affect.7          7.000                           
##     affect.8          8.000                           
##     affect.9          9.000                           
##     affect.10         0.000                           
##     affect.11         0.000                           
##     affect.12         0.000                           
##     affect.13         0.000                           
##     affect.14         0.000                           
##     affect.15         0.000                           
##   slope3_affect =~                                    
##     affect.1          0.000                           
##     affect.2          0.000                           
##     affect.3          0.000                           
##     affect.4          0.000                           
##     affect.5          0.000                           
##     affect.6          0.000                           
##     affect.7          0.000                           
##     affect.8          0.000                           
##     affect.9          0.000                           
##     affect.10        10.000                           
##     affect.11        11.000                           
##     affect.12        12.000                           
##     affect.13        13.000                           
##     affect.14        14.000                           
##     affect.15        15.000                           
## 
## Covariances:
##                    Estimate  Std.Err  z-value  P(>|z|)
##   level1_affect ~~                                    
##     level2_affect     1.620    0.250    6.470    0.000
##     level3_affect     0.282    0.292    0.965    0.335
##     slope1_affect     0.252    0.101    2.491    0.013
##     slope2_affect     0.092    0.091    1.012    0.312
##     slope3_affect     0.570    0.099    5.753    0.000
##   level2_affect ~~                                    
##     level3_affect     0.776    0.416    1.867    0.062
##     slope1_affect     0.277    0.141    1.962    0.050
##     slope2_affect     0.187    0.134    1.402    0.161
##     slope3_affect     0.526    0.137    3.834    0.000
##   level3_affect ~~                                    
##     slope1_affect    -0.202    0.174   -1.163    0.245
##     slope2_affect     0.264    0.160    1.651    0.099
##     slope3_affect    -0.147    0.174   -0.850    0.396
##   slope1_affect ~~                                    
##     slope2_affect     0.216    0.055    3.905    0.000
##     slope3_affect     0.388    0.060    6.511    0.000
##   slope2_affect ~~                                    
##     slope3_affect    -0.329    0.054   -6.052    0.000
## 
## Intercepts:
##                    Estimate  Std.Err  z-value  P(>|z|)
##    .affect.1          0.000                           
##    .affect.2          0.000                           
##    .affect.3          0.000                           
##    .affect.4          0.000                           
##    .affect.5          0.000                           
##    .affect.6          0.000                           
##    .affect.7          0.000                           
##    .affect.8          0.000                           
##    .affect.9          0.000                           
##    .affect.10         0.000                           
##    .affect.11         0.000                           
##    .affect.12         0.000                           
##    .affect.13         0.000                           
##    .affect.14         0.000                           
##    .affect.15         0.000                           
##     level1_affect     3.816    0.091   41.746    0.000
##     level2_affect     3.780    0.130   29.166    0.000
##     level3_affect     3.997    0.160   25.043    0.000
##     slope1_affect    -0.485    0.054   -8.932    0.000
##     slope2_affect     0.532    0.050   10.659    0.000
##     slope3_affect    -0.478    0.052   -9.204    0.000
## 
## Variances:
##                    Estimate  Std.Err  z-value  P(>|z|)
##     lvl1_ff           1.839    0.241    7.645    0.000
##     lvl2_ff           1.608    0.499    3.222    0.001
##     lvl3_ff           1.140    0.795    1.434    0.152
##     slp1_ff           0.977    0.084   11.705    0.000
##     slp2_ff           0.896    0.071   12.708    0.000
##     slp3_ff           1.023    0.076   13.395    0.000
##    .affct.1 (rs_v)    1.002    0.030   33.357    0.000
##    .affct.2 (rs_v)    1.002    0.030   33.357    0.000
##    .affct.3 (rs_v)    1.002    0.030   33.357    0.000
##    .affct.4 (rs_v)    1.002    0.030   33.357    0.000
##    .affct.5 (rs_v)    1.002    0.030   33.357    0.000
##    .affct.6 (rs_v)    1.002    0.030   33.357    0.000
##    .affct.7 (rs_v)    1.002    0.030   33.357    0.000
##    .affct.8 (rs_v)    1.002    0.030   33.357    0.000
##    .affct.9 (rs_v)    1.002    0.030   33.357    0.000
##    .affc.10 (rs_v)    1.002    0.030   33.357    0.000
##    .affc.11           1.041    0.087   11.980    0.000
##    .affc.12           1.090    0.090   12.089    0.000
##    .affc.13           1.070    0.091   11.715    0.000
##    .affc.14           0.982    0.090   10.952    0.000
##    .affc.15           0.953    0.100    9.494    0.000

Again, the basis coefficients are the important piece here:

'


# latent slope for first third basis coefficients

slope1_affect =~ 1*affect.1 + 2*affect.2 + 3*affect.3 + 4*affect.4 + 0*affect.5 + 0*affect.6 + 0*affect.7 + 0*affect.8 + 0*affect.9 + 0*affect.10 + 0*affect.11 + 0*affect.12 + 0*affect.13 + 0*affect.14 + 0*affect.15

# latent slope for second third basis coefficients

slope2_affect =~ 0*affect.1 + 0*affect.2 + 0*affect.3 + 0*affect.4 + 5*affect.5 + 6*affect.6 + 7*affect.7 + 8*affect.8 + 9*affect.9 + 0*affect.10 + 0*affect.11 + 0*affect.12 + 0*affect.13 + 0*affect.14 + 0*affect.15

# latent slope for final third basis coefficients

slope3_affect =~ 0*affect.1 + 0*affect.2 + 0*affect.3 + 0*affect.4 + 0*affect.5 + 0*affect.6 + 0*affect.7 + 0*affect.8 + 0*affect.9 + 10*affect.10 + 11*affect.11 + 12*affect.12 + 13*affect.13 + 14*affect.14 + 15*affect.15



'
## [1] "\n\n\n# latent slope for first third basis coefficients\n\nslope1_affect =~ 1*affect.1 + 2*affect.2 + 3*affect.3 + 4*affect.4 + 0*affect.5 + 0*affect.6 + 0*affect.7 + 0*affect.8 + 0*affect.9 + 0*affect.10 + 0*affect.11 + 0*affect.12 + 0*affect.13 + 0*affect.14 + 0*affect.15\n\n# latent slope for second third basis coefficients\n\nslope2_affect =~ 0*affect.1 + 0*affect.2 + 0*affect.3 + 0*affect.4 + 5*affect.5 + 6*affect.6 + 7*affect.7 + 8*affect.8 + 9*affect.9 + 0*affect.10 + 0*affect.11 + 0*affect.12 + 0*affect.13 + 0*affect.14 + 0*affect.15\n\n# latent slope for final third basis coefficients\n\nslope3_affect =~ 0*affect.1 + 0*affect.2 + 0*affect.3 + 0*affect.4 + 0*affect.5 + 0*affect.6 + 0*affect.7 + 0*affect.8 + 0*affect.9 + 10*affect.10 + 11*affect.11 + 12*affect.12 + 13*affect.13 + 14*affect.14 + 15*affect.15\n\n\n\n"

Bo\(^2\)m