One-sample/paired t-tests
One-sample t-tests are already covered in the basic
usage vignette, so the example presented here is simply a repeat of
that example. However, we’ll change the order of presentation around a
bit so that the relationship between the one-sample and independent
samples models are a little clearer.
We’ll start with an example from Rouder et al.
(2009), in which they analyse the results if a one sample
t-test. Rouder et al. (2009)
report a t statistic of 2.03, from a sample size of 80. The
sampling distribution of the t statistic, when the null
hypothesis is false is the noncentral
t-distribution distribution. Therefore, we can use this
fact to construct a likelihood function for our t
value from the noncentral t-distribution. For this, we’ll need
two parameters—the t value itself, and the degrees of
freedom. In the one-sample case, the degrees of freedom
will just be N - 1.
t <- 2.03
n <- 80
data_model <- likelihood("noncentral_t", t = t, df = n - 1)
plot(data_model)
From plotting our likelihood function we can see that values
of t at our observation (t = 2.03) and most consistent
without our observation while more extreme values of t (e.g.,
t = -2, or t = 6) are less consistent without our
observed data. This is of course, as we would expect.
Now that we’ve defined our likelihood, the next thing to think about
is the prior. Rouder et al.
(2009) provide a more detailed justification for their prior, so
we won’t go into that here; however, their prior is based on the
Cauchy distribution. A Cauchy distribution is like a very fat
tailed t distribution (in fact, it is a very
fat tailed t distribution). We can see an example of a
standard Cauchy distribution below.
plot(prior("cauchy", location = 0, scale = 1))
The Cauchy distribution can be scaled so that it is
wider or narrow. You might want to choose a scale parameter
based on the range of t values that your theory prediction. For
example, here is a slightly narrower Cauchy
distribution.
plot(prior("cauchy", location = 0, scale = .707))
Although scaling the prior distribution by our theoretical
predictions makes sense, we also need to factor in another
thing—sample size. Because of how the t statistic is calculated
(\(t = \frac{\mu}{\sigma/\sqrt{n}}\)),
we can see that t statistics are also dependent on sample size.
That is, for any given effect size (standardized mean
difference), the corresponding t value will be a function
of sample size—or more, specifically, the square root of the sample
size. That is, for a given underlying true effect size, when we
have a large sample size we can expect to see larger values of
t that if we have a small sample size. Therefore, we should
also scale our prior by the square root of the sample size. In the
specification below, we’re specifying a standard Cauchy, but
further scaling this by \(\sqrt{n}\).
alt_prior <- prior("cauchy", location = 0, scale = 1 * sqrt(n))
plot(alt_prior)
For our null we’ll use a point located at 0.
null_prior <- prior("point", point = 0)
plot(null_prior)
Now we can proceed to compute the Bayes factor in the standard
way.
bf_onesample_1 <- integral(data_model * alt_prior) / integral(data_model * null_prior)
summary(bf_onesample_1)
#> Bayes factor
#> Using the levels from Wagenmakers et al (2017)
#> A BF of 0.6421 indicates:
#> Anecdotal evidence
In the preceding example, we’ve taken account of the fact that the
same underlying effect size will produce different values of
t depending on the sample, and we’ve scaled our prior
accordingly. However, we could also apply the scaling at the other end
by re-scaling our likelihood. The bayesplay
package contains two additional noncentral t likelihoods that
have been rescaled. The first of the these is the
noncentral_d likelihood. This is a likelihood based on the
sample distribution of the one-sample/paired samples Cohen’s d.
This is calculated simply as \(d =
\frac{\mu_{\mathrm{diff}}}{\sigma_{\mathrm{diff}}}\).
Alternatively, we can conver the observed t value by dividing
it by \(\sqrt{n}\). The
noncentral_d likelihood just takes this effect size and the
sample size are parameters.
d <- t / sqrt(n)
data_model2 <- likelihood("noncentral_d", d = d, n = n)
plot(data_model2)
Now that we’ve applied our scaling to the likelihood, we
don’t need to apply this scaling to the prior. Therefore, if we
wanted to use the same prior as the preceding example, we’ll have a
scale value of 1, rather than
1 * sqrt(n).
alt_prior2 <- prior("cauchy", location = 0, scale = 1)
plot(alt_prior2)
We can re-use our null prior from before and calculate the Bayes
factor the same way as before.
bf_onesample_2 <- integral(data_model2 * alt_prior2) / integral(data_model2 * null_prior)
summary(bf_onesample_2)
#> Bayes factor
#> Using the levels from Wagenmakers et al (2017)
#> A BF of 0.6421 indicates:
#> Anecdotal evidence
As expected, the two results are identical.
Independent samples t-tests
Rouder et al. (2009) also provide an
extension of their method to the two sample case, although they do not
provide a worked example. Instead, we can generate our own example and
directly compare the results from the bayesplay package
with the results from the BayesFactor
package.
For this example, we’ll start by generating some data from an
independent samples design.
set.seed(2125519)
group1 <- 25 + scale(rnorm(n = 15)) * 15
group2 <- 35 + scale(rnorm(n = 16)) * 16
First, let us see the results from the BayesFactor
package.
BayesFactor::ttestBF(x = group1, y = group2, paired = FALSE, rscale = 1)
Bayes factor analysis
--------------
[1] Alt., r=1 : 0.9709424 ±0%
Against denominator:
Null, mu1-mu2 = 0
---
Bayes factor type: BFindepSample, JZS
As with the one-sample case, we can run the analysis in the
bayesplay package using either the t statistic or
the Cohen’s d. We’ll start by running the analysis using the
t statistic. The easiest way to do this, is to simply use the
t.test() function in R.
t_result <- t.test(x = group1, y = group2, paired = FALSE, var.equal = TRUE)
t_result
#>
#> Two Sample t-test
#>
#> data: group1 and group2
#> t = -1.7922, df = 29, p-value = 0.08354
#> alternative hypothesis: true difference in means is not equal to 0
#> 95 percent confidence interval:
#> -21.411872 1.411872
#> sample estimates:
#> mean of x mean of y
#> 25 35
From this output we need the t statistic itself, and the
degrees of freedom.
t <- t_result$statistic
t
#> t
#> -1.792195
df <- t_result$parameter
df
#> df
#> 29
With the t, and df value in hand, we can specify
our likelihood using the same noncentral t distribution as the
one sample case.
data_model3 <- likelihood("noncentral_t", t = t, df = df)
As with the one-sample case, a Cauchy prior is used for the
alternative hypothesis. Again, this will need to be appropriately
scaled. In one sample case we scaled it by \(\sqrt{n}\). In the two-sample case,
however, we’ll scale it by \(\sqrt{\frac{n_1
\times n_2}{n_1 + n_2}}\)
n1 <- length(group1)
n2 <- length(group2)
alt_prior3 <- prior("cauchy", location = 0, scale = 1 * sqrt((n1 * n2) / (n1 + n2)))
plot(alt_prior3)
We’ll use the same point null prior as before, and then compute the
Bayes factor in the usual way.
bf_independent_1 <- integral(data_model3 * alt_prior3) / integral(data_model3 * null_prior)
summary(bf_independent_1)
#> Bayes factor
#> Using the levels from Wagenmakers et al (2017)
#> A BF of 0.9709 indicates:
#> Anecdotal evidence
Appropriately scaling our Cauchy prior can be tricky, so an
alternative is instead, as before, to scale our likelihood. The
bayesplay package contains a likelihood that is appropriate
for independent samples Cohen’s d, the
noncentral_d2 likelihood. To use this, we’ll need the
Cohen’s d value, and the two sample sizes.
For independent samples designs the Cohen’s d is calculated
as follows:
\[d = \frac{m_1 -
m2}{s_\mathrm{pooled}},\]
where \(s_\mathrm{pooled}\) is given
as follows:
\[s_\mathrm{pooled} = \sqrt{\frac{(n_1 -
1)s^2_1 + (n_2 - 1)s^2_2}{n_1 + n_2 -2}}\]
This is fairly straightforward to calculate, as shown below.
m1 <- mean(group1)
m2 <- mean(group2)
s1 <- sd(group1)
s2 <- sd(group2)
md_diff <- m1 - m2
sd_pooled <- sqrt((((n1 - 1) * s1^2) + ((n2 - 1) * s2^2)) / (n1 + n2 - 2))
d <- md_diff / sd_pooled
d
#> [1] -0.6441105
However, it can also be obtained from the effsize
package using the following syntax.
effsize::cohen.d(group1, group2, paired = FALSE, hedged.correction = FALSE)
With the d value in hand, we can how specify a new
likelihood.
data_model4 <- likelihood("noncentral_d2", d = d, n1 = n1, n2 = n2)
data_model4
#> Likelihood
#> Family
#> noncentral_d2
#> Parameters
#> d: -0.644110547740848
#> n1: 15
#> n2: 16
#>
Because we’ve used the appropriately scaled noncentral t
likelihood, the noncentral_d2, we no longer need to scale
the Cauchy prior.
alt_prior4 <- prior("cauchy", location = 0, scale = 1)
plot(alt_prior4)
And we can now calculate the Bayes factor in the usual way.
bf_independent_2 <- integral(data_model4 * alt_prior4) / integral(data_model4 * null_prior)
summary(bf_independent_2)
#> Bayes factor
#> Using the levels from Wagenmakers et al (2017)
#> A BF of 0.9709 indicates:
#> Anecdotal evidence
Again, we obtain the same result as the BayesFactor
package.