st: RE: RE: RE: two-parts model bootstrap using XTGEE

[tshmak <tshmak@hku.hk>]

Hi Zhonghe, 

Sorry maybe I was wrong in saying you can run -margins- to get the effects you wanted, cos it's a combination of predictions from two models. But from your original post it seems that the problem you had was that convergence failed in one of the bootstrap samples for one of your xtgee models. The best thing is probably to work out exactly why the convergence failed. -bootstrap- comes with the noisily and trace subcommand, so you can use these to help you. Typical reasons for failure in convergence are perhaps: you encountered a flat region. This could be because you have variables that are nearly perfectly collinear. Bear in mind that even if your original data is fine, your bootstrap samples (or at least some of them) could have problems. Stata also automatically removes perfectly collinear variables. So it could be that in one of your bootstrap samples, one of your i.post, i.pilot, or i.post_pilot got removed. 

I think by default -bootstrap- continues even if the sample results in a missing value. So, try running it with the noisily option, and see what error message it gives you. 

Hope it helps. 

Tim




-----Original Message-----
From: owner-statalist@hsphsun2.harvard.edu [mailto:owner-statalist@hsphsun2.harvard.edu] On Behalf Of Li, Elena (Zhonghe)
Sent: 04 December 2012 12:53
To: statalist@hsphsun2.harvard.edu
Subject: st: RE: RE: two-parts model bootstrap using XTGEE

Thank you Tim.

If it were bootstrapping from one model, I know how to get the marginal effect - I assume that was what you meant?
And I think I have the formula to get the marginal effects for the two parts - but we need the p value for this difference in difference of p0*p1 that I don't know how to get without bootstrapping.

Best,
Zhonghe
________________________________________
From: owner-statalist@hsphsun2.harvard.edu [owner-statalist@hsphsun2.harvard.edu] on behalf of tshmak [tshmak@hku.hk]
Sent: Monday, December 03, 2012 10:11 PM
To: statalist@hsphsun2.harvard.edu
Subject: st: RE: two-parts model bootstrap using XTGEE

Dear Zhonghe,

Instead of writing the program the way you wrote below, I think you should consider using the command -margins-. That would simplify your program a lot and let's see if you still get the same problem.

Best,

Tim

-----Original Message-----
From: owner-statalist@hsphsun2.harvard.edu [mailto:owner-statalist@hsphsun2.harvard.edu] On Behalf Of Li, Elena (Zhonghe)
Sent: 04 December 2012 03:12
To: statalist@hsphsun2.harvard.edu
Subject: st: two-parts model bootstrap using XTGEE

Hi
We are fitting a two-parts model on health care cost using longitudinal data.
Part1 is modeling whether there is any cost, and
Part2 is modeling the cost when if there is a cost.

The goal is to bootstrap the difference in difference (post_pilot) in terms of the p0*p1 (the predicted probability of having cost * the predicted cost).
I am a SAS user and tried to put the following Stata program together.  The models converged the first time but when the replication starts the modeling does not converge:
"convergence not achieved
an error occurred when bootstrap executed tempa, posting missing values"

Can you please help?  Thank you.
Zhonghe Li

_____________________________________________________________
sort member qtr
by member qtr: drop if _n~=1
egen memberid=group(member)
egen practice=group(pcp)

gen newid=memberid
xtset newid qtr
generate sample=(1-missing(havecost, logcost, female, aet, uhc, ant,  Qtr_months, age1-age4, qtr, health,  post, pilot, post_pilot)) keep if sample

matrix define results = (1,2,3)


capture program drop tempa

program define tempa, rclass


xtgee havecost i.female aet uhc ant  Qtr_months age1-age4 qtr health  i.post i.pilot i.post_pilot,  family(binomial) link(logit) corr(ar 1)


        gen truepost=post
        gen truepilot=pilot
        gen true_post_pilot=post_pilot

        replace post=0
        replace pilot=0
        replace post_pilot=0
        predict pre_control_p1

        replace pilot=1
        predict pre_pilot_p1

        replace post=1
        replace post_pilot=1
        predict post_pilot_p1

        replace pilot=0
        replace post_pilot=0
        predict post_control_p1

        replace post=truepost
        replace pilot=truepilot
        replace post_pilot=true_post_pilot

xtgee logcost female aet uhc ant  Qtr_months age1-age4 qtr health  post pilot post_pilot if logcost>0,  family(gamma) link(log) corr(ar 1) vce(robust)


        replace post=0
        replace pilot=0
        replace post_pilot=0
        predict pre_control_p2, mu

        replace pilot=1
        predict pre_pilot_p2, mu

        replace post=1
        replace post_pilot=1
        predict post_pilot_p2, mu

        replace pilot=0
        replace post_pilot=0
        predict post_control_p2, mu

        drop post pilot post_pilot
        rename truepost post
        rename truepilot pilot
        rename true_post_pilot post_pilot


gen phat_pre_control=pre_control_p1*pre_control_p2
gen phat_pre_pilot=pre_pilot_p1*pre_pilot_p2
gen phat_post_pilot=post_pilot_p1*post_pilot_p2
gen phat_post_control=post_control_p1*post_control_p2

gen diff_pilot=phat_post_pilot-phat_pre_pilot
gen diff_control=phat_post_control-phat_pre_control
gen diff_in_diff=diff_pilot-diff_control

summ diff_pilot
return scalar m_diff_pilot = r(mean)
local var1 = r(mean)

summ diff_control
return scalar m_diff_control = r(mean)
local var2 = r(mean)

summ diff_in_diff
return scalar m_diff_in_diff = r(mean)
local var3 = r(mean)

matrix define temp=(`var1',`var2',`var3') matrix define results = results \ temp

matrix list results

quietly: drop diff*  pre_* post_control_p1 post_control_p2 post_pilot_p1 post_pilot_p2 phat*
*quietly: scalar drop _all

xtset,clear

end



bootstrap mean_diff_pilot=r(m_diff_pilot) mean_diff_control=r(m_diff_control) mean_diff_in_diff=r(m_diff_in_diff), noisily strata(pilot) cluster(memberid) idcluster(newid) group(qtr) reps(100) seed(1234567) : tempa

svmat results, names(diff)
outsheet twopartdiff* using diff.csv, comma replace

*
*   For searches and help try:
*   http://www.stata.com/help.cgi?search
*   http://www.stata.com/support/faqs/resources/statalist-faq/
*   http://www.ats.ucla.edu/stat/stata/

*
*   For searches and help try:
*   http://www.stata.com/help.cgi?search
*   http://www.stata.com/support/faqs/resources/statalist-faq/
*   http://www.ats.ucla.edu/stat/stata/

*
*   For searches and help try:
*   http://www.stata.com/help.cgi?search
*   http://www.stata.com/support/faqs/resources/statalist-faq/
*   http://www.ats.ucla.edu/stat/stata/

*
*   For searches and help try:
*   http://www.stata.com/help.cgi?search
*   http://www.stata.com/support/faqs/resources/statalist-faq/
*   http://www.ats.ucla.edu/stat/stata/