Bootstrap bias corrected confidence intervals in multiply imputed weighted multinomial regression models

[Darina Peycheva]

Hi all,

I'd appreciate your help with this. I'm trying to obtain bootstrap bias corrected confidence intervals from multinomial regression model run on imputed data, also weighted. This is my code for the multinomial model:

mi estimate, eform: mlogit outcome exposure covariates [pweight], rrr b(1)

I'm trying the following code with added bootstrap bias corrected confidence interval:

program define mymlogit
mlogit outcome exposure covariates [pweight], rrr b(1)
end

use "mydata.dta", clear

*mi set mlong
mi convert mlong

mi xeq: bootstrap, reps(100) bca: mymlogit

estat bootstrap, bca

And I'm getting the following error message:

'insufficient observations to compute jackknife standard errors'

But I also think the issue may be with the code and not just the small number of replications / small sample size. My dataset contains 1,158 cases; no rare outcomes.

Thanks for your help.

Darina

[Felix Bittmann]

First, are you sure you need BCa instead of only BC?
Second, combining bootstrapping and imputation is not straightforward and your setup will probably not work. I have outlined two options to do this in Stata, which is described here: https://www.preprints.org/manuscript...13/v1/download

[Darina Peycheva]

Thanks very much, this was very helpful!

Yes, sorry, this was BCCI. I managed to adapt an old code of mine where I used the same regression model as part of mediation analysis. The code now looks as follows and runs smoothly:

capture program drop mymlogit
program mymlogit, rclass

*Estimate TE
mi estimate, saving(miest02,replace) eform: mlogit outcome exposure covariates [pweight=r3lhwtresp], b(1) rrr

*matrix bb_total = e(b_mi)
scalar b_total1 = bb_total[1,27]
*scalar list b_total1
return scalar b_total1 = bb_total[1,27]

*return list
*display r(b_total1)
*matrix list bb_total

scalar b_total2 = bb_total[1,53]
*scalar list b_total2
return scalar b_total2 = bb_total[1,53]

end

bootstrap exp(r(b_total1)) exp(r(b_total2)) , seed(12345) reps(5) noisily: mymlogit
estat bootstrap,all