Pairwise Treatment Comparisons in Bayesian Analysis (bayesmh)

[Mollie Paynee]

Hello,

I am new to Bayesian analysis in Stata and would appreciate some guidance.

I am simulating data with four treatment groups, one control group, and a continuous outcome. I use a Bayesian regression model to estimate the treatment effects. So far, I have used bayestest interval to compute the probability that each treatment effect exceeds a pre-specified cutoff. However, I would like to estimate pairwise probabilities, such as:

• Treatment 1 having a greater effect than Treatment 2
• Treatment 1 having a greater effect than Treatment 3
• ...etc.

I initially thought of specifying something like this:

 bayestest interval ({y:1.r}, upper({y:2.r}))

But this is not allowed. Is there an alternative way to obtain these pairwise comparisons in Stata’s Bayesian framework?

Thank you! Below is my current simplified simulation code:

clear    
set obs 200
    
gen r = floor(5 * runiform())  // Assigns 0-4
label define groups 0 "Control" 1 "Treatment 1" 2 "Treatment 2" 3 "Treatment 3" 4 "Treatment 4"
label values r groups


gen bd = cond(r == 1, -5, -1) + rnormal(0, 1)  // Treatment effects
gen e = rnormal(0, 2)  // Random error
gen y = bd + e  // Outcome variable

//Bayesian Analysis
bayesmh y i.r, ///
likelihood(normal({var})) ///
prior({y:_cons} {y:1.r} {y:2.r} {y:3.r} {y:4.r}, flat) ///
prior({var}, jeffreys) ///
rseed(12345)

bayesstats summary _all
//Which has the highest probability of having a treatment effect
bayestest interval ({y:1.r}, upper(-3)) ({y:2.r}, upper(-3)) ({y:3.r}, upper(-3)) ({y:4.r} , upper(-3))

[Joseph Coveney]

. . . I would like to estimate pairwise probabilities, such as:

• Treatment 1 having a greater effect than Treatment 2
• Treatment 1 having a greater effect than Treatment 3
• ...etc.

I initially thought of specifying something like this:

 bayestest interval ({y:1.r}, upper({y:2.r}))

But this is not allowed. Is there an alternative way to obtain these pairwise comparisons in Stata’s Bayesian framework?

You can specify scalar expressions of model parameters with bayestest interval.

So, you can get your pairwise probabilities for "A greater (lower) than B" kind of questions with syntax along these lines.

bayestest interval ({y:1.r}-{y:2.r}), upper(0)

// and

bayestest interval ({y:1.r}-{y:2.r}), lower(0)

[Mollie Paynee]

You can specify scalar expressions of model parameters with bayestest interval.

So, you can get your pairwise probabilities for "A greater (lower) than B" kind of questions with syntax along these lines.

bayestest interval ({y:1.r}-{y:2.r}), upper(0)

// and

bayestest interval ({y:1.r}-{y:2.r}), lower(0)

Thank you, Joseph! This worked perfectly, thank you so much.