# Using a simulation to show sceptics that t-test gives the right answer, and how to do the same for a bayesian inference case

I want to demonstrate to sceptical colleagues and friends that statistical analysis works. And there’s a follow-on question regarding doing the same with bayesian inference.

Suppose that I want to demonstrate to someone that a two sample t-test at a given level of significance yields the correct results. That is, it correctly tells me that two populations which have no difference between them are shown to not be significantly different, and that two populations which do have a significant difference are successfully detected.

To do this, I’m going to create a simulation containing three populations of people who all suffer from high blood pressure,and we want to use statistics to examine how effective different drugs are in treating it. Population A receives drug A, population B receives drug B, population C (for ‘control’) receives a placebo. We’re going to compare systolic blood pressures in every case. I set up the populations in my simulation so that A has a reduced systolic blood pressure (i.e. the drug is working), whereas B and C have the same elevated level. If the demonstration works successfully, we expect to find, as we take samples from the population and apply the two sample t-test to these sampled data, that there is no significant difference between B and C data, but there is a difference between A and C.

And of course, because this is a statistical process, I need to repeat my test many times, ie. perform many iterations of the same test, to show that the statistical nature of the conclusion is correct.

For the first iteration I select 40 people at random from each of populations A, B and C, compute the mean and standard deviation for each sample set and perform the t-test for the B and C comparison, and also the A and C comparison, i.e. does the value of t exhibit a significant difference or not? In hypothesis language, does the value of t support the null hypothesis or not?

Now I do the same thing for another 99 iterations. What I expect to see in the one hundred A/C comparisons (where we know that there is a reduction in blood pressure in the people taking the drug) is that if I have tested against a t value corresponding to 5% significance, about 19 in every 20 iterations should show a significant difference, i.e. the t result is larger than the value in the table for alpha = 0.05, n = 40.

Similarly, the one hundred results for the B/C comparisons should show fewer than 1 in 20 give a significant difference.

So my first question is, is this a correct premise? Is my proposed method of obtaining the conclusion valid?

Secondly, how would we do a similar thing for a ‘drug trial’ using bayesian inference to find whether drug A or B is better than nothing?

I want to demonstrate to sceptical colleagues and friends that statistical analysis works.

It "works" in the sense that it does what it is designed to do, under the circumstances for which it is designed.

Suppose that I want to demonstrate to someone that a two sample t-test at a given level of significance yields the correct results. That is, it correctly tells me that two populations which have no difference between them are shown to not be significantly different, and that two populations which do have a significant difference are successfully detected.

It definitely can't do that, at least not with the certainty that second sentence seems to imply.

two populations which have no difference between them are shown to not be significantly different

When the null hypothesis is true it should not be rejected at a higher proportion than the significance level, $$\alpha$$. So, for example, if you set your significance level at say 2%, over many such experiments without a difference in the mean (given the remaining assumptions), no more than 2% of them should result in rejection.

The rate of rejection when the null is true is the Type I error rate.

two populations which do have a significant difference are successfully detected

Take care -- this is a mistatement, because 'significant' has a specific meaning. The populations are either different or they are not. "Significance" is a property of hypothesis tests performed on samples, not a property of populations. If you mean a substantial difference or a meaningful difference, you might want to look at what minimum effect size you want to have a good chance to detect.

If you specify such an effect size, then each sample size (or set of sample sizes with multiple groups) will correspond to some probability of rejecting the null (i.e. the power) and some complementary probability of failing to reject it (the type II error rate for that effect size and sample size(s)).

To do this, I’m going to create a simulation containing three populations of people who all suffer from high blood pressure, and we want to use statistics to examine how effective different drugs are in treating it. Population A receives drug A, population B receives drug B, population C (for ‘control’) receives a placebo.

If you have distinct populations and you assign each of them a different drug you have confounded differences in drug effect with differences in population characteristics. To compare drugs what you would likely aim to do is randomly sample some single population and assign the elements of your sample randomly to treatment groups. If you have distinct sub-populations, you'd perhaps want to sample each of them at random and assign each of them randomly to treatment groups (and if at all possible, use their sub-population membership as a covariate)

I'll assume hereafter you mean that A, B and C are distinct treatment groups, rather than distinct target populations.

We’re going to compare systolic blood pressures in every case. I set up the populations in my simulation so that A has a reduced systolic blood pressure (i.e. the drug is working), whereas B and C have the same elevated level. If the demonstration works successfully, we expect to find, as we take samples from the population and apply the two sample t-test to these sampled data, that there is no significant difference between B and C data, but there is a difference between A and C.

If I understand correctly, your intention is that drug B is completely ineffective, while drug A is effective at reducing blood pressure.

And of course, because this is a statistical process, I need to repeat my test many times, ie. perform many iterations of the same test, to show that the statistical nature of the conclusion is correct.

This is simply recognizing that the probability of rejection in the two tests are not 1 and 0 respectively, but something different from those values; as a result we need to sample enough to obtain good estimates of the rejection probability for each (and check that it matches the a priori calculations we get from choosing a significance level and doing a power calculation).

For the first iteration I select 40 people at random from each of populations A, B and C,

Why 40? How do you know that will be enough to detect whatever the effect size was that you chose?

Now I do the same thing for another 99 iterations.

Why only 100 such sets of samples? You could do a hundred thousand in a matter of seconds, or a million in seconds if you do it efficiently; this will reduce the simulation error substantially. The most time consuming part is setting it up, even though the calculations are only a few lines of code (e.g. in R you can make use of replicate).

What I expect to see in the one hundred A/C comparisons (where we know that there is a reduction in blood pressure in the people taking the drug) is that if I have tested against a t value corresponding to 5% significance, about 19 in every 20 iterations should show a significant difference, i.e. the t result is larger than the value in the table for alpha = 0.05, n = 40.

Why do you think the power will be 95% in this instance? I see nothing to indicate that this will be the case.

Similarly, the one hundred results for the B/C comparisons should show fewer than 1 in 20 give a significant difference.

No, the number of rejections in 100 trials will be binomial(100,0.05); the chance you'll observe fewer than 5 rejections (i.e. fewer than 1 in 20) is only about 44%.

So my first question is, is this a correct premise? Is my proposed method of obtaining the conclusion valid?

No, as noted above. You'll need to investigate the power for your effect size and sample size, and then show that in the simulation you attain the advertized power. You will need a much larger simulation size (repeats of the experiment) so that your results are not driven by randomness (i.e. so that you have accurate estimates of the type I and type II rejection rates)

Secondly, how would we do a similar thing for a ‘drug trial’ using bayesian inference to find whether drug A or B is better than nothing?

You can't expect Bayesian calculations to have a set of frequentist properties that they don't attempt to achieve.*

The point of a Bayesian calculation is slightly different.

It depends on what the precise Bayesian calculations are that you intend to perform. For example, if you're calculating a posterior distribution on the effect of the drug you could attempt to demonstrate that the probability statements associated with the posterior for that parameter would be approximately correct. (One difficulty is that the posteriors depend on the priors.)

Or if you're invoking Bayesian decision theory then you can simulate to show that the resulting decisions have the properties they're supposed to have.

It's a little difficult to say much without a clearly expressed procedure. Your discussion of the t-test was very specific so it was easy to discuss.

*(Of course, in straightforward cases, Bayesian calculations often do attain useful frequentist properties but to expect them to do so as a matter of course would be to miss the point)

• Hi @Glen_b, and thank you for that response. If I were instead to say that the expected number of the B/C comparisons showing a significant difference is 1 in 20, is that closer to the mark, with less sloppy terminology? And also, am I now saying the same thing you said in the discussion of the xkcd jelly bean strip regarding expected number of significant results when nothing is really going on? stats.stackexchange.com/questions/88065/… – D R Ball Dec 30 '18 at 22:18
• 1. Yes, with a t-test (given the assumptions) you should attain the nominal type I error rate, so the expected proportion would be 1/20. The issue was two things going on at once -- in the long run you should have no more than a proportion $α$ of rejections, but in 100 such tests the proportion of rejections may exceed its expectation. $\:$ 2. Yes, the green jelly bean thing does make that point about multiple testing but its underlying points are also about problems with reporting of science in the media (and perhaps implicitly, publication bias). – Glen_b -Reinstate Monica Dec 31 '18 at 2:25
• Thanks @Glen_b. So we're saying that the simulation should reflect the type I error rate correctly. Can we make a corresponding statement for the type II error rate? In the A/C case, a type II error corresponds to the calculations showing the null hypothesis to be true. I stated in my OP that I thought this would also occur 5% of the time, but if I read your response correctly, you say this is wrong. What is the subtlety I'm missing please? Why is the probability of type II error not the same as the significance level I'm using. Is it just the way? – D R Ball Dec 31 '18 at 16:27
• The calculations would lead you not to reject the null; they can't show the null to be true. You're correct just above - the type II error is not equal to the significance level in general. You seem to feel it should be the same as the type I error rate but that's almost never the case, unless by design. By way of illustration, let me pose this hypothetical: imagine the treatment applied to group A is effective and so groups A and C differ in their average blood pressure. ... ctd – Glen_b -Reinstate Monica Jan 1 '19 at 1:27
• ctd ... Imagine you and I both conduct similar experiments. You use some sample size $m$ for both the treatment and control group and imagine it does turn out that the power is 95%. Now I run the experiment independently of you but my sample sizes are each $n$, where $n$ is much smaller than $m$ (e.g. $n=5$ and $m=500$ say). You estimated both means much more precisely than I did. Ponder what sort of confidence interval we might each obtain for the difference in means. Should my rejection rate be the same as yours? – Glen_b -Reinstate Monica Jan 1 '19 at 1:27