0.1 Overview

This script shares The Cortisol Data, an N = 34, T = 9 data set we have used to illustrate a variety of growth modeling and mixture modeling methods. We describe the data and then walk through R-based implementations of the models covered in

Ram, N., & Grimm, K. (2007). Using simple and complex growth models to articulate developmental change: Matching theory to method. International Journal of Behavioral Development, 31, 303-316.

The data are shared with intent that others may find them useful for learning about growth modeling or developing new methods for analysis of change. New publications based on these data require citation and acknowledgement of the full set of papers that have used the data, inlcuding the above paper, and …

Ram, N., & Grimm, K. (2009). Growth mixture modeling: A method for identifying differences in longitudinal change among unobserved groups. International Journal of Behavioral Development, 33, 565-576.
Ram, N., Grimm, K., Gatzke-Kopp, L. & Molenaar, P.C.M. (2011). Longitudinal mixture models and the identification of archetypes: Action-adventure, mystery, science fiction, fantasy, or romance? In B. Laursen, T. Little, & N. Card (Eds.) Handbook of Developmental Research Methods (pp. 481-500). New York: Guilford.
Grimm, K.J., Steele, J.S., Ram, N., & Nesselroade, J.R. (2013). Exploratory latent growth models in the structural equation modeling framework. Structural Equation Modeling, 20, 568-591.

Thank you to John Nesselroade, Teresa Seeman, and Marilyn Albert for data collection and inspiration.

Please contact us with questions about the data and their use.

Nilam & Kevin

0.1.1 Outline

This script covers …

A. Describing The Cortisol Data
B. Fitting the series of growth models from Ram & Grimm (2007): Linear, Quadratic, Latent Basis, Exponential, Multiphase

0.1.1.1 Preliminaries

Loading libraries used in this script.

library(psych)
library(ggplot2)
## Warning: package 'ggplot2' was built under R version 3.4.4
library(nlme)
library(lme4)
## Warning: package 'lme4' was built under R version 3.4.4

0.1.1.2 Reading in the data

Loading the public data

#set filepath for data file
filepath <- "https://quantdev.ssri.psu.edu/sites/qdev/files/TheCortisolData.csv"
#read in the .csv file using the url() function
cortisol_wide <- read.csv(file=url(filepath),header=TRUE)

Looking at the top few rows of the wide data.

head(cortisol_wide,6)
id cort_0 cort_1 cort_2 cort_3 cort_4 cort_5 cort_6 cort_7 cort_8
1 4.2 4.1 9.7 14.0 19.0 18.0 20.0 23.0 24.0
2 5.5 5.6 14.0 16.0 19.0 17.0 18.0 20.0 19.0
3 4.0 3.8 7.5 12.0 14.0 13.0 9.1 8.2 7.9
4 6.1 5.6 14.0 20.0 26.0 23.0 26.0 25.0 26.0
5 4.6 4.4 7.2 12.3 15.8 16.1 17.0 17.8 19.1
6 6.8 9.5 14.2 19.6 19.0 13.9 13.4 12.5 11.7

0.1.1.3 Reshaping the data

Two main data schema are used to accommodate repeated measures data - “Wide Format” and “Long Format”. Different functions work with different kinds of data input. We already have the wide format data. We make a set of long format data.

Reshape from wide to long

#reshaping wide to long
cortisol_long <- reshape(data=cortisol_wide, 
                         timevar=c("time"), 
                         idvar="id",
                         varying=c("cort_0","cort_1","cort_2","cort_3",
                                   "cort_4","cort_5","cort_6","cort_7","cort_8"),
                         direction="long", sep="_")
#sorting for easy viewing
# order by id and time
cortisol_long <- cortisol_long[order(cortisol_long$id,cortisol_long$time), ]

To match the scaling of time used in some of the papers we add an additional time variable that runs from 0 to 1.

cortisol_long$timescaled <- (cortisol_long$time - 0)/8

Looking at the top few rows of the long data.

head(cortisol_long,18)
id time cort timescaled
1.0 1 0 4.2 0.000
1.1 1 1 4.1 0.125
1.2 1 2 9.7 0.250
1.3 1 3 14.0 0.375
1.4 1 4 19.0 0.500
1.5 1 5 18.0 0.625
1.6 1 6 20.0 0.750
1.7 1 7 23.0 0.875
1.8 1 8 24.0 1.000
2.0 2 0 5.5 0.000
2.1 2 1 5.6 0.125
2.2 2 2 14.0 0.250
2.3 2 3 16.0 0.375
2.4 2 4 19.0 0.500
2.5 2 5 17.0 0.625
2.6 2 6 18.0 0.750
2.7 2 7 20.0 0.875
2.8 2 8 19.0 1.000

0.2 A. Describing The Cortisol Data

Basic descriptives of the 9-occasion data.

#means, sd, etc.
describe(cortisol_wide)
vars n mean sd median trimmed mad min max range skew kurtosis se
id 1 34 17.500000 9.958246 17.50 17.500000 12.60210 1 34.0 33.0 0.0000000 -1.3062827 1.7078251
cort_0 2 34 5.247059 2.454462 5.00 4.903571 2.29803 2 12.0 10.0 1.1784853 0.9152094 0.4209368
cort_1 3 34 5.064706 2.593685 4.50 4.782143 2.22390 2 11.3 9.3 1.0159850 -0.0131327 0.4448133
cort_2 4 34 10.847059 3.110124 10.05 10.789286 3.85476 6 16.7 10.7 0.0864489 -1.1724863 0.5333818
cort_3 5 34 17.355882 3.685614 17.35 17.392857 3.48411 7 24.9 17.9 -0.2934763 0.3862287 0.6320776
cort_4 6 34 19.479412 3.722314 19.10 19.492857 2.81694 12 27.1 15.1 -0.0439142 -0.5599728 0.6383716
cort_5 7 34 16.461765 4.180838 16.95 16.414286 4.52193 8 26.0 18.0 0.0779079 -0.5824582 0.7170077
cort_6 8 34 15.029412 4.569366 14.65 14.900000 4.15128 6 26.0 20.0 0.2614218 -0.3075741 0.7836399
cort_7 9 34 15.620588 4.809446 16.15 15.603571 5.55975 7 25.0 18.0 -0.0237796 -0.9863200 0.8248133
cort_8 10 34 15.864706 5.609455 16.00 15.814286 6.07866 6 26.0 20.0 0.1359885 -1.0784551 0.9620137
#boxplot
#boxplot by day
ggplot(data=cortisol_long, aes(x=factor(time), y=cort)) + 
  geom_boxplot(notch = TRUE) +
  stat_summary(fun.y="mean", geom="point", shape=23, size=3, fill="white") +
  labs(x = "Time", y = "Cortisol")

#correlations
# Correlations 
round(cor(cortisol_wide[,c("cort_0","cort_1","cort_2","cort_3",
                           "cort_4","cort_5","cort_6","cort_7","cort_8")], 
          use="complete.obs",method="spearman"),2)
##        cort_0 cort_1 cort_2 cort_3 cort_4 cort_5 cort_6 cort_7 cort_8
## cort_0   1.00   0.97   0.45   0.27   0.23   0.17   0.28   0.22   0.20
## cort_1   0.97   1.00   0.48   0.27   0.18   0.14   0.23   0.19   0.16
## cort_2   0.45   0.48   1.00   0.71   0.66   0.47   0.42   0.39   0.34
## cort_3   0.27   0.27   0.71   1.00   0.71   0.37   0.16   0.10   0.14
## cort_4   0.23   0.18   0.66   0.71   1.00   0.81   0.58   0.50   0.52
## cort_5   0.17   0.14   0.47   0.37   0.81   1.00   0.88   0.84   0.82
## cort_6   0.28   0.23   0.42   0.16   0.58   0.88   1.00   0.97   0.93
## cort_7   0.22   0.19   0.39   0.10   0.50   0.84   0.97   1.00   0.94
## cort_8   0.20   0.16   0.34   0.14   0.52   0.82   0.93   0.94   1.00
#pairs plot from the psych library
pairs.panels(cortisol_wide[,c("cort_0","cort_1","cort_2","cort_3",
                              "cort_4","cort_5","cort_6","cort_7","cort_8"),])