use "INSERT PATH TO IMPUTED DATA FILE HERE/GSSimputed.dta", clear
*mi import flong, m(imp) id(id) imp(sibs-SAHM_childXFamilyIncome)
mi import flong, m(imp) id(id) imp(polviews fund reliten class comprend ///
FamilyIncome income_norinc real_satjob fematt fematt_mean ///
year1Xpolviews year1XFamilyIncome SAHM_childXyear1Xfematt ///
SAHM_childXfematt year1Xfematt SAHM_childXFamilyIncome ///
SAHM_childXyear1XFamilyIncome year1Xincome_norinc)
/*check import*/
mi describe
*mi varying
drop imp
/*rename variables so theirs names are not too long (i.e., less than 33 characters)*/
mi rename SAHM_childXyear1XFamilyIncome SAHMXyearXFIncome
/*reshape to normal MI format*/
mi convert wide, clear
/*add in attitudes towards being a SAHM from WVS*/
/*those who strongly agree that being a housewife is just as fulfililng as working for pay*/
/*numbers from online data analysis tool of US data on WVS website, missing excluded*/
/*years are 1990, 1995, 1999, 2006, 2011*/
/*when a WVS year falls between two GSS years, assign the WVS values to
both GSS years*/
gen housewife_fulfil = .
replace housewife_fulfil = 21 if year1==16
replace housewife_fulfil = 28 if year1==20
replace housewife_fulfil = 26 if year1==24
replace housewife_fulfil = 32 if year1==32
replace housewife_fulfil = 24 if year1==36
/*collapse bottom two categories of real_satjob*/
gen satjob3cat = 1 if real_satjob < 3
replace satjob3cat = 2 if real_satjob==3
replace satjob3cat = 3 if real_satjob==4
label define years 0 "1974" 1 "1975" 2 "1976" ///
3 "1977" 4 "1978" 6 "1980" 8 "1982" 9 "1983" 10 "1984" 11 "1985" ///
12 "1986" 13 "1987" 14 "1988" 15 "1989" 16 "1990" 17 "1991" ///
19 "1993" 20 "1994" 22 "1996" 24 "1998" 26 "2000" 28 "2002" ///
30 "2004" 32 "2006" 34 "2008" 36 "2010" 38 "2012" 40 "2014", replace
label values year1 years
label define moms 0 "employed mother" 1 "SAHM"
label values SAHM_child moms
label define edlevels 0 "high school or less" 1 "associates degree or higher"
label values highed edlevels
/*drop anyone originally missing on the outcome*/
gen dropobs = missing(real_satjob)
drop if dropobs==1
mi update
/*are response patterns similar for SAHMs who have and have not ever worked?*/
tab satjob3cat everwork if SAHM_child==1 [iw=wtssall], chi2
reg satjob3cat everwork if SAHM_child==1 [pw=wtssall] /*use reg not ttest so can use weights*/
/*not using robust standard errors - logit models are biased if there are
omitted variables anyhow (even if the omitted variable is not correlated
with the predictor of interest), so not much reason to put a robust SE
on a coefficient that might be wrong - see Greene's (2008) discussion of this*/
/*using iweights because for some reason you can't use pweights with
the vce(oim) option*/
/* Table 1, model 1 - basic trend, no controls */
mi est, dots post: ologit satjob3cat i.SAHM_child##c.year1 [iw=wtssall], vce(oim)
est store mbase
/* Figure 1 - predicted probabilities */
/*note, the 84% CIs are to give a better approximation of where the differences
between the two lines are actually significantly different from one another
- you can use this R script to calculate it, but it returns a CI of 83.4
#figure out how big the SE bars should be using the Goldstein and Healy 1995 method
#figure out appropriate z-value comparisons between all the means, based on SEs of the predictions
#so enter all the SEs of predictions
sigma=scan()
#this function gets the SE ratio - see equation 2 in Goldstein and Healy 1995
adj.comp=function(x){
return((x[1]+x[2])/sqrt(x[1]^2 + x[2]^2))
}
#get ratios for all two-way combinations
adj=combn(sigma, 2, adj.comp)
#go with the highest ratio so that you are being conservative for all comparisons
ratio=max(adj)
z=seq(from=1, to=2, by=.001)
#find the z value that corresponds to 95% confidence
findZ=function(z.vec){
prob.level=c()
for(i in 1:length(z.vec)){
prob.level=c(prob.level, 2*(1-pnorm(z.vec[i]*ratio)))
}
return(prob.level)
}
pvals=findZ(z)
ind=which(abs(.05-pvals)==min(abs(.05-pvals)))
z[ind]
1-2*pnorm(z[ind], lower.tail=F)
*/
mimrgns SAHM_child, at(year1=(0(1)4 6 8(1)17 19 20(2)38)) predict(outcome(3)) ///
post saving(satwork3, replace)
est store out3
marginsplot, title("SAHM and Employed Mothers:") ///
sub("Predicted Probability of being Very Satisfied with Work") xtitle("Year") ///
ytitle("Predicted Probability") level(84) ///
xlabel(,angle(vertical)) recastci(rarea) ///
plot1opts(lpattern(solid)) plot2opts(lpattern(longdash)) ///
ci1opts(lpattern(shortdash) fi(20)) ci2opts(lpattern(shortdash) fi(20)) ///
ysc(range(.1 .6)) ylabel(#6) plotopts(msymbol(none))
est restore mbase
mimrgns SAHM_child, at(year1=(0(1)4 6 8(1)17 19 20(2)38)) predict(outcome(2)) ///
post saving(satwork2, replace)
est store out2
marginsplot, title("SAHM and Employed Mothers:") ///
sub("Predicted Probability of being Moderately Satisfied with Work") xtitle("Year") ///
ytitle("Predicted Probability") level(84) ///
xlabel(,angle(vertical)) recastci(rarea) ///
plot1opts(lpattern(solid)) plot2opts(lpattern(longdash)) ///
ci1opts(lpattern(shortdash) fi(20)) ci2opts(lpattern(shortdash) fi(20)) ///
ysc(range(.1 .6)) ylabel(#6) plotopts(msymbol(none))
est restore mbase
mimrgns SAHM_child, at(year1=(0(1)4 6 8(1)17 19 20(2)38)) predict(outcome(1)) ///
post saving(satwork1, replace)
est store out1
marginsplot, title("SAHM and Employed Mothers:") ///
sub("Predicted Probability of being A Little or Very Dissatisfied with Work") xtitle("Year") ///
ytitle("Predicted Probability") level(84) ///
xlabel(,angle(vertical)) recastci(rarea) ///
plot1opts(lpattern(solid)) plot2opts(lpattern(longdash)) ///
ci1opts(lpattern(shortdash) fi(20)) ci2opts(lpattern(shortdash) fi(20)) ///
ysc(range(.1 .6)) ylabel(#6) plotopts(msymbol(none))
/*look at proportions very or moderately satisfied*/
est restore out3
matrix define A = e(b)
est restore out2
matrix define B = e(b)
matrix C = A + B
matrix list C
/*marginal effects for base model*/
est restore mbase
mimrgns, dydx(SAHM_child) at(year1=(0(1)4 6 8(1)17 19 20(2)38)) predict(outcome(3)) ///
post
est store meffs_base
/* Table 2- selection into SAHM status */
/*interactions for characteristics that you think might be becoming more or less
salient over the years */
/*mean levels of feminist norms don't predict, so rather than try to justify their
inclusion only to then take them out, just don't put them in*/
mi est, dots post: logit SAHM_child year1 polviews CP highed income_norinc ///
childhome_num babies married divorced Age10 south Black OtherRace fematt ///
year1Xpolviews year1XCP year1Xhighed year1Xincome_norinc ///
year1Xmarried year1Xdivorced year1XAge10 year1Xfematt [iweight=wtssall], vce(oim)
est store selectmod
/*test to make sure the story is the same if you limit it to just
the significant predictors (the ones you will carry forward and use as controls
in the next analysis)*/
mi est, dots post: logit SAHM_child year1 polviews CP highed income_norinc ///
childhome_num babies married divorced Age10 south Black OtherRace fematt ///
year1Xpolviews year1XCP year1Xhighed year1Xincome_norinc ///
year1Xmarried year1Xdivorced year1XAge10 year1Xfematt [iweight=wtssall], vce(oim)
est store selectmod2
mi test Black OtherRace south CP
mi est, dots post: logit SAHM_child year1 polviews highed income_norinc ///
childhome_num babies married divorced Age10 fematt ///
year1Xpolviews year1XCP year1Xhighed year1Xincome_norinc ///
year1Xmarried year1Xdivorced year1XAge10 year1Xfematt [iweight=wtssall], vce(oim)
est store selectmod3
mi test year1Xpolviews year1Xhighed year1XAge10 year1Xmarried year1Xdivorced
mi est, dots post: logit SAHM_child year1 polviews highed income_norinc ///
childhome_num babies married divorced Age10 fematt ///
year1Xincome_norinc year1Xfematt [iweight=wtssall], vce(oim)
est store selectmod4
mi est, dots post: logit SAHM_child year1 polviews highed income_norinc ///
childhome_num babies divorced Age10 fematt ///
year1Xincome_norinc year1Xfematt [iweight=wtssall], vce(oim)
est store selectmod5
/*so the story is the same even in a trimmed model, so these are the predictors
to use as controls*/
/* Table 2, model 2 - model with controls for selection into SAHM status */
mi est, dots post: ologit satjob3cat i.SAHM_child##c.year1 ///
polviews highed income_norinc childhome_num babies divorced Age10 fematt ///
year1Xincome_norinc year1Xfematt [iw=wtssall], vce(oim)
est store mselect
/* Figure 2 - marginal effects */
est restore mselect
mimrgns, dydx(SAHM_child) at(year1=(0(1)4 6 8(1)17 19 20(2)38)) predict(outcome(3)) post
est store meffs
/*95% CIs since comparing to 0, not another set of marginal means*/
marginsplot, title("Marginal Effect of SAHM status on Being Very Satisfied with Work") ///
xtitle("Year") ytitle("Marginal Effect") ///
xlabel(,angle(vertical)) ///
recastci(rarea) ci1opts(lpattern(shortdash) fcolor(gold) lcolor(gold) fintensity(20)) ///
plotopts(msymbol(none)) yline(0)
/*does the trend exist for just years for which we have the WVS housewife
attitudes data? - yes*/
mi est, dots post: ologit satjob3cat i.SAHM_child##c.year1 ///
polviews highed income_norinc childhome_num babies divorced Age10 fematt ///
year1Xincome_norinc year1Xfematt [iw=wtssall] if !missing(housewife_fulfil), vce(oim)
est store mselect_hwatt
/* Table 3 - this model was mistakenly reported for Table 3 in the paper -
the substantive conclusions are the same, but it uses married as a control
instead of divorced */
mi est, dots post: ologit satjob3cat i.SAHM_child##c.housewife_fulfil year1 ///
polviews highed income_norinc married childhome_num babies Age10 fematt ///
year1Xincome_norinc year1Xfematt ///
[iw=wtssall] if !missing(housewife_fulfil), vce(oim)
est store mselect_hwatt2_reported
/* this is the actual model with the correct control, and the one used to generate
the numbers and marginal effects reported in the text (not table) */
mi est, dots post: ologit satjob3cat i.SAHM_child##c.housewife_fulfil year1 ///
polviews highed income_norinc childhome_num babies divorced Age10 fematt ///
year1Xincome_norinc year1Xfematt ///
[iw=wtssall] if !missing(housewife_fulfil), vce(oim)
est store mselect_hwatt2
mimrgns, dydx(SAHM) predict(outcome(3)) at(housewife_fulfil=(21 32)) post
est store meffs_hwatt
/*does the trend differ by class?*/
/* Table 1, model 3 - family income */
mi est, dots post: ologit satjob3cat i.SAHM_child##c.year1 ///
polviews highed income_norinc childhome_num babies divorced Age10 fematt FamilyIncome ///
year1Xincome_norinc year1Xfematt year1XFamilyIncome ///
SAHM_childXFamilyIncome SAHMXyearXFIncome ///
[iw=wtssall], vce(oim)
est store mincome
/*how much does this change if I use a non-imputed interaction (so I can
use margins?) - some differences, but substantively similar, though results
are cleaner and stronger when using imputed interactions, as expected - see footnote
8 in published article*/
mi est, dots post: ologit satjob3cat i.SAHM_child##c.year1##c.FamilyIncome ///
polviews highed income_norinc childhome_num babies divorced Age10 fematt ///
year1Xincome_norinc year1Xfematt ///
[iw=wtssall], vce(oim)
est store mincome2
est restore mincome2
mimrgns, dydx(SAHM_child) at(year1=(0(1)4 6 8(1)17 19 20(2)38) (p25) FamilyIncome) ///
predict(outcome(3)) post saving(meff_income_lo, replace)
est store meffs_educ1
marginsplot, title("Marginal Effects of SAHM status on Being Very Satisfied with Work") ///
sub("for women at 1st and 3rd Quartiles of Household Income") xtitle("Year") ///
ytitle("Marginal Effect") xlabel(,angle(vertical)) recastci(rarea) ///
ci1opts(lpattern(shortdash) fcolor(red) lcolor(red) fintensity(20)) ///
plotopts(msymbol(none)) yline(0) level(84) ///
plot1opts(lpattern(solid))
est restore mincome2
mimrgns, dydx(SAHM_child) at(year1=(0(1)4 6 8(1)17 19 20(2)38) (p75) FamilyIncome) ///
predict(outcome(3)) post saving(meff_income_hi, replace)
est store meffs_educ2
marginsplot, title("Marginal Effects of SAHM status on Being Very Satisfied with Work") ///
sub("for women at 1st and 3rd Quartiles of Household Income") xtitle("Year") ///
ytitle("Marginal Effect") xlabel(,angle(vertical)) recastci(rarea) ///
ci1opts(lpattern(shortdash) fcolor(blue) lcolor(blue) fintensity(20)) ///
plotopts(msymbol(none)) yline(0) level(84) ///
plot1opts(lpattern(solid))
/*Table 1, model 4 - education */
/*no missing on interaction terms, so use Stata interaction lingo
so that I can use margins*/
mi est, dots post: ologit satjob3cat i.SAHM_child##c.year1##i.highed ///
polviews income_norinc childhome_num babies divorced Age10 fematt ///
year1Xincome_norinc year1Xfematt [iw=wtssall], vce(oim)
est store meduc
est restore meduc
mimrgns highed, dydx(SAHM_child) at(year1=(0(1)4 6 8(1)17 19 20(2)38)) predict(outcome(3)) post
est store meffs_educ
marginsplot, title("Marginal Effects of SAHM status on Being Very Satisfied with Work") ///
sub("for women with and without a college degree") xtitle("Year") ///
ytitle("Marginal Effect") xlabel(,angle(vertical)) recastci(rarea) ///
ci1opts(lpattern(shortdash) fcolor(blue) lcolor(blue) fintensity(20)) ///
ci2opts(lpattern(shortdash) fcolor(red) lcolor(red) fintensity(20)) ///
plotopts(msymbol(none)) yline(0) level(84) ///
plot1opts(lpattern(solid)) plot2opts(lpattern(longdash))
/*so yes, trend is more positive for those in higher classes, as predicted */
/*college educated SAHMs seem to be more satisfied than employed women in the 2000s -
is this just an artifact of using a straight line, which might increase the extremity
of predicitons at the ends of the line?*/
gen SAHM_highed = SAHM_child==1 & highed==1
gen SAHM_lowed = SAHM_child==1 & highed==0
mi update
/*examine the average of the last two years from the trend model vs.
a model that uses just those two years*/
est replay meffs_educ
di (.1084367 +.0942774 )/2
mi est, dots post: ologit satjob3cat i.SAHM_child year1 ///
polviews highed income_norinc childhome_num babies divorced Age10 fematt ///
year1Xincome_norinc year1Xfematt [iw=wtssall] if year1 >= 36 & highed==1, vce(oim)
est store meduc_hi
mimrgns, dydx(SAHM) predict(outcome(3)) post
est store meffs_meduc_hi
/*effect size is very similar, so the trend line prediction is plausible
not statistically significant, but that is not surprising given how few
college-educated SAHMs there are in just those two years*/
/***SENSITIVITY CHECKS***/
/*** Table 4 - same results if we control for possible misunderstanding of the job satisfaction
question? ***/
/*limit sample to those who have never worked - more likely to interpret the
question as having to do with housework*/
/*include all working mothers, and those SAHMs who have never worked*/
mi est, dots post: ologit satjob3cat i.SAHM_child##c.year1 ///
polviews highed income_norinc childhome_num babies divorced Age10 fematt ///
year1Xincome_norinc year1Xfematt [iw=wtssall] if SAHM_child==0 | everwork==0, vce(oim)
est store mselect2
mi est, dots post: ologit satjob3cat i.SAHM_child##c.year1 ///
polviews highed income_norinc childhome_num babies divorced Age10 fematt FamilyIncome ///
year1Xincome_norinc year1Xfematt year1XFamilyIncome SAHM_childXFamilyIncome ///
SAHMXyearXFIncome [iw=wtssall] if SAHM_child==0 | everwork==0, vce(oim)
est store mincome2
mi est, dots post: ologit satjob3cat i.SAHM_child##c.year1##i.highed ///
polviews income_norinc childhome_num babies divorced Age10 fematt ///
year1Xincome_norinc year1Xfematt [iw=wtssall] if SAHM_child==0 | everwork==0, vce(oim)
est store meduc2
/*sensitivity analysis not reported in the paper*/
/*same results with OLS models with robust SE's?*/
/*this is to check whether the bias from omitted variables in the logit
models is too great to warrant using logit models*/
mi est, dots post: reg satjob3cat i.SAHM_child##c.year1 ///
polviews CP highed FamilyIncome childhome_num married ///
Age10 babies fematt fematt_mean year1Xpolviews year1XFamilyIncome ///
year1Xfematt [iw=wtssall], vce(robust)
est store mselect_ols
mi est, dots post: reg satjob3cat i.SAHM_child##c.year1 ///
polviews CP highed FamilyIncome childhome_num married ///
Age10 babies fematt fematt_mean year1Xpolviews year1XFamilyIncome ///
year1Xfematt SAHM_childXFamilyIncome ///
SAHMXyearXFIncome [iw=wtssall], vce(robust)
est store mincome_ols
mi est, dots post: reg satjob3cat i.SAHM_child##c.year1##i.highed ///
polviews CP highed FamilyIncome childhome_num married ///
Age10 babies fematt fematt_mean year1Xpolviews year1XFamilyIncome ///
year1Xfematt [iw=wtssall], vce(robust)
est store meduc_ols
/*all very similar*/
/*Appendix - descriptives - full sample*/
proportion satjob3cat [pw=wtssall] /*no missing so don't use mi prefix*/
mi est, dots post: mean SAHM_child polviews CP highed FamilyIncome income_norinc ///
childhome_num babies married divorced Age10 south Black OtherRace fematt [pw=wtssall]
/*descriptives - by employment status*/
proportion satjob3cat [pw=wtssall], over(SAHM_child)
mi est, dots post: mean polviews CP highed FamilyIncome income_norinc ///
childhome_num babies married divorced Age10 south Black OtherRace ///
fematt [pw=wtssall], over(SAHM_child)
/*percentage missing*/
misstable sum satjob3cat SAHM_child polviews CP highed FamilyIncome income_norinc ///
childhome_num babies married divorced Age10 south Black OtherRace fematt