Code for bootstrapped CI for inference about mediation

[Jon Hoefer]

Dear Statalists,

below I have run a Monte Carlo simulation in order to obtain CI's which I can use to infer wether or not the indirect effect of my mediation analysis is significant. However, I would prefer doing the same with bootstrapped CI's. Could anyone be so kind to provide a code / hint how I can code it for the same procedure:
The coefficients b4 and b5 are the coefficients from the two regressions conducted "along the indirect path" for the mediation. Standard errors are marked in red and coefficients in blue, to make it clearer to understand what I have coded:

. set obs 5000
number of observations (_N) was 2,651, now 5,000

. gen b4 = (invnorm(uniform())*.0329509)+-.0263693

. gen b5 = (invnorm(uniform())*.1144345)+-.0705507

. gen CI_1D = b4*b5

. centile CI_1D, centile(2.5 97.5)

                                                       -- Binom. Interp. --
    Variable |       Obs  Percentile    Centile        [95% Conf. Interval]
-------------+-------------------------------------------------------------
       CI_1D |     5,000        2.5   -.0072381       -.0077505   -.0067966
             |                 97.5     .015011        .0142704    .0158532

Thank you in advance for your suggestions

[Mike Lacy]

The conventional way to bootstrap a sample estimate is to write a program that calculates that estimate, puts it in either the return list or the ereturn list, after which you feed that program to -bootstrap-. It appears that you are doing something different. Is there a reason why the conventional method doesn't work for you here? By "conventional," I'd mean something like this:

cap prog drop mymediation
prog mymediation, rclass
...run regression models to calculate a mediation coefficient by some method
return b_med = ( ... whatever coefficient you calculated)
end
bootstrap b_med = r(b_med), reps(5000): mymediation

This presumes also that the mediation calculation you want can't be done by one of the user-written mediation commands, such as -khb- (ssc describe khb).

[Jon Hoefer]

Thanks for your response Mike.

Frankly I have not worked with the khb tool before, as I would prefer using a code for my thesis that does not rely on non-default commands. Also for me the intuition behind the "manual" calculation of the bootstrapped CI is very intuitive and hence I decided to go with it over some of the individual tools. When considering your code, some questions arise:
1. When I code :

 cap prog drop mymediation

prog mymediation, rclass
  1. run reg ROA Premium1DW Family_Firm_Identifier Founder_Identity FirmAge FirmSize Indebtedness i.SIC_2 i.Fiscal_Year, vce(cluster id)

it gives me

option vce() not allowed

How come and what could I do about it? I need clustered errors for my estimation
2. the indirect effect is the product of two separate regression's coefficients, how can I get the two into the code you provided? Is there a way of just manually plugging them in into the code after running the two regressions? If so, this might be easier for someone with my Stata programming knowledge.

Edit: Another thing I have tried to get to work is the SEM equation builder. However, I wasn't sure how to include the control variables and so again decided to go with the "manual" way

[Jon Hoefer]

Dear Mike,

I have tried to get the khb to work and in fact replicated the authors code from their example, but do not get it to work:

khb regress ROA Family_Firm_Identifier  || Acquisition_Frequency , concomitant(Founder_Identity FirmAge FirmSize Growthopp Indebtedness i.SIC_2 i.Fiscal_Year) cluster(GKVEY) summary
invalid stub() option;
stub name too long
r(198);

do you have a clue about what is going on / wrong. I even tried the code with any options (that is without any control variable), but even then the error remains the same.
Thank you in advance for your help!
EDIT: Fixed it myself, apparently variable names are were too long. After trimming Family_Firm_Identifier and Acquisition_Frequency I got it to work. Attached you find the output.
. khb regress ROA FFI || ACF , concomitant(Founder_Identity FirmAge FirmSize Growthopp Indebtedness i.SIC_2 i.Fiscal_Year) sum
> mary vce(cluster GVKEY)

Decomposition using Linear Probability Models

Model-Type: regress Number of obs = 1802
Variables of Interest: FFI R-squared = 0.25
Z-variable(s): ACF
Concomitant: Founder_Identity FirmAge FirmSize Growthopp Indebtedness i.SIC_2 i.Fiscal_Year
------------------------------------------------------------------------------
| Robust
ROA | Coef. Std. Err. z P>|z| [95% Conf. Interval]
-------------+----------------------------------------------------------------
FFI |
Reduced | .0286796 .0129134 2.22 0.026 .0033698 .0539894
Full | .0287867 .012922 2.23 0.026 .0034601 .0541133
Diff | -.0001071 .0002396 -0.45 0.655 -.0005768 .0003626
------------------------------------------------------------------------------

Summary of confounding

Variable | Conf_ratio Conf_Pct Resc_Fact
-------------+-------------------------------------
FFI | .99627919 -0.37 1
---------------------------------------------------
Two things: First as Diff p>>0.05 there is no indirect effect?! Second: Can someone maybe in a few words elaborate what exactly the estimator does (as I will have to describe it in my thesis?) Does it use bootstrapping, too?

[Mike Lacy]

If I recall correctly, for a regular regression model, I believe that -khb- gives the same results as other older mediation techniques. As the help file shows, there are quite a few articles associated with -khb-, so I think you'll just have to dive into that literature and the documentation for -khb-.

Relatively few estimation methods in Stata incorporate bootstrapping into the estimation command itself. Per my comments above, the conventional way to bootstrap in Stata is to use your program with Stata's -bootstrap- prefix, which will work with almost any estimation program in Stata, user-written or built-in. You'll have to extract the coefficient you want from the ereturn list of -khb-, probably from e(b). See the help for khb and bootstrap.

[Mike Zyphur]

Hi Joe
Just came across your post after replying to someone on a similar topic: the medsem package might be useful for you. It's a Monte Carlo (parametric) bootstrap approach that is faster than non-parametric bootstrapping, which can be a significant time saver with complex models. Perhaps it can cover your specific use case?

In case it's relevant, Instats is offering a series of seminars on Stata's sem command including bootstrapped CIs for indirect effects and conditional indirect effects. All seminars are combined in a Structured Course in case it's of interest.

Best wishes
Mike