Hi,
I'd be grateful if anyone can help me figure out where the problem stands with the my attempt to obtain bootstrap CIs for total, direct, and indirect effects (for the association between smoking during pregnancy (of mothers) and time to natural menopause (of daughters)).
I'm trying to do mediation analysis using inverse odds weighting. More about the method can be found here: 10.1093/aje/kwu278. I've adapted the code used here: dx.doi.org/10.1136/jech-2016-208671; dx.doi.org/10.1136/bmjopen-2018-026258.
The exposure (smoking during pregnancy) is coded 0/1. The outcome (experience of natural menopause) is coded 0/1. Time to event is in years since age 11. The joint mediating role of birthweight standardized for gestation (z-score) and breastfeeding (0/1) is of interest. Covariates include maternal age at birth (in years), parity (continuous), and social class (manual/nonmanual coded 0/1).
This is what my code looks like:
*performs IOW with stcox and applies weightscapture program drop IOW
capture program drop IOWMI
program IOWMI, rclass
capture drop logodds predprob inverseodds wt_iow smokingm_mi smokingm_new
*Run exposure model
mi estimate, saving(miest01,replace): logit smokingm stdbweight i.breastfed ib3.manual mage parity
*Create inverse odds weights
mi predict logodds using miest01
mi xeq: gen predprob=exp(logodds)/(1+exp(logodds))
mi xeq: gen inverseodds = ((1-predprob)/predprob)
*This deals with missingness in the original 'smokingm' variable
order _1_smokingm _2_smoking _3_smoking _4_smoking _5_smoking
egen smokingm_mi=rowmean( _1_smokingm- _5_smokingm)
*mean(smokingm_mi)
*list _1_smokingm- _20_smokingm smokingm_mi in 1/100 if smokingm==.
*tabstat smokingm_mi,stat(mean sd p50 min max)
*histogram smokingm_mi
*tab smokingm_mi
gen smokingm_new=1 if smokingm_mi>=0.5
replace smokingm_new=0 if smokingm_mi<0.5
*tab smokingm_new
*tab smokingm_new if smokingm==.
*tab smokingm_new if smokingm!=.
mi xeq: gen wt_iow = 1 if smokingm_new==0
mi xeq: replace wt_iow = inverseodds if smokingm_new==1
*Estimate TE
mi stset timeout, failure(meno) enter(timein) origin(timein)
mi estimate, saving(miest02,replace): stcox i.smokingm ib3.manual mage parity
matrix bb_total = e(b_mi)
scalar b_total = bb_total[1,2]
return scalar b_total = bb_total[1,2]
*Estimate NDE
mi stset timeout [pweight=wt_iow], failure(meno) enter(timein) origin(timein)
mi estimate, saving(miest03,replace): stcox i.smokingm ib3.manual mage parity
matrix bb_direct = e(b_mi)
scalar b_direct = bb_direct[1,2]
return scalar b_direct = bb_direct[1,2]
*Estimate NIE
return scalar b_indirect = b_total-b_direct
end
* Request bootstrapped estimates of indirect, direct and total effects.
bootstrap r(b_indirect) r(b_direct) r(b_total), seed(12345) reps(5) noisily: IOWMI
The error that I get is as follows:
. IOWMI
variable wt_iow not found
Your mi data are stset and some of the variables previously declared by stset are not in the dataset. mi verifies that none of the stset variables are also
registered as imputed or passive. Type mi stset, clear to clear old no-longer-valid settings.
an error occurred when bootstrap executed IOWMI, posting missing values
I've read several posts about issues with MI and bootstrap, survival and bootstrap, weighting and bootstrap, but have not seen an example where all this complexity is applied simultaneously.
I'd appreciate your help.
Darina
---------------------------------------------------------------
This is the output from dataex if this could help:
As I cannot share an example of the imputed data as there are too many varibales. The following should be run:
mi set wide
mi register imputed nmeno manual mage parity smokingm stdbweight breastfed
mi set M=5
mi impute chained (logit) nmeno smokingm manual breastfed (regress) stdbweight mage parity, replace rseed(678910) noisily
*run to enable setup of survival analysis
*add date of birth
gen dob_m=3
gen dob_y=1970
gen dob=ym(dob_y,dob_m)
format dob %tm
*add 11th birthday (origin, start of follow-up)
gen age11b_m=3
gen age11b_y=1981
gen age11b=ym(age11b_y,age11b_m)
format age11b %tm
gen age11=(age11b-dob)/12
*add age at event (end of follow-up)
gen event_age=(event_date_ym - dob)/12
gen meno=nmeno
gen timein=age11
gen timeout=event_age
I'd be grateful if anyone can help me figure out where the problem stands with the my attempt to obtain bootstrap CIs for total, direct, and indirect effects (for the association between smoking during pregnancy (of mothers) and time to natural menopause (of daughters)).
I'm trying to do mediation analysis using inverse odds weighting. More about the method can be found here: 10.1093/aje/kwu278. I've adapted the code used here: dx.doi.org/10.1136/jech-2016-208671; dx.doi.org/10.1136/bmjopen-2018-026258.
The exposure (smoking during pregnancy) is coded 0/1. The outcome (experience of natural menopause) is coded 0/1. Time to event is in years since age 11. The joint mediating role of birthweight standardized for gestation (z-score) and breastfeeding (0/1) is of interest. Covariates include maternal age at birth (in years), parity (continuous), and social class (manual/nonmanual coded 0/1).
This is what my code looks like:
*performs IOW with stcox and applies weightscapture program drop IOW
capture program drop IOWMI
program IOWMI, rclass
capture drop logodds predprob inverseodds wt_iow smokingm_mi smokingm_new
*Run exposure model
mi estimate, saving(miest01,replace): logit smokingm stdbweight i.breastfed ib3.manual mage parity
*Create inverse odds weights
mi predict logodds using miest01
mi xeq: gen predprob=exp(logodds)/(1+exp(logodds))
mi xeq: gen inverseodds = ((1-predprob)/predprob)
*This deals with missingness in the original 'smokingm' variable
order _1_smokingm _2_smoking _3_smoking _4_smoking _5_smoking
egen smokingm_mi=rowmean( _1_smokingm- _5_smokingm)
*mean(smokingm_mi)
*list _1_smokingm- _20_smokingm smokingm_mi in 1/100 if smokingm==.
*tabstat smokingm_mi,stat(mean sd p50 min max)
*histogram smokingm_mi
*tab smokingm_mi
gen smokingm_new=1 if smokingm_mi>=0.5
replace smokingm_new=0 if smokingm_mi<0.5
*tab smokingm_new
*tab smokingm_new if smokingm==.
*tab smokingm_new if smokingm!=.
mi xeq: gen wt_iow = 1 if smokingm_new==0
mi xeq: replace wt_iow = inverseodds if smokingm_new==1
*Estimate TE
mi stset timeout, failure(meno) enter(timein) origin(timein)
mi estimate, saving(miest02,replace): stcox i.smokingm ib3.manual mage parity
matrix bb_total = e(b_mi)
scalar b_total = bb_total[1,2]
return scalar b_total = bb_total[1,2]
*Estimate NDE
mi stset timeout [pweight=wt_iow], failure(meno) enter(timein) origin(timein)
mi estimate, saving(miest03,replace): stcox i.smokingm ib3.manual mage parity
matrix bb_direct = e(b_mi)
scalar b_direct = bb_direct[1,2]
return scalar b_direct = bb_direct[1,2]
*Estimate NIE
return scalar b_indirect = b_total-b_direct
end
* Request bootstrapped estimates of indirect, direct and total effects.
bootstrap r(b_indirect) r(b_direct) r(b_total), seed(12345) reps(5) noisily: IOWMI
The error that I get is as follows:
. IOWMI
variable wt_iow not found
Your mi data are stset and some of the variables previously declared by stset are not in the dataset. mi verifies that none of the stset variables are also
registered as imputed or passive. Type mi stset, clear to clear old no-longer-valid settings.
an error occurred when bootstrap executed IOWMI, posting missing values
I've read several posts about issues with MI and bootstrap, survival and bootstrap, weighting and bootstrap, but have not seen an example where all this complexity is applied simultaneously.
I'd appreciate your help.
Darina
---------------------------------------------------------------
This is the output from dataex if this could help:
Code:
* Example generated by -dataex-. For more info, type help dataex clear input float(event_date_ym nmeno) byte mage float stdbweight byte(manual smokingm parity breastfed) 702 0 34 -.210943 0 0 4 0 598 0 25 1.0130603 0 0 1 0 611 0 . . . . . . 702 0 18 1.4068865 1 0 0 . 702 0 30 -1.0437096 1 1 1 0 702 0 18 -2.1919973 1 0 0 . 702 0 26 -1.048419 0 0 0 1 677 0 22 .3473743 1 1 1 1 702 0 20 -.9111431 1 1 0 0 702 0 24 -.654446 1 0 5 . 702 0 29 -2.2724223 1 1 0 0 702 0 24 -.14652179 1 1 0 0 619 1 31 -.710833 1 1 1 1 702 0 31 .3924279 1 0 1 1 702 0 22 .454743 0 1 0 0 593 1 . . . . . . 702 0 26 -.29078162 1 1 2 . 702 0 26 1.4280024 0 1 2 0 702 0 29 . 1 0 2 0 702 0 36 .0467419 1 1 4 1 end format %tm event_date_ym label values nmeno nmeno label def nmeno 0 "no", modify label def nmeno 1 "yes", modify label values mage VAR00007 label values manual manual label def manual 0 "non-manual", modify label def manual 1 "manual / other", modify label values parity A0166 label values breastfed breastfed label def breastfed 0 "<1m", modify label def breastfed 1 "1+ months", modify
mi set wide
mi register imputed nmeno manual mage parity smokingm stdbweight breastfed
mi set M=5
mi impute chained (logit) nmeno smokingm manual breastfed (regress) stdbweight mage parity, replace rseed(678910) noisily
*run to enable setup of survival analysis
*add date of birth
gen dob_m=3
gen dob_y=1970
gen dob=ym(dob_y,dob_m)
format dob %tm
*add 11th birthday (origin, start of follow-up)
gen age11b_m=3
gen age11b_y=1981
gen age11b=ym(age11b_y,age11b_m)
format age11b %tm
gen age11=(age11b-dob)/12
*add age at event (end of follow-up)
gen event_age=(event_date_ym - dob)/12
gen meno=nmeno
gen timein=age11
gen timeout=event_age
Comment