Hi everyone,
I've made some initial ventures into Stata's bayesmh command, and I'm tying to get to grips with the commands potential, and limitations.
I have been trying to replicate a panel regression I ran, in a Bayesian framework. However, the syntax is really quite different, and I'm not sure I can pull out exactly what I'm looking for from the model. If I run through the example:
I have a simple random-effects regression which tests for treatment effects between 2 groups of individuals (treatgroup, either 0 or 1), with 2 follow-ups (visit, either 1 or 2):
Code:
. xtreg x visit##treatgroup, i(id) re
Random-effects GLS regression Number of obs = 200
Group variable: id Number of groups = 100
R-sq: Obs per group:
within = 0.9386 min = 2
between = 0.7622 avg = 2.0
overall = 0.8649 max = 2
Wald chi2(3) = 1811.87
corr(u_i, X) = 0 (assumed) Prob > chi2 = 0.0000
----------------------------------------------------------------------------------
x | Coef. Std. Err. z P>|z| [95% Conf. Interval]
-----------------+----------------------------------------------------------------
2.visit | -.2580899 .1966148 -1.31 0.189 -.6434478 .127268
1.treatgroup | .1748543 .2701736 0.65 0.518 -.3546761 .7043848
|
visit#treatgroup |
2 1 | 7.86268 .2780553 28.28 0.000 7.317701 8.407658
|
_cons | -.2941914 .1910416 -1.54 0.124 -.668626 .0802432
-----------------+----------------------------------------------------------------
sigma_u | .92650385
sigma_e | .98307397
rho | .47040156 (fraction of variance due to u_i)
----------------------------------------------------------------------------------
The value I want to extract from the model above is the change in each treatment group over time, which for the first group is _b[2.visit], and I can obtain the other group's one using lincom 1.treatgroup + 2.visit#1.treatgroup. Alternatively, I know I can also use margins visit, dydx(treatgroup) to get this same information:
Code:
. margins visit, dydx(treatgroup)
Conditional marginal effects Number of obs = 200
Model VCE : Conventional
Expression : Linear prediction, predict()
dy/dx w.r.t. : 1.treatgroup
-------------------------------------------------------------------------------
| Delta-method
| dy/dx Std. Err. z P>|z| [95% Conf. Interval]
--------------+----------------------------------------------------------------
1.treatgroup |
visit |
1 | .1748543 .2701736 0.65 0.518 -.3546761 .7043848
2 | 8.037534 .2701736 29.75 0.000 7.508004 8.567064
-------------------------------------------------------------------------------
Note: dy/dx for factor levels is the discrete change from the base level.
However, I'm not aware of any way that I can get this information, or indeed any linear combination of coefficients, from a bayesmh command. Does anyone know of a way that this is possible, ideally where not only the linear combination of posterior means is available, but also the credible interval for this combination?
The model I have so far is the following:
Code:
. bayesmh x visit##treatgroup, likelihood(normal({var})) prior({x:}, flat) prior({var}, jeffreys)
Burn-in ...
Simulation ...
Model summary
----------------------------------------------------------------------------------
Likelihood:
x ~ normal(xb_x,{var})
Priors:
{x:2.visit 1.treatgroup 2.visit#1.treatgroup _cons} ~ 1 (flat) (1)
{var} ~ jeffreys
----------------------------------------------------------------------------------
(1) Parameters are elements of the linear form xb_x.
Bayesian normal regression MCMC iterations = 12,500
Random-walk Metropolis-Hastings sampling Burn-in = 2,500
MCMC sample size = 10,000
Number of obs = 200
Acceptance rate = .2615
Efficiency: min = .04421
avg = .05227
Log marginal likelihood = -346.23983 max = .06428
----------------------------------------------------------------------------------
| Equal-tailed
| Mean Std. Dev. MCSE Median [95% Cred. Interval]
-----------------+----------------------------------------------------------------
x |
2.visit | -.2687678 .2686743 .011508 -.2702313 -.7965471 .2631774
1.treatgroup | .1680425 .2703249 .012708 .1667856 -.3614423 .6942513
|
visit#treatgroup |
2 1 | 7.884192 .3770318 .01636 7.884781 7.15229 8.630273
|
_cons | -.2924941 .1875156 .008919 -.2888126 -.6504271 .0803263
-----------------+----------------------------------------------------------------
var | 1.842442 .181575 .007162 1.840625 1.503308 2.2076
----------------------------------------------------------------------------------
Any help would be most appreciated.
Thanks in advance.
Matt
Comment