computing SE with survey weights for the Arhomme command

[Luis Mijares Castaneda]

Hello, I have the following problem, i want to use the survey stratification and psu using the -arhomme- command, I first tried the following code and received the following error, "arhomme is not supported by svy with vce(bootstrap); see help svy estimation for a list of Stata estimation commands that are supported by svy r(322);" I then tried writing the program in the second code block but for some reason that program does not compile, any help for how to use svyset with arhomme would be greatly appreciated.

svyset raehsamp [pweight=new_weight], strata (raestrat)
bsweights bs_, n(-1) reps(100)seed(4881269) 
svyset [pw=new_weight], bsrw(bs_*)
xi: svy bootstrap, nodrop _b: arhomme log_avrg_cost i.inc_d endentulism race age_cat ///
   male education veteran mothered wealth smoke_now ///
        chronicdisease, ///
        select(r11dentst = dentalinsurance_w1 endentulism ///
        inc_d race age_cat male education veteran mothered wealth ///
        smoke_now chronicdisease) quantiles(0.5) taupoints(20) rhopoints(49) ///
        meshsize(1) graph nostderrors gaussian

arhomme is not supported by svy with vce(bootstrap); see help svy estimation for a list of Stata estimation commands that are supported by svy
r(322);

 


cap program drop boot_arhomme
program define boot_arhomme, eclass
    preserve
    * Resample data while keeping PSU structure (survey design)
    bsample, cluster(raehsamp) strata(raestrat)  

    * Run arhomme with probability weights
    quietly xi:arhomme log_avrg_cost i.inc_d i.endentulism i.race i.age_cat ///
        i.male i.education i.veteran i.mothered i.wealth i.smoke_now ///
        chronicdisease [pw=new_weight], ///
        select(r11dentst = dentalinsurance_w1 endentulism ///
        inc_d race age_cat male education veteran mothered wealth ///
        smoke_now chronicdisease) quantiles(0.5) taupoints(20) rhopoints(49) ///
        meshsize(1) graph nostderrors gaussian

    * Save bootstrapped coefficients
    return scalar b_inc_d = _b[inc_d]
    return scalar b_race  = _b[race]
    return scalar b_edu   = _b[education]

    restore
end


* Run bootstrap with 1000 replications
simulate b_inc_d=r(b_inc_d) b_race=r(b_race) b_edu=r(b_edu), reps(1000) seed(12345): boot_arhomme

* Compute bootstrapped standard errors
summarize b_inc_d b_race b_edu

* Compute bootstrapped 95% confidence intervals
centile b_inc_d b_race b_edu, centile(2.5 97.5)