Chapter 26 Week 5 - Class
The effect of professional child care on the cognitive, emotional, and social development of children is a much debated toping among policy makers as well as parents with young children. In particular the effect of child care during the first year of life has been of interest.
The dataset WorkingMom_2014.sav describes a study in which the researchers are interested in the relationship between mothers’ work status during the first year of life and their children’s cognitive development by age 4.5 years, measured with the Woodcock Johnson Achievement and Cognitive Batteries (WCJ). Other variables that are considered of interest and which were included are: Mother’s earnings at age 4.5; Home environment quality; Mother’s sensitivity; and Mother’s depression.
Mother’s work status has three levels:
- full time working moms
- part time working moms, and
- stay-at-home moms.
26.0.1 The analysis
If you were to analyze these data, you could do an ANCOVA analysis, with Work (mother’s working status) as factor, and the other independent variables as covariates. Here you find the results of the ANCOVA analysis (assumptions of the ANCOVA were checked and okay).
## Analysis of Variance Table
##
## Response: WCJ
## Df Sum Sq Mean Sq F value Pr(>F)
## Work 2 766 382.8 1.7807 0.1691
## MothEar 1 369 368.9 1.7160 0.1906
## HomeEnv 1 4982 4982.0 23.1749 1.743e-06 ***
## MothSen 1 5023 5022.6 23.3638 1.584e-06 ***
## MothDep 1 4 4.0 0.0185 0.8918
## Residuals 872 187456 215.0
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1These are the parameter estimates for that model:
##
## Call:
## lm(formula = WCJ ~ Work + MothEar + HomeEnv + MothSen + MothDep,
## data = data)
##
## Residuals:
## Min 1Q Median 3Q Max
## -42.903 -9.835 -0.324 9.759 47.868
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 79.39820 5.52998 14.358 < 2e-16 ***
## WorkPart time -0.84041 1.24223 -0.677 0.4989
## Workno work 2.69758 1.33528 2.020 0.0437 *
## MothEar -0.48422 0.50179 -0.965 0.3348
## HomeEnv 0.58431 0.12485 4.680 3.33e-06 ***
## MothSen 6.02656 1.24888 4.826 1.65e-06 ***
## MothDep -0.01705 0.12525 -0.136 0.8918
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 14.66 on 872 degrees of freedom
## Multiple R-squared: 0.05611, Adjusted R-squared: 0.04961
## F-statistic: 8.639 on 6 and 872 DF, p-value: 3.789e-0926.0.2 Question 1
Based on the results, what would be your advise for women with children under one year of age, and why?
26.0.3 Dummy coding
Consider the other variables that were included. We may expect that these are actually influenced by mother’s work status. For instance, mother’s earnings by age 4.5 may depend on whether or not the mother was working during the first year, and also on whether she was working full time or part time. This implies there may be mediated or indirect effects of work status on the outcome variable.
To investigate this, we can specify a SEM model in which the effects of work status are mediated through the four other variables. Note that work status is a categorical (ordinal) variable, with three categories. To include this variable as a predictor, we need to dummy-code it. lm() does this automatically; lavaan doesn’t. The function dummy.code in library(psych) can help us:
library(psych)
dummy.code(data$Work)## Full time Part time no work
## [1,] 1 0 0
## [2,] 1 0 0
## [3,] 1 0 0
## [4,] 1 0 0
## [5,] 1 0 0
## [6,] 1 0 0We need to add these variables to our data. Also, the names of the dummies include spaces, so we should rename them. One easy way to do this is to change the levels of the original variable. Then, we add the dummies to data:
levels(data$Work)## [1] "Full time" "Part time" "no work"levels(data$Work) <- c("workfull", "workpart", "workno")
data <- data.frame(data, dummy.code(data$Work))If we want to use the stay-at-home moms as the reference group, we just use the other two dummies as predictors. Additionally, if we want to know the intercept of WCJ for the reference category, we need to add the argument meanstructure = TRUE to the sem() function call.
26.0.4 Question 2
We now want to specify a SEM model in which the effects of work status are mediated through the four other variables. What is the syntax of this model?
Click for explanation
model <- '
WCJ ~ MothEar + HomeEnv + MothSen + MothDep + workfull + workpart
MothEar ~ workfull + workpart
HomeEnv ~ workfull + workpart
MothSen ~ workfull + workpart
MothDep ~ workfull + workpart
'26.0.5 Question 3
If you had no access to lavaan, what regression analyses could you run in order to piece together this mediated path model? What else would we need to do?
Click for explanation
You would estimate these regression analyses:
fit_y <- lm(WCJ ~ MothEar + HomeEnv + MothSen + MothDep + workfull + workpart, data)
fit_m1 <- lm(MothEar ~ workfull + workpart, data)
fit_m2 <- lm(HomeEnv ~ workfull + workpart, data)
fit_m3 <- lm(MothSen ~ workfull + workpart, data)
fit_m4 <- lm(MothDep ~ workfull + workpart, data)Additionally, we would have to multiply the effects of the dummies on the mediators, with the effects of the mediators on WCJ, in order to calculate the indirect effects.
26.0.6 Question 4
Add syntax to compute the indirect effects to your model from question 2, and run the analysis. Apply an appropriate solution to test whether the indirect effects are significant or not. Discuss the indirect effects of work status on the outcome variable. Taking the direct effects into account, can you explain the differences between the two dummy variables?
Click for explanation
First, label all regressions, and add syntax for the indirect effects:
model <- '
WCJ ~ ear * MothEar + env * HomeEnv + sen * MothSen + dep * MothDep + workfull + workpart
MothEar ~ earfull * workfull + earpart * workpart
HomeEnv ~ envfull * workfull + envpart * workpart
MothSen ~ senfull * workfull + senpart * workpart
MothDep ~ depfull * workfull + deppart * workpart
i_earfull := ear * earfull
i_envfull := env * envfull
i_senfull := sen * senfull
i_depfull := dep * depfull
i_earpart := ear * earpart
i_envpart := env * envpart
i_senpart := sen * senpart
i_deppart := dep * deppart
'Then, estimate the model. Specify bootstrapped SEs and 1000 bootstrap samples. Remember that, to get the intercepts of dependent variables for the reference category, we should specify meanstructure = TRUE:
library(lavaan)
fit <- sem(model,
data,
meanstructure = TRUE,
se = "bootstrap",
bootstrap = 1000)To obtain the 95% CIs for the indirect effects, we can use the following syntax. First, I obtain the parameterestimates() and put them in an object called pars. Then, I ask for only the rows that contain my defined parameters (pars$op == ":="), and only the columns c("label", "est", "ci.upper", "ci.lower"):
pars <- parameterestimates(fit, boot.ci.type = "bca.simple")
pars[pars$op == ":=", c("label", "est", "ci.upper", "ci.lower")]## label est ci.upper ci.lower
## i_earfull -0.323 0.267 -0.972
## i_envfull 0.290 0.782 -0.086
## i_senfull 0.436 0.971 0.089
## i_depfull -0.014 0.222 -0.283
## i_earpart -0.111 0.064 -0.459
## i_envpart 0.478 1.136 0.064
## i_senpart 0.818 1.487 0.320
## i_deppart 0.001 0.163 -0.10426.0.7 Question 5
If we want to determine whether working or not is harmful for children’s cognitive development, what should we consider: the direct effects, the indirect effects, or the total effects, and why?
26.0.8 Question 6
Compare the effect of the dummy variables working full time and working part time in a regression model with all control variables, to the parameters in the lavaan model. Which parameter estimates from the regression correspond to your lavaan model, and why?
Click for explanation
First, run the requested regression:
fit_reg <- lm(WCJ ~ workfull + workpart + MothEar + HomeEnv + MothSen + MothDep, data)
summary(fit_reg)##
## Call:
## lm(formula = WCJ ~ workfull + workpart + MothEar + HomeEnv +
## MothSen + MothDep, data = data)
##
## Residuals:
## Min 1Q Median 3Q Max
## -42.903 -9.835 -0.324 9.759 47.868
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 82.09578 5.47602 14.992 < 2e-16 ***
## workfull -2.69758 1.33528 -2.020 0.0437 *
## workpart -3.53799 1.51825 -2.330 0.0200 *
## MothEar -0.48422 0.50179 -0.965 0.3348
## HomeEnv 0.58431 0.12485 4.680 3.33e-06 ***
## MothSen 6.02656 1.24888 4.826 1.65e-06 ***
## MothDep -0.01705 0.12525 -0.136 0.8918
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 14.66 on 872 degrees of freedom
## Multiple R-squared: 0.05611, Adjusted R-squared: 0.04961
## F-statistic: 8.639 on 6 and 872 DF, p-value: 3.789e-09These coefficients are all identical to the direct effects on WCJ from your lavaan model.
Differences are that your lavaan model ADDITIONALLY estimates the effects of the dummies on the mediators, AND computes the indirect effects. The lavaan model also uses bootstrapped standard errors. As a result, some of the conclusions about significance may differ, and our confidence intervals for the parameters differ.
26.0.9 Question 7
What are your most important conclusions about the effect of working status on cognitive development? You can use your answers to the previous questions. But a model with many variables such as this one has many possible interesting options to look at. Just look around for any interesting results!