Back to schedule


Growth Curve Modeling

Example Data Set: Kashy

Outcome: satisf_A Satisfaction
Predictor Variables:
time: time in days (there are 14 days with 0 = study midpoint)
genderE: Gender effects coded
man, woman: Gender dummy coded
Moderators: Cavoid_A, Cavoid_P Grand mean centered attachment avoidance (Actor and Partner)

First, read in the new dataset. It’s already in the person-period pairwise sturcture.

library(tidyr)
library(dplyr)
library(nlme)

kashy_ppp <- read.csv(file.choose(), header=TRUE)

head(kashy_ppp)
##   dyadid person Day dayOfWeek time genderE genderS conflict_A conflict_P
## 1      1      1   1         1 -6.5       1       M          1          2
## 2      1      1   2         2 -5.5       1       M          4          3
## 3      1      1   3         3 -4.5       1       M          3          3
## 4      1      1   4         4 -3.5       1       M          2          2
## 5      1      1   5         5 -2.5       1       M          2          2
## 6      1      1   6         6 -1.5       1       M          1          2
##   Cconflict_A Cconflict_P   Cavoid_A Cavoid_A.1 woman man satisf_A
## 1       -0.99        0.01 -0.8877427  -1.137743     0   1      7.0
## 2        2.01        1.01 -0.8877427  -1.137743     0   1      7.0
## 3        1.01        1.01 -0.8877427  -1.137743     0   1      6.5
## 4        0.01        0.01 -0.8877427  -1.137743     0   1      7.0
## 5        0.01        0.01 -0.8877427  -1.137743     0   1      7.0
## 6       -0.99        0.01 -0.8877427  -1.137743     0   1      7.0
##   satisf_P Csatisf_A  Csatis_P
## 1     7.00 0.6756513 0.6756513
## 2     6.75 0.6756513 0.4256513
## 3     6.50 0.1756513 0.1756513
## 4     7.00 0.6756513 0.6756513
## 5     7.00 0.6756513 0.6756513
## 6     7.00 0.6756513 0.6756513

Basic Growth Curve Model for Men Only

First, for illustration purposes, we want to run a growth curve model for men only. We can use the Kashy data but select only men with this syntax:

kashy_men <- kashy_ppp %>%
  filter(genderE == 1)

Instead of gls(), we use the lme() in the nlme package, but now we have occasion nested within person and no need for correlated errors. Now we are using traditional MLM and will make use of the random option instead of the correlation option in the syntax below. The random option asks lme() to estimate variance in the intercepts of satisfaction (time = 0, the study midpoint), and variance in the slopes (change in satisfaction overtime). The subject is dyadid because this variable is serving as our participant ID-recall that we are going to use men only and this is a heterosexual sample.

GC_men <- lme(satisf_A ~ 1 + time,
              data = kashy_men, 
              random = ~ 1 + time|dyadid,
              na.action = na.omit)
  
summary(GC_men)
## Linear mixed-effects model fit by REML
##  Data: kashy_men 
##        AIC      BIC    logLik
##   2913.754 2945.342 -1450.877
## 
## Random effects:
##  Formula: ~1 + time | dyadid
##  Structure: General positive-definite, Log-Cholesky parametrization
##             StdDev     Corr  
## (Intercept) 0.68969184 (Intr)
## time        0.04991018 -0.055
## Residual    0.57494853       
## 
## Fixed effects: satisf_A ~ 1 + time 
##                Value  Std.Error   DF t-value p-value
## (Intercept) 6.259438 0.06964341 1327 89.8784  0.0000
## time        0.019287 0.00621176 1327  3.1049  0.0019
##  Correlation: 
##      (Intr)
## time -0.043
## 
## Standardized Within-Group Residuals:
##        Min         Q1        Med         Q3        Max 
## -6.6635944 -0.2279288  0.0955475  0.3916416  4.6840918 
## 
## Number of Observations: 1431
## Number of Groups: 103

The final fixed equation is:

\[Satisfaction = 6.26 + .019(Time)\]

The Intercept = 6.26, which is interpreted as the average level of satisfaction at time = 0 (the study midpoint). The Coefficient for Time = .019. That is, over time, satisfaction increases .019 units a day. The slope is small, although it is statistically significant-thus, there is some evidence of an average increase in satisfaction over time for men.

As for the random effects, the three important estimates are:

  1. Standard deviation of Intercepts: (Intercept) = .670
  2. Standard deviation of Slopes: time = .050
  3. Correlation between Int/Slope: Corr = -.055

There is positive variance in the intercepts-some men were more satisfied than others at the midpoint. There is also positive variance in the slopes-some men are changing in satisfaction more than others. The slope-intercept covariance is in the negative direction. We would interpret this as, men with higher values at time 0 change more slowly than those with lower values.

Next, we can include the men’s attachment avoidance as a moderator of the level of satisfaction at the study midpoint (main effect of avoidance), and as a moderator of the change in satisfaction over time (time by avoidance interaction). This is the addition of two fixed effects-the random effects remain unchanged. The interpretation of the random effects are as before, but they are interpreted with the effect of avoidance partialled. Cavoid_A is grand mean centered avoidance of the actor.

GC_men_mod <- lme(satisf_A ~ 1 + time + Cavoid_A + time*Cavoid_A,
              data = kashy_men, 
              random = ~ 1 + time|dyadid,
              na.action = na.omit)
  
summary(GC_men_mod)
## Linear mixed-effects model fit by REML
##  Data: kashy_men 
##        AIC      BIC    logLik
##   2908.895 2950.922 -1446.447
## 
## Random effects:
##  Formula: ~1 + time | dyadid
##  Structure: General positive-definite, Log-Cholesky parametrization
##             StdDev     Corr  
## (Intercept) 0.67848447 (Intr)
## time        0.05026826 -0.032
## Residual    0.57764139       
## 
## Fixed effects: satisf_A ~ 1 + time + Cavoid_A + time * Cavoid_A 
##                   Value  Std.Error   DF  t-value p-value
## (Intercept)    6.260975 0.06904110 1313 90.68475  0.0000
## time           0.019266 0.00629264 1313  3.06163  0.0022
## Cavoid_A      -0.139771 0.06878183  100 -2.03210  0.0448
## time:Cavoid_A  0.004503 0.00626001 1313  0.71928  0.4721
##  Correlation: 
##               (Intr) time   Cavd_A
## time          -0.025              
## Cavoid_A      -0.060  0.003       
## time:Cavoid_A  0.003 -0.060 -0.025
## 
## Standardized Within-Group Residuals:
##         Min          Q1         Med          Q3         Max 
## -6.63893351 -0.23027180  0.08512544  0.38955936  4.61792496 
## 
## Number of Observations: 1417
## Number of Groups: 102

The final fixed effects growth equation with attachment avoidance as a moderator is:

\[Satisfaction = 6.26 + .019(Time) -.140(Avoid) + .0045(Time*Avoid)\]

The \(intercept = 6.26\), the predicted satisfaction for men at the study midpoint who are at the mean on avoidance. There is a significant main effect of time, \(b = 0.019\), that is, satisfaction increases by .019 each day for men who are at the mean on avoidance. There is a significant main effect of avoidance, \(b = -.140\), such that men who are higher in avoidance are less satisfied at the study midpoint than men who are lower in avoidance. There is no significant Time X avoidance interaction, \(b = .0045\), but to interpret the direction of this interaction we could calculate the slope for men who are high and low in avoidance:

\(Men's Avoidance SD = 0.9499\)
\(High Avoidance (+1sd) Slope = .019 + .0045* 0.9499 = .0233\)
\(Low Avoidance (-1sd) Slope = .019 + .0045*- 0.9499 = .0147\)

(Alternatively, we could re-center Avoidance to obtain these two values.) Men high on avoidance become more satisfied over the course of the study than men lower on avoidance, but this difference is not significant.

Dyadic Growth Curve Modeling

For dyadic growth curve modeling we are going to start with a two intercept model. Note how man and woman is included below as well as -1. This will give us separate intercepts for women and men. We again use the lme() procedure, but now we need a random = statement as well as a correlation = statement:
The random = statement estimates separate intercepts and slopes for men and women and all of the within-and-between person correlations (see slides for more description of all of these random effects).

The correlation = Statement specifies the variances and covariances between the residuals. weights = varIdent(form = ~1|genderS) asks R to estimate different residual variances for men and women (set to be the same at each time point). There is one correlation between the residuals (i.e., the time-specific correlation between satisfaction scores will be Rho in the output).

kashy_ppp <- kashy_ppp %>%
  mutate(slope_m = man*(time), slope_w = woman*(time), obsid = Day+14*(dyadid-1))
dyadGC_di_two  <- lme(satisf_A ~ man + woman + slope_m + slope_w  - 1, 
                    data = kashy_ppp, 
                    random = ~ man + woman  + slope_m + slope_w - 1|dyadid,
                    correlation = corCompSymm(form = ~1|dyadid/obsid),  
                    weights = varIdent(form = ~1|genderS),
                    na.action = na.omit)

summary(dyadGC_di_two)
## Linear mixed-effects model fit by REML
##  Data: kashy_ppp 
##        AIC      BIC    logLik
##   5679.997 5781.305 -2822.998
## 
## Random effects:
##  Formula: ~man + woman + slope_m + slope_w - 1 | dyadid
##  Structure: General positive-definite, Log-Cholesky parametrization
##          StdDev     Corr                
## man      0.68939425 man    woman  slop_m
## woman    0.55533836  0.811              
## slope_m  0.04997165 -0.058  0.053       
## slope_w  0.04786064 -0.017  0.149  0.462
## Residual 0.57490536                     
## 
## Correlation Structure: Compound symmetry
##  Formula: ~1 | dyadid/obsid 
##  Parameter estimate(s):
##       Rho 
## 0.4325675 
## Variance function:
##  Structure: Different standard deviations per stratum
##  Formula: ~1 | genderS 
##  Parameter estimates:
##        M        F 
## 1.000000 1.106536 
## Fixed effects: satisf_A ~ man + woman + slope_m + slope_w - 1 
##            Value  Std.Error   DF   t-value p-value
## man     6.259147 0.06961211 2760  89.91462   0.000
## woman   6.388408 0.05724112 2760 111.60522   0.000
## slope_m 0.019185 0.00621354 2760   3.08757   0.002
## slope_w 0.009855 0.00630161 2760   1.56392   0.118
##  Correlation: 
##         man    woman  slop_m
## woman    0.784              
## slope_m -0.045  0.040       
## slope_w -0.013  0.106  0.448
## 
## Standardized Within-Group Residuals:
##        Min         Q1        Med         Q3        Max 
## -6.7757804 -0.2417459  0.1127611  0.4210600  4.7334050 
## 
## Number of Observations: 2866
## Number of Groups: 103

Finally, as a reminder, one should include gender in the model as a moderator. This is referred to as the “gender moderation model.” The dyadic growth curve model for the gender moderation model is as follows.

dyadGC_di_int  <- lme(satisf_A ~ genderE + time + genderE*time, 
                    data = kashy_ppp, 
                    random = ~ man + woman  + slope_m + slope_w - 1|dyadid,
                    correlation = corCompSymm(form = ~1|dyadid/obsid),  
                    weights = varIdent(form = ~1|genderS),
                    na.action = na.omit)

summary(dyadGC_di_int)
## Linear mixed-effects model fit by REML
##  Data: kashy_ppp 
##        AIC      BIC    logLik
##   5682.769 5784.077 -2824.385
## 
## Random effects:
##  Formula: ~man + woman + slope_m + slope_w - 1 | dyadid
##  Structure: General positive-definite, Log-Cholesky parametrization
##          StdDev     Corr                
## man      0.68939513 man    woman  slop_m
## woman    0.55533714  0.811              
## slope_m  0.04997161 -0.058  0.053       
## slope_w  0.04786057 -0.017  0.149  0.462
## Residual 0.57490548                     
## 
## Correlation Structure: Compound symmetry
##  Formula: ~1 | dyadid/obsid 
##  Parameter estimate(s):
##       Rho 
## 0.4325681 
## Variance function:
##  Structure: Different standard deviations per stratum
##  Formula: ~1 | genderS 
##  Parameter estimates:
##        M        F 
## 1.000000 1.106536 
## Fixed effects: satisf_A ~ genderE + time + genderE * time 
##                  Value  Std.Error   DF   t-value p-value
## (Intercept)   6.323777 0.05993571 2760 105.50935  0.0000
## genderE      -0.064631 0.02165425 2760  -2.98466  0.0029
## time          0.014520 0.00532387 2760   2.72734  0.0064
## genderE:time  0.004665 0.00328870 2760   1.41841  0.1562
##  Correlation: 
##              (Intr) gendrE time  
## genderE       0.302              
## time          0.022 -0.168       
## genderE:time -0.048  0.036 -0.016
## 
## Standardized Within-Group Residuals:
##        Min         Q1        Med         Q3        Max 
## -6.7757776 -0.2417450  0.1127613  0.4210644  4.7334063 
## 
## Number of Observations: 2866
## Number of Groups: 103

Back to schedule