ℹ️ Skipped - page is already crawled
| Filter | Status | Condition | Details |
|---|---|---|---|
| HTTP status | PASS | download_http_code = 200 | HTTP 200 |
| Age cutoff | PASS | download_stamp > now() - 6 MONTH | 0.2 months ago |
| History drop | PASS | isNull(history_drop_reason) | No drop reason |
| Spam/ban | PASS | fh_dont_index != 1 AND ml_spam_score = 0 | ml_spam_score=0 |
| Canonical | PASS | meta_canonical IS NULL OR = '' OR = src_unparsed | Not set |
| Property | Value |
|---|---|
| URL | https://learningstatisticswithr.com/book/bayes.html |
| Last Crawled | 2026-04-13 09:57:09 (5 days ago) |
| First Indexed | 2019-01-08 19:25:01 (7 years ago) |
| HTTP Status Code | 200 |
| Meta Title | Chapter 17 Bayesian statistics | Learning statistics with R: A tutorial for psychology students and other beginners. (Version 0.6.1) |
| Meta Description | Learning Statistics with R covers the contents of an introductory statistics class, as typically taught to undergraduate psychology students, focusing on the use of the R statistical software. |
| Meta Canonical | null |
| Boilerpipe Text | In our reasonings concerning matter of fact, there are all imaginable degrees of assurance, from the highest certainty to the lowest species of moral evidence. A wise man, therefore, proportions his belief to the evidence.
– David Hume
254
.
The ideas I’ve presented to you in this book describe inferential statistics from the frequentist perspective. I’m not alone in doing this. In fact, almost every textbook given to undergraduate psychology students presents the opinions of the frequentist statistician as
the
theory of inferential statistics, the one true way to do things. I have taught this way for practical reasons. The frequentist view of statistics dominated the academic field of statistics for most of the 20th century, and this dominance is even more extreme among applied scientists. It was and is current practice among psychologists to use frequentist methods. Because frequentist methods are ubiquitous in scientific papers, every student of statistics needs to understand those methods, otherwise they will be unable to make sense of what those papers are saying! Unfortunately – in my opinion at least – the current practice in psychology is often misguided, and the reliance on frequentist methods is partly to blame. In this chapter I explain why I think this, and provide an introduction to Bayesian statistics, an approach that I think is generally superior to the orthodox approach.
This chapter comes in two parts. In Sections
17.1
through
17.3
I talk about what Bayesian statistics are all about, covering the basic mathematical rules for how it works as well as an explanation for why I think the Bayesian approach is so useful. Afterwards, I provide a brief overview of how you can do Bayesian versions of chi-square tests (Section
17.6
),
\(t\)
-tests (Section
17.7
), regression (Section
17.8
) and ANOVA (Section
17.9
).
Probabilistic reasoning by rational agents
From a Bayesian perspective, statistical inference is all about
belief revision
. I start out with a set of candidate hypotheses
\(h\)
about the world. I don’t know which of these hypotheses is true, but do I have some beliefs about which hypotheses are plausible and which are not. When I observe the data
\(d\)
, I have to revise those beliefs. If the data are consistent with a hypothesis, my belief in that hypothesis is strengthened. If the data inconsistent with the hypothesis, my belief in that hypothesis is weakened. That’s it! At the end of this section I’ll give a precise description of how Bayesian reasoning works, but first I want to work through a simple example in order to introduce the key ideas. Consider the following reasoning problem:
I’m carrying an umbrella. Do you think it will rain?
In this problem, I have presented you with a single piece of data (
\(d =\)
I’m carrying the umbrella), and I’m asking you to tell me your beliefs about whether it’s raining. You have two possible
hypotheses
,
\(h\)
: either it rains today or it does not. How should you solve this problem?
Priors: what you believed before
The first thing you need to do ignore what I told you about the umbrella, and write down your pre-existing beliefs about rain. This is important: if you want to be honest about how your beliefs have been revised in the light of new evidence, then you
must
say something about what you believed before those data appeared! So, what might you believe about whether it will rain today? You probably know that I live in Australia, and that much of Australia is hot and dry. And in fact you’re right: the city of Adelaide where I live has a Mediterranean climate, very similar to southern California, southern Europe or northern Africa. I’m writing this in January, and so you can assume it’s the middle of summer. In fact, you might have decided to take a quick look on Wikipedia
255
and discovered that Adelaide gets an average of 4.4 days of rain across the 31 days of January. Without knowing anything else, you might conclude that the probability of January rain in Adelaide is about 15%, and the probability of a dry day is 85%. If this is really what you believe about Adelaide rainfall (and now that I’ve told it to you, I’m betting that this really
is
what you believe) then what I have written here is your
prior distribution
, written
\(P(h)\)
:
Hypothesis
Degree of Belief
Rainy day
0.15
Dry day
0.85
Likelihoods: theories about the data
To solve the reasoning problem, you need a theory about my behaviour. When does Dan carry an umbrella? You might guess that I’m not a complete idiot,
256
and I try to carry umbrellas only on rainy days. On the other hand, you also know that I have young kids, and you wouldn’t be all that surprised to know that I’m pretty forgetful about this sort of thing. Let’s suppose that on rainy days I remember my umbrella about 30% of the time (I really am awful at this). But let’s say that on dry days I’m only about 5% likely to be carrying an umbrella. So you might write out a little table like this:
Hypothesis
Umbrella
No umbrella
Rainy day
0.30
0.70
Dry day
0.05
0.95
It’s important to remember that each cell in this table describes your beliefs about what data
\(d\)
will be observed,
given
the truth of a particular hypothesis
\(h\)
. This “conditional probability” is written
\(P(d|h)\)
, which you can read as “the probability of
\(d\)
given
\(h\)
”. In Bayesian statistics, this is referred to as
likelihood
of data
\(d\)
given hypothesis
\(h\)
.
257
The joint probability of data and hypothesis
At this point, all the elements are in place. Having written down the priors and the likelihood, you have all the information you need to do Bayesian reasoning. The question now becomes,
how
do we use this information? As it turns out, there’s a very simple equation that we can use here, but it’s important that you understand why we use it, so I’m going to try to build it up from more basic ideas.
Let’s start out with one of the rules of probability theory. I listed it way back in Table
9.1
, but I didn’t make a big deal out of it at the time and you probably ignored it. The rule in question is the one that talks about the probability that
two
things are true. In our example, you might want to calculate the probability that today is rainy (i.e., hypothesis
\(h\)
is true)
and
I’m carrying an umbrella (i.e., data
\(d\)
is observed). The
joint probability
of the hypothesis and the data is written
\(P(d,h)\)
, and you can calculate it by multiplying the prior
\(P(h)\)
by the likelihood
\(P(d|h)\)
. Mathematically, we say that:
\[
P(d,h) = P(d|h) P(h)
\]
So, what is the probability that today is a rainy day
and
I remember to carry an umbrella? As we discussed earlier, the prior tells us that the probability of a rainy day is 15%, and the likelihood tells us that the probability of me remembering my umbrella on a rainy day is 30%. So the probability that both of these things are true is calculated by multiplying the two:
\[
\begin{array}
P(\mbox{rainy}, \mbox{umbrella}) & = & P(\mbox{umbrella} | \mbox{rainy}) \times P(\mbox{rainy}) \\
& = & 0.30 \times 0.15 \\
& = & 0.045
\end{array}
\]
In other words, before being told anything about what actually happened, you think that there is a 4.5% probability that today will be a rainy day and that I will remember an umbrella. However, there are of course
four
possible things that could happen, right? So let’s repeat the exercise for all four. If we do that, we end up with the following table:
Umbrella
No-umbrella
Rainy
0.045
0.105
Dry
0.0425
0.8075
This table captures all the information about which of the four possibilities are likely. To really get the full picture, though, it helps to add the row totals and column totals. That gives us this table:
Umbrella
No-umbrella
Total
Rainy
0.0450
0.1050
0.15
Dry
0.0425
0.8075
0.85
Total
0.0875
0.9125
1
This is a very useful table, so it’s worth taking a moment to think about what all these numbers are telling us. First, notice that the row sums aren’t telling us anything new at all. For example, the first row tells us that if we ignore all this umbrella business, the chance that today will be a rainy day is 15%. That’s not surprising, of course: that’s our prior. The important thing isn’t the number itself: rather, the important thing is that it gives us some confidence that our calculations are sensible! Now take a look at the column sums, and notice that they tell us something that we haven’t explicitly stated yet. In the same way that the row sums tell us the probability of rain, the column sums tell us the probability of me carrying an umbrella. Specifically, the first column tells us that on average (i.e., ignoring whether it’s a rainy day or not), the probability of me carrying an umbrella is 8.75%. Finally, notice that when we sum across all four logically-possible events, everything adds up to 1. In other words, what we have written down is a proper probability distribution defined over all possible combinations of data and hypothesis.
Now, because this table is so useful, I want to make sure you understand what all the elements correspond to, and how they written:
Umbrella
No-umbrella
Rainy
\(P\)
(Umbrella, Rainy)
\(P\)
(No-umbrella, Rainy)
\(P\)
(Rainy)
Dry
\(P\)
(Umbrella, Dry)
\(P\)
(No-umbrella, Dry)
\(P\)
(Dry)
\(P\)
(Umbrella)
\(P\)
(No-umbrella)
Finally, let’s use “proper” statistical notation. In the rainy day problem, the data corresponds to the observation that I do or do not have an umbrella. So we’ll let
\(d_1\)
refer to the possibility that you observe me carrying an umbrella, and
\(d_2\)
refers to you observing me not carrying one. Similarly,
\(h_1\)
is your hypothesis that today is rainy, and
\(h_2\)
is the hypothesis that it is not. Using this notation, the table looks like this:
\(d_1\)
\(d_2\)
\(h_1\)
\(P(h_1, d_1)\)
\(P(h_1, d_2)\)
\(P(h_1)\)
\(h_2\)
\(P(h_2, d_1)\)
\(P(h_2, d_2)\)
\(P(h_2)\)
\(P(d_1)\)
\(P(d_2)\)
Updating beliefs using Bayes’ rule
The table we laid out in the last section is a very powerful tool for solving the rainy day problem, because it considers all four logical possibilities and states exactly how confident you are in each of them before being given any data. It’s now time to consider what happens to our beliefs when we are actually given the data. In the rainy day problem, you are told that I really
am
carrying an umbrella. This is something of a surprising event: according to our table, the probability of me carrying an umbrella is only 8.75%. But that makes sense, right? A guy carrying an umbrella on a summer day in a hot dry city is pretty unusual, and so you really weren’t expecting that. Nevertheless, the problem tells you that it is true. No matter how unlikely you thought it was, you must now adjust your beliefs to accommodate the fact that you now
know
that I have an umbrella.
258
To reflect this new knowledge, our
revised
table must have the following numbers:
Umbrella
No-umbrella
Rainy
0
Dry
0
Total
1
0
In other words, the facts have eliminated any possibility of “no umbrella”, so we have to put zeros into any cell in the table that implies that I’m not carrying an umbrella. Also, you know for a fact that I am carrying an umbrella, so the column sum on the left must be 1 to correctly describe the fact that
\(P(\mbox{umbrella})=1\)
.
What two numbers should we put in the empty cells? Again, let’s not worry about the maths, and instead think about our intuitions. When we wrote out our table the first time, it turned out that those two cells had almost identical numbers, right? We worked out that the joint probability of “rain and umbrella” was 4.5%, and the joint probability of “dry and umbrella” was 4.25%. In other words, before I told you that I am in fact carrying an umbrella, you’d have said that these two events were almost identical in probability, yes? But notice that
both
of these possibilities are consistent with the fact that I actually am carrying an umbrella. From the perspective of these two possibilities, very little has changed. I hope you’d agree that it’s
still
true that these two possibilities are equally plausible. So what we expect to see in our final table is some numbers that preserve the fact that “rain and umbrella” is
slightly
more plausible than “dry and umbrella”, while still ensuring that numbers in the table add up. Something like this, perhaps?
Umbrella
No-umbrella
Rainy
0.514
0
Dry
0.486
0
Total
1
0
What this table is telling you is that, after being told that I’m carrying an umbrella, you believe that there’s a 51.4% chance that today will be a rainy day, and a 48.6% chance that it won’t. That’s the answer to our problem! The
posterior probability
of rain
\(P(h|d)\)
given that I am carrying an umbrella is 51.4%
How did I calculate these numbers? You can probably guess. To work out that there was a 0.514 probability of “rain”, all I did was take the 0.045 probability of “rain and umbrella” and divide it by the 0.0875 chance of “umbrella”. This produces a table that satisfies our need to have everything sum to 1, and our need not to interfere with the relative plausibility of the two events that are actually consistent with the data. To say the same thing using fancy statistical jargon, what I’ve done here is divide the joint probability of the hypothesis and the data
\(P(d,h)\)
by the
marginal probability
of the data
\(P(d)\)
, and this is what gives us the posterior probability of the hypothesis
given
that we know the data have been observed. To write this as an equation:
259
\[
P(h | d) = \frac{P(d,h)}{P(d)}
\]
However, remember what I said at the start of the last section, namely that the joint probability
\(P(d,h)\)
is calculated by multiplying the prior
\(P(h)\)
by the likelihood
\(P(d|h)\)
. In real life, the things we actually know how to write down are the priors and the likelihood, so let’s substitute those back into the equation. This gives us the following formula for the posterior probability:
\[
P(h | d) = \frac{P(d|h) P(h)}{P(d)}
\]
And this formula, folks, is known as
Bayes’ rule
. It describes how a learner starts out with prior beliefs about the plausibility of different hypotheses, and tells you how those beliefs should be revised in the face of data. In the Bayesian paradigm, all statistical inference flows from this one simple rule.
Bayesian hypothesis tests
In Chapter
11
I described the orthodox approach to hypothesis testing. It took an entire chapter to describe, because null hypothesis testing is a very elaborate contraption that people find very hard to make sense of. In contrast, the Bayesian approach to hypothesis testing is incredibly simple. Let’s pick a setting that is closely analogous to the orthodox scenario. There are two hypotheses that we want to compare, a null hypothesis
\(h_0\)
and an alternative hypothesis
\(h_1\)
. Prior to running the experiment we have some beliefs
\(P(h)\)
about which hypotheses are true. We run an experiment and obtain data
\(d\)
. Unlike frequentist statistics Bayesian statistics does allow to talk about the probability that the null hypothesis is true. Better yet, it allows us to calculate the
posterior probability of the null hypothesis
, using Bayes’ rule:
\[
P(h_0 | d) = \frac{P(d|h_0) P(h_0)}{P(d)}
\]
This formula tells us exactly how much belief we should have in the null hypothesis after having observed the data
\(d\)
. Similarly, we can work out how much belief to place in the alternative hypothesis using essentially the same equation. All we do is change the subscript:
\[
P(h_1 | d) = \frac{P(d|h_1) P(h_1)}{P(d)}
\]
It’s all so simple that I feel like an idiot even bothering to write these equations down, since all I’m doing is copying Bayes rule from the previous section.
260
The Bayes factor
In practice, most Bayesian data analysts tend not to talk in terms of the raw posterior probabilities
\(P(h_0|d)\)
and
\(P(h_1|d)\)
. Instead, we tend to talk in terms of the
posterior odds
ratio. Think of it like betting. Suppose, for instance, the posterior probability of the null hypothesis is 25%, and the posterior probability of the alternative is 75%. The alternative hypothesis is three times as probable as the null, so we say that the
odds
are 3:1 in favour of the alternative. Mathematically, all we have to do to calculate the posterior odds is divide one posterior probability by the other:
\[
\frac{P(h_1 | d)}{P(h_0 | d)} = \frac{0.75}{0.25} = 3
\]
Or, to write the same thing in terms of the equations above:
\[
\frac{P(h_1 | d)}{P(h_0 | d)} = \frac{P(d|h_1)}{P(d|h_0)} \times \frac{P(h_1)}{P(h_0)}
\]
Actually, this equation is worth expanding on. There are three different terms here that you should know. On the left hand side, we have the posterior odds, which tells you what you believe about the relative plausibilty of the null hypothesis and the alternative hypothesis
after
seeing the data. On the right hand side, we have the
prior odds
, which indicates what you thought
before
seeing the data. In the middle, we have the
Bayes factor
, which describes the amount of evidence provided by the data:
\[
\begin{array}{ccccc}\displaystyle
\frac{P(h_1 | d)}{P(h_0 | d)} &=& \displaystyle\frac{P(d|h_1)}{P(d|h_0)} &\times& \displaystyle\frac{P(h_1)}{P(h_0)} \\[6pt] \\[-2pt]
\uparrow && \uparrow && \uparrow \\[6pt]
\mbox{Posterior odds} && \mbox{Bayes factor} && \mbox{Prior odds}
\end{array}
\]
The Bayes factor (sometimes abbreviated as
BF
) has a special place in the Bayesian hypothesis testing, because it serves a similar role to the
\(p\)
-value in orthodox hypothesis testing: it quantifies the strength of evidence provided by the data, and as such it is the Bayes factor that people tend to report when running a Bayesian hypothesis test. The reason for reporting Bayes factors rather than posterior odds is that different researchers will have different priors. Some people might have a strong bias to believe the null hypothesis is true, others might have a strong bias to believe it is false. Because of this, the polite thing for an applied researcher to do is report the Bayes factor. That way, anyone reading the paper can multiply the Bayes factor by their own
personal
prior odds, and they can work out for themselves what the posterior odds would be. In any case, by convention we like to pretend that we give equal consideration to both the null hypothesis and the alternative, in which case the prior odds equals 1, and the posterior odds becomes the same as the Bayes factor.
Interpreting Bayes factors
One of the really nice things about the Bayes factor is the numbers are inherently meaningful. If you run an experiment and you compute a Bayes factor of 4, it means that the evidence provided by your data corresponds to betting odds of 4:1 in favour of the alternative. However, there have been some attempts to quantify the standards of evidence that would be considered meaningful in a scientific context. The two most widely used are from
Jeffreys (
1961
)
and
Kass and Raftery (
1995
)
. Of the two, I tend to prefer the
Kass and Raftery (
1995
)
table because it’s a bit more conservative. So here it is:
Bayes factor
Interpretation
1 - 3
Negligible evidence
3 - 20
Positive evidence
20 - 150
Strong evidence
$>$150
Very strong evidence
And to be perfectly honest, I think that even the Kass and Raftery standards are being a bit charitable. If it were up to me, I’d have called the “positive evidence” category “weak evidence”. To me, anything in the range 3:1 to 20:1 is “weak” or “modest” evidence at best. But there are no hard and fast rules here: what counts as strong or weak evidence depends entirely on how conservative you are, and upon the standards that your community insists upon before it is willing to label a finding as “true”.
In any case, note that all the numbers listed above make sense if the Bayes factor is greater than 1 (i.e., the evidence favours the alternative hypothesis). However, one big practical advantage of the Bayesian approach relative to the orthodox approach is that it also allows you to quantify evidence
for
the null. When that happens, the Bayes factor will be less than 1. You can choose to report a Bayes factor less than 1, but to be honest I find it confusing. For example, suppose that the likelihood of the data under the null hypothesis
\(P(d|h_0)\)
is equal to 0.2, and the corresponding likelihood
\(P(d|h_0)\)
under the alternative hypothesis is 0.1. Using the equations given above, Bayes factor here would be:
\[
\mbox{BF} = \frac{P(d|h_1)}{P(d|h_0)} = \frac{0.1}{0.2} = 0.5
\]
Read literally, this result tells is that the evidence in favour of the alternative is 0.5 to 1. I find this hard to understand. To me, it makes a lot more sense to turn the equation “upside down”, and report the amount op evidence in favour of the
null
. In other words, what we calculate is this:
\[
\mbox{BF}^\prime = \frac{P(d|h_0)}{P(d|h_1)} = \frac{0.2}{0.1} = 2
\]
And what we would report is a Bayes factor of 2:1 in favour of the null. Much easier to understand, and you can interpret this using the table above.
Why be a Bayesian?
Up to this point I’ve focused exclusively on the logic underpinning Bayesian statistics. We’ve talked about the idea of “probability as a degree of belief”, and what it implies about how a rational agent should reason about the world. The question that you have to answer for yourself is this: how do
you
want to do your statistics? Do you want to be an orthodox statistician, relying on sampling distributions and
\(p\)
-values to guide your decisions? Or do you want to be a Bayesian, relying on Bayes factors and the rules for rational belief revision? And to be perfectly honest, I can’t answer this question for you. Ultimately it depends on what you think is right. It’s your call, and your call alone. That being said, I can talk a little about why
I
prefer the Bayesian approach.
Statistics that mean what you think they mean
You keep using that word. I do not think it means what you think it means
– Inigo Montoya, The Princess Bride
261
To me, one of the biggest advantages to the Bayesian approach is that it answers the right questions. Within the Bayesian framework, it is perfectly sensible and allowable to refer to “the probability that a hypothesis is true”. You can even try to calculate this probability. Ultimately, isn’t that what you
want
your statistical tests to tell you? To an actual human being, this would seem to be the whole
point
of doing statistics: to determine what is true and what isn’t. Any time that you aren’t exactly sure about what the truth is, you should use the language of probability theory to say things like “there is an 80% chance that Theory A is true, but a 20% chance that Theory B is true instead”.
This seems so obvious to a human, yet it is explicitly forbidden within the orthodox framework. To a frequentist, such statements are a nonsense because “the theory is true” is not a repeatable event. A theory is true or it is not, and no probabilistic statements are allowed, no matter how much you might want to make them. There’s a reason why, back in Section
11.5
, I repeatedly warned you
not
to interpret the
\(p\)
-value as the probability of that the null hypothesis is true. There’s a reason why almost every textbook on statstics is forced to repeat that warning. It’s because people desperately
want
that to be the correct interpretation. Frequentist dogma notwithstanding, a lifetime of experience of teaching undergraduates and of doing data analysis on a daily basis suggests to me that most actual humans thing that “the probability that the hypothesis is true” is not only meaningful, it’s the thing we care
most
about. It’s such an appealing idea that even trained statisticians fall prey to the mistake of trying to interpret a
\(p\)
-value this way. For example, here is a quote from an official Newspoll report in 2013, explaining how to interpret their (frequentist) data analysis:
262
Throughout the report, where relevant, statistically significant changes have been noted. All significance tests have been based on the 95 percent level of confidence.
This means that if a change is noted as being statistically significant, there is a 95 percent probability that a real change has occurred
, and is not simply due to chance variation. (emphasis added)
Nope! That’s
not
what
\(p<.05\)
means. That’s
not
what 95% confidence means to a frequentist statistician. The bolded section is just plain wrong. Orthodox methods cannot tell you that “there is a 95% chance that a real change has occurred”, because this is not the kind of event to which frequentist probabilities may be assigned. To an ideological frequentist, this sentence should be meaningless. Even if you’re a more pragmatic frequentist, it’s still the wrong definition of a
\(p\)
-value. It is simply not an allowed or correct thing to say if you want to rely on orthodox statistical tools.
On the other hand, let’s suppose you are a Bayesian. Although the bolded passage is the wrong definition of a
\(p\)
-value, it’s pretty much exactly what a Bayesian means when they say that the posterior probability of the alternative hypothesis is greater than 95%. And here’s the thing. If the Bayesian posterior is actually thing you
want
to report, why are you even trying to use orthodox methods? If you want to make Bayesian claims, all you have to do is be a Bayesian and use Bayesian tools.
Speaking for myself, I found this to be a the most liberating thing about switching to the Bayesian view. Once you’ve made the jump, you no longer have to wrap your head around counterinuitive definitions of
\(p\)
-values. You don’t have to bother remembering why you can’t say that you’re 95% confident that the true mean lies within some interval. All you have to do is be honest about what you believed before you ran the study, and then report what you learned from doing it. Sounds nice, doesn’t it? To me, this is the big promise of the Bayesian approach: you do the analysis you really want to do, and express what you really believe the data are telling you.
Evidentiary standards you can believe
If [
\(p\)
] is below .02 it is strongly indicated that the [null] hypothesis fails to account for the whole of the facts. We shall not often be astray if we draw a conventional line at .05 and consider that [smaller values of
\(p\)
] indicate a real discrepancy.
– Sir Ronald
Fisher (
1925
)
Consider the quote above by Sir Ronald Fisher, one of the founders of what has become the orthodox approach to statistics. If anyone has ever been entitled to express an opinion about the intended function of
\(p\)
-values, it’s Fisher. In this passage, taken from his classic guide
Statistical Methods for Research Workers
, he’s pretty clear about what it means to reject a null hypothesis at
\(p<.05\)
. In his opinion, if we take
\(p<.05\)
to mean there is “a real effect”, then “we shall not often be astray”. This view is hardly unusual: in my experience, most practitioners express views very similar to Fisher’s. In essence, the
\(p<.05\)
convention is assumed to represent a fairly stringent evidentiary standard.
Well, how true is that? One way to approach this question is to try to convert
\(p\)
-values to Bayes factors, and see how the two compare. It’s not an easy thing to do because a
\(p\)
-value is a fundamentally different kind of calculation to a Bayes factor, and they don’t measure the same thing. However, there have been some attempts to work out the relationship between the two, and it’s somewhat surprising. For example,
Johnson (
2013
)
presents a pretty compelling case that (for
\(t\)
-tests at least) the
\(p<.05\)
threshold corresponds roughly to a Bayes factor of somewhere between 3:1 and 5:1 in favour of the alternative. If that’s right, then Fisher’s claim is a bit of a stretch. Let’s suppose that the null hypothesis is true about half the time (i.e., the prior probability of
\(H_0\)
is 0.5), and we use those numbers to work out the posterior probability of the null hypothesis given that it has been rejected at
\(p<.05\)
. Using the data from
Johnson (
2013
)
, we see that if you reject the null at
\(p<.05\)
, you’ll be correct about 80% of the time. I don’t know about you, but in my opinion an evidentiary standard that ensures you’ll be wrong on 20% of your decisions isn’t good enough. The fact remains that, quite contrary to Fisher’s claim, if you reject at
\(p<.05\)
you shall quite often go astray. It’s not a very stringent evidentiary threshold at all.
The
\(p\)
-value is a lie.
The cake is a lie.
The cake is a lie.
The cake is a lie.
The cake is a lie.
– Portal
263
Okay, at this point you might be thinking that the real problem is not with orthodox statistics, just the
\(p<.05\)
standard. In one sense, that’s true. The recommendation that
Johnson (
2013
)
gives is not that “everyone must be a Bayesian now”. Instead, the suggestion is that it would be wiser to shift the conventional standard to something like a
\(p<.01\)
level. That’s not an unreasonable view to take, but in my view the problem is a little more severe than that. In my opinion, there’s a fairly big problem built into the way most (but not all) orthodox hypothesis tests are constructed. They are grossly naive about how humans actually do research, and because of this most
\(p\)
-values are wrong.
Sounds like an absurd claim, right? Well, consider the following scenario. You’ve come up with a really exciting research hypothesis and you design a study to test it. You’re very diligent, so you run a power analysis to work out what your sample size should be, and you run the study. You run your hypothesis test and out pops a
\(p\)
-value of 0.072. Really bloody annoying, right?
What should you do? Here are some possibilities:
You conclude that there is no effect, and try to publish it as a null result
You guess that there might be an effect, and try to publish it as a “borderline significant” result
You give up and try a new study
You collect some more data to see if the
\(p\)
value goes up or (preferably!) drops below the “magic” criterion of
\(p<.05\)
Which would
you
choose? Before reading any further, I urge you to take some time to think about it. Be honest with yourself. But don’t stress about it too much, because you’re screwed no matter what you choose. Based on my own experiences as an author, reviewer and editor, as well as stories I’ve heard from others, here’s what will happen in each case:
Let’s start with option 1. If you try to publish it as a null result, the paper will struggle to be published. Some reviewers will think that
\(p=.072\)
is not really a null result. They’ll argue it’s borderline significant. Other reviewers will agree it’s a null result, but will claim that even though some null results
are
publishable, yours isn’t. One or two reviewers might even be on your side, but you’ll be fighting an uphill battle to get it through.
Okay, let’s think about option number 2. Suppose you try to publish it as a borderline significant result. Some reviewers will claim that it’s a null result and should not be published. Others will claim that the evidence is ambiguous, and that you should collect more data until you get a clear significant result. Again, the publication process does not favour you.
Given the difficulties in publishing an “ambiguous” result like
\(p=.072\)
, option number 3 might seem tempting: give up and do something else. But that’s a recipe for career suicide. If you give up and try a new project else every time you find yourself faced with ambiguity, your work will never be published. And if you’re in academia without a publication record you can lose your job. So that option is out.
It looks like you’re stuck with option 4. You don’t have conclusive results, so you decide to collect some more data and re-run the analysis. Seems sensible, but unfortunately for you, if you do this all of your
\(p\)
-values are now incorrect.
All
of them. Not just the
\(p\)
-values that you calculated for
this
study. All of them. All the
\(p\)
-values you calculated in the past and all the
\(p\)
-values you will calculate in the future. Fortunately, no-one will notice. You’ll get published, and you’ll have lied.
Wait, what? How can that last part be true? I mean, it sounds like a perfectly reasonable strategy doesn’t it? You collected some data, the results weren’t conclusive, so now what you want to do is collect more data until the the results
are
conclusive. What’s wrong with that?
Honestly, there’s nothing wrong with it. It’s a reasonable, sensible and rational thing to do. In real life, this is exactly what every researcher does. Unfortunately, the theory of null hypothesis testing as I described it in Chapter
11
forbids
you from doing this.
264
The reason is that the theory assumes that the experiment is finished and all the data are in. And because it assumes the experiment is over, it only considers
two
possible decisions. If you’re using the conventional
\(p<.05\)
threshold, those decisions are:
Outcome
Action
\(p\)
less than .05
Reject the null
\(p\)
greater than .05
Retain the null
What
you’re
doing is adding a third possible action to the decision making problem. Specifically, what you’re doing is using the
\(p\)
-value itself as a reason to justify continuing the experiment. And as a consequence you’ve transformed the decision-making procedure into one that looks more like this:
Outcome
Action
\(p\)
less than .05
Stop the experiment and reject the null
\(p\)
between .05 and .1
Continue the experiment
\(p\)
greater than .1
Stop the experiment and retain the null
The “basic” theory of null hypothesis testing isn’t built to handle this sort of thing, not in the form I described back in Chapter
11
. If you’re the kind of person who would choose to “collect more data” in real life, it implies that you are
not
making decisions in accordance with the rules of null hypothesis testing. Even if you happen to arrive at the same decision as the hypothesis test, you aren’t following the decision
process
it implies, and it’s this failure to follow the process that is causing the problem.
265
Your
\(p\)
-values are a lie.
Worse yet, they’re a lie in a dangerous way, because they’re all
too small
. To give you a sense of just how bad it can be, consider the following (worst case) scenario. Imagine you’re a really super-enthusiastic researcher on a tight budget who didn’t pay any attention to my warnings above. You design a study comparing two groups. You desperately want to see a significant result at the
\(p<.05\)
level, but you really don’t want to collect any more data than you have to (because it’s expensive). In order to cut costs, you start collecting data, but every time a new observation arrives you run a
\(t\)
-test on your data. If the
\(t\)
-tests says
\(p<.05\)
then you stop the experiment and report a significant result. If not, you keep collecting data. You keep doing this until you reach your pre-defined spending limit for this experiment. Let’s say that limit kicks in at
\(N=1000\)
observations. As it turns out, the truth of the matter is that there is no real effect to be found: the null hypothesis is true. So, what’s the chance that you’ll make it to the end of the experiment and (correctly) conclude that there is no effect? In an ideal world, the answer here should be 95%. After all, the whole
point
of the
\(p<.05\)
criterion is to control the Type I error rate at 5%, so what we’d hope is that there’s only a 5% chance of falsely rejecting the null hypothesis in this situation. However, there’s no guarantee that will be true. You’re breaking the rules: you’re running tests repeatedly, “peeking” at your data to see if you’ve gotten a significant result, and all bets are off.
Figure 17.1: How badly can things go wrong if you re-run your tests every time new data arrive? If you are a frequentist, the answer is “very wrong”.
So how bad is it? The answer is shown as the solid black line in Figure
17.1
, and it’s
astoundingly
bad. If you peek at your data after every single observation, there is a 49% chance that you will make a Type I error. That’s, um, quite a bit bigger than the 5% that it’s supposed to be. By way of comparison, imagine that you had used the following strategy. Start collecting data. Every single time an observation arrives, run a
Bayesian
\(t\)
-test (Section
17.7
and look at the Bayes factor. I’ll assume that
Johnson (
2013
)
is right, and I’ll treat a Bayes factor of 3:1 as roughly equivalent to a
\(p\)
-value of .05.
266
This time around, our trigger happy researcher uses the following procedure: if the Bayes factor is 3:1 or more in favour of the null, stop the experiment and retain the null. If it is 3:1 or more in favour of the alternative, stop the experiment and reject the null. Otherwise continue testing. Now, just like last time, let’s assume that the null hypothesis is true. What happens? As it happens, I ran the simulations for this scenario too, and the results are shown as the dashed line in Figure
17.1
. It turns out that the Type I error rate is much much lower than the 49% rate that we were getting by using the orthodox
\(t\)
-test.
In some ways, this is remarkable. The entire
point
of orthodox null hypothesis testing is to control the Type I error rate. Bayesian methods aren’t actually designed to do this at all. Yet, as it turns out, when faced with a “trigger happy” researcher who keeps running hypothesis tests as the data come in, the Bayesian approach is much more effective. Even the 3:1 standard, which most Bayesians would consider unacceptably lax, is much safer than the
\(p<.05\)
rule.
Is it really this bad?
The example I gave in the previous section is a pretty extreme situation. In real life, people don’t run hypothesis tests every time a new observation arrives. So it’s not fair to say that the
\(p<.05\)
threshold “really” corresponds to a 49% Type I error rate (i.e.,
\(p=.49\)
). But the fact remains that if you want your
\(p\)
-values to be honest, then you either have to switch to a completely different way of doing hypothesis tests, or you must enforce a strict rule:
no peeking
. You are
not
allowed to use the data to decide when to terminate the experiment. You are
not
allowed to look at a “borderline”
\(p\)
-value and decide to collect more data. You aren’t even allowed to change your data analyis strategy after looking at data. You are strictly required to follow these rules, otherwise the
\(p\)
-values you calculate will be nonsense.
And yes, these rules are surprisingly strict. As a class exercise a couple of years back, I asked students to think about this scenario. Suppose you started running your study with the intention of collecting
\(N=80\)
people. When the study starts out you follow the rules, refusing to look at the data or run any tests. But when you reach
\(N=50\)
your willpower gives in… and you take a peek. Guess what? You’ve got a significant result! Now, sure, you know you
said
that you’d keep running the study out to a sample size of
\(N=80\)
, but it seems sort of pointless now, right? The result is significant with a sample size of
\(N=50\)
, so wouldn’t it be wasteful and inefficient to keep collecting data? Aren’t you tempted to stop? Just a little? Well, keep in mind that if you do, your Type I error rate at
\(p<.05\)
just ballooned out to 8%. When you report
\(p<.05\)
in your paper, what you’re
really
saying is
\(p<.08\)
. That’s how bad the consequences of “just one peek” can be.
Now consider this … the scientific literature is filled with
\(t\)
-tests, ANOVAs, regressions and chi-square tests. When I wrote this book I didn’t pick these tests arbitrarily. The reason why these four tools appear in most introductory statistics texts is that these are the bread and butter tools of science. None of these tools include a correction to deal with “data peeking”: they all assume that you’re not doing it. But how realistic is that assumption? In real life, how many people do you think have “peeked” at their data before the experiment was finished and adapted their subsequent behaviour after seeing what the data looked like? Except when the sampling procedure is fixed by an external constraint, I’m guessing the answer is “most people have done it”. If that has happened, you can infer that the reported
\(p\)
-values are wrong. Worse yet, because we don’t know what decision process they actually followed, we have no way to know what the
\(p\)
-values
should
have been. You can’t compute a
\(p\)
-value when you don’t know the decision making procedure that the researcher used. And so the reported
\(p\)
-value remains a lie.
Given all of the above, what is the take home message? It’s not that Bayesian methods are foolproof. If a researcher is determined to cheat, they can always do so. Bayes’ rule cannot stop people from lying, nor can it stop them from rigging an experiment. That’s not my point here. My point is the same one I made at the very beginning of the book in Section
1.1
: the reason why we run statistical tests is to protect us from ourselves. And the reason why “data peeking” is such a concern is that it’s so tempting,
even for honest researchers
. A theory for statistical inference has to acknowledge this. Yes, you might try to defend
\(p\)
-values by saying that it’s the fault of the researcher for not using them properly. But to my mind that misses the point. A theory of statistical inference that is so completely naive about humans that it doesn’t even consider the possibility that the researcher might
look at their own data
isn’t a theory worth having. In essence, my point is this:
Good laws have their origins in bad morals.
– Ambrosius Macrobius
267
Good rules for statistical testing have to acknowledge human frailty. None of us are without sin. None of us are beyond temptation. A good system for statistical inference should still work even when it is used by actual humans. Orthodox null hypothesis testing does not.
268
Bayesian analysis of contingency tables
Time to change gears. Up to this point I’ve been talking about what Bayesian inference is and why you might consider using it. I now want to briefly describe how to do Bayesian versions of various statistical tests. The discussions in the next few sections are not as detailed as I’d like, but I hope they’re enough to help you get started. So let’s begin.
The first kind of statistical inference problem I discussed in this book appeared in Chapter
12
, in which we discussed categorical data analysis problems. In that chapter I talked about several different statistical problems that you might be interested in, but the one that appears most often in real life is the analysis of
contingency tables
. In this kind of data analysis situation, we have a cross-tabulation of one variable against another one, and the goal is to find out if there is some
association
between these variables. The data set I used to illustrate this problem is found in the
chapek9.Rdata
file, and it contains a single data frame
chapek9
load
(
file.path
(projecthome,
"data"
,
"chapek9.Rdata"
))
head
(chapek9)
## species choice
## 1 robot flower
## 2 human data
## 3 human data
## 4 human data
## 5 robot data
## 6 human flower
In this data set, we supposedly sampled 180 beings and measured two things. First, we checked whether they were humans or robots, as captured by the
species
variable. Second, we asked them to nominate whether they most preferred flowers, puppies, or data. When we produce the cross-tabulation, we get this as the results:
crosstab <-
xtabs
(
~
species
+
choice, chapek9 )
crosstab
## choice
## species puppy flower data
## robot 13 30 44
## human 15 13 65
Surprisingly, the humans seemed to show a much stronger preference for data than the robots did. At the time we speculated that this might have been because the questioner was a large robot carrying a gun, and the humans might have been scared.
The orthodox text
Just to refresh your memory, here’s how we analysed these data back in
Chapter@refch
:chisquare. Because we want to determine if there is some
association
between
species
and
choice
, we used the
associationTest()
function in the
lsr
package to run a chi-square test of association. The results looked like this:
library
(lsr)
associationTest
(
~
species
+
choice, chapek9 )
##
## Chi-square test of categorical association
##
## Variables: species, choice
##
## Hypotheses:
## null: variables are independent of one another
## alternative: some contingency exists between variables
##
## Observed contingency table:
## choice
## species puppy flower data
## robot 13 30 44
## human 15 13 65
##
## Expected contingency table under the null hypothesis:
## choice
## species puppy flower data
## robot 13.5 20.8 52.7
## human 14.5 22.2 56.3
##
## Test results:
## X-squared statistic: 10.722
## degrees of freedom: 2
## p-value: 0.005
##
## Other information:
## estimated effect size (Cramer's v): 0.244
Because we found a small
\(p\)
value (in this case
\(p<.01\)
), we concluded that the data are inconsistent with the null hypothesis of no association, and we rejected it.
The Bayesian test
How do we run an equivalent test as a Bayesian? Well, like every other bloody thing in statistics, there’s a lot of different ways you
could
do it. However, for the sake of everyone’s sanity, throughout this chapter I’ve decided to rely on one R package to do the work. Specifically, I’m going to use the
BayesFactor
package written by Jeff Rouder and Rich Morey, which as of this writing is in version 0.9.10.
For the analysis of contingency tables, the
BayesFactor
package contains a function called
contingencyTableBF()
. The data that you need to give to this function is the contingency table itself (i.e., the
crosstab
variable above), so you might be expecting to use a command like this:
library
( BayesFactor )
# ...because we have to load the package
contingencyTableBF
( crosstab )
# ...because that makes sense, right?
However, if you try this you’ll get an error message. This is because the
contingencyTestBF()
function needs one other piece of information from you: it needs to know what
sampling plan
you used to run your experiment. You can specify the sampling plan using the
sampleType
argument. So I should probably tell you what your options are! The
contingencyTableBF()
function distinguishes between four different types of experiment:
Fixed sample size
. Suppose that in our
chapek9
example, our experiment was designed like this: we deliberately set out to test 180 people, but we didn’t try to control the number of humans or robots, nor did we try to control the choices they made. In this design, the total number of observations
\(N\)
is fixed, but everything else is random. This is referred to as “joint multinomial” sampling, and if that’s what you did you should specify
sampleType = "jointMulti"
. In the case of the
chapek9
data, that’s actually what I had in mind when I invented the data set.
Fixed row (or column) totals
. A different kind of design might work like this. We decide ahead of time that we want 180 people, but we try to be a little more systematic about it. Specifically, the
experimenter
constrains it so that we get a predetermined number of humans and robots (e.g., 90 of each). In this design,
either
the row totals or the column totals are fixed, but not both. This is referred to as “independent multinomial” sampling, and if that’s what you did you should specify
sampleType = "indepMulti"
.
Both row and column totals fixed
. Another logical possibility is that you designed the experiment so that
both
the row totals and the column totals are fixed. This doesn’t make any sense at all in the
chapek9
example, but there are other deisgns that can work this way. Suppose that I show you a collection of 20 toys, and then given them 10 stickers that say
boy
and another 10 that say
girl
. I then give them 10
blue
stickers and 10
pink
stickers. I then ask you to put the stickers on the 20 toys such that every toy has a colour and every toy has a gender. No matter how you assign the stickers, the total number of pink and blue toys will be 10, as will the number of boys and girls. In this design
both
the rows and columns of the contingency table are fixed. This is referred to as “hypergeometric” sampling, and if that’s what you’ve done you should specify
sampleType = "hypergeom"
.
Nothing is fixed
. Finally, it might be the case that
nothing
is fixed. Not the row columns, not the column totals, and not the total sample size either. For instance, in the
chapek9
scenario, suppose what I’d done is run the study for a fixed length of
time
. By chance, it turned out that I got 180 people to turn up to study, but it could easily have been something else. This is referred to as “Poisson” sampling, and if that’s what you’ve done you should specify
sampleType="poisson"
.
Okay, so now we have enough knowledge to actually run a test. For the
chapek9
data, I implied that we designed the study such that the total sample size
\(N\)
was fixed, so we should set
sampleType = "jointMulti"
. The command that we need is,
library
( BayesFactor )
contingencyTableBF
( crosstab,
sampleType =
"jointMulti"
)
## Bayes factor analysis
## --------------
## [1] Non-indep. (a=1) : 15.92684 ±0%
##
## Against denominator:
## Null, independence, a = 1
## ---
## Bayes factor type: BFcontingencyTable, joint multinomial
As with most R commands, the output initially looks suspiciously similar to utter gibberish. Fortunately, it’s actually pretty simple once you get past the initial impression. Firstly, note that the stuff at the top and bottom are irrelevant fluff. You already know that you’re doing a Bayes factor analysis. You already know that you’re analysing a contingency table, and you already know that you specified a joint multinomial sampling plan. So let’s strip that out and take a look at what’s left over:
[1] Non-indep. (a=1) : 15.92684 @plusorminus0%
Against denominator:
Null, independence, a = 1
Let’s also ignore those two
a=1
bits, since they’re technical details that you don’t need to know about at this stage.
269
The rest of the output is actually pretty straightforward. At the bottom, the output defines the null hypothesis for you: in this case, the null hypothesis is that there is no relationship between
species
and
choice
. Or, to put it another way, the null hypothesis is that these two variables are
independent
. Now if you look at the line above it, you might (correctly) guess that the
Non-indep.
part refers to the
alternative
hypothesis. In this case, the alternative is that there
is
a relationship between
species
and
choice
: that is, they are not independent. So the only thing left in the output is the bit that reads
15.92684 @plusorminus0%
The 15.9 part is the Bayes factor, and it’s telling you that the odds for the alternative hypothesis against the null are about 16:1. The
\(\pm0\%\)
part is not very interesting: essentially, all it’s telling you is that R has calculated an exact Bayes factor, so the uncertainty about the Bayes factor is 0%.
270
In any case, the data are telling us that we have moderate evidence for the alternative hypothesis.
Writing up the results
When writing up the results, my experience has been that there aren’t quite so many “rules” for how you “should” report Bayesian hypothesis tests. That might change in the future if Bayesian methods become standard and some task force starts writing up style guides, but in the meantime I would suggest using some common sense. For example, I would avoid writing this:
A Bayesian test of association found a significant result (BF=15.92)
To my mind, this write up is unclear. Even assuming that you’ve already reported the relevant descriptive statistics, there are a number of things I am unhappy with. First, the concept of “statistical significance” is pretty closely tied with
\(p\)
-values, so it reads slightly strangely. Second, the “BF=15.92” part will only make sense to people who already understand Bayesian methods, and not everyone does. Third, it is somewhat unclear exactly which test was run and what software was used to do so.
On the other hand, unless precision is
extremely
important, I think that this is taking things a step too far:
We ran a Bayesian test of association using version 0.9.10-1 of the BayesFactor package using default priors and a joint multinomial sampling plan. The resulting Bayes factor of 15.92 to 1 in favour of the alternative hypothesis indicates that there is moderately strong evidence for the non-independence of species and choice.
Everything about that passage is correct, of course.
Morey and Rouder (
2015
)
built their Bayesian tests of association using the paper by
Gunel and Dickey (
1974
)
, the specific test we used assumes that the experiment relied on a joint multinomial sampling plan, and indeed the Bayes factor of 15.92 is moderately strong evidence. It’s just far too wordy.
In most situations you just don’t need that much information. My preference is usually to go for something a little briefer. First, if you’re reporting multiple Bayes factor analyses in your write up, then somewhere you only need to cite the software once, at the beginning of the results section. So you might have one sentence like this:
All analyses were conducted using the BayesFactor package in R , and unless otherwise stated default parameter values were used
Notice that I don’t bother including the version number? That’s because the citation itself includes that information (go check my reference list if you don’t believe me). There’s no need to clutter up your results with redundant information that almost no-one will actually need. When you get to the actual test you can get away with this:
A test of association produced a Bayes factor of 16:1 in favour of a relationship between species and choice.
Short and sweet. I’ve rounded 15.92 to 16, because there’s not really any important difference between 15.92:1 and 16:1. I spelled out “Bayes factor” rather than truncating it to “BF” because not everyone knows the abbreviation. I indicated exactly what the effect is (i.e., “a relationship between species and choice”) and how strong the evidence was. I
didn’t
bother indicating whether this was “moderate” evidence or “strong” evidence, because the odds themselves tell you! There’s nothing stopping you from including that information, and I’ve done so myself on occasions, but you don’t strictly need it. Similarly, I didn’t bother to indicate that I ran the “joint multinomial” sampling plan, because I’m assuming that the method section of my write up would make clear how the experiment was designed. (I might change my mind about that if the method section was ambiguous.) Neither did I bother indicating that this was a
Bayesian
test of association: if your reader can’t work that out from the fact that you’re reporting a Bayes factor and the fact that you’re citing the
BayesFactor
package for all your analyses, then there’s no chance they’ll understand anything you’ve written. Besides, if you keep writing the word “Bayes” over and over again it starts to look stupid. Bayes Bayes Bayes Bayes Bayes. See?
Other sampling plans
Up to this point all I’ve shown you is how to use the
contingencyTableBF()
function for the joint multinomial sampling plan (i.e., when the total sample size
\(N\)
is fixed, but nothing else is). For the Poisson sampling plan (i.e., nothing fixed), the command you need is identical except for the
sampleType
argument:
contingencyTableBF
(crosstab,
sampleType =
"poisson"
)
## Bayes factor analysis
## --------------
## [1] Non-indep. (a=1) : 28.20757 ±0%
##
## Against denominator:
## Null, independence, a = 1
## ---
## Bayes factor type: BFcontingencyTable, poisson
Notice that the Bayes factor of 28:1 here is
not
the identical to the Bayes factor of 16:1 that we obtained from the last test. The sampling plan actually does matter.
What about the design in which the row columns (or column totals) are fixed? As I mentioned earlier, this corresponds to the “independent multinomial” sampling plan. Again, you need to specify the
sampleType
argument, but this time you need to specify whether you fixed the rows or the columns. For example, suppose I deliberately sampled 87 humans and 93 robots, then I would need to indicate that the
fixedMargin
of the contingency table is the
"rows"
. So the command I would use is:
contingencyTableBF
(crosstab,
sampleType =
"indepMulti"
,
fixedMargin=
"rows"
)
## Bayes factor analysis
## --------------
## [1] Non-indep. (a=1) : 8.605897 ±0%
##
## Against denominator:
## Null, independence, a = 1
## ---
## Bayes factor type: BFcontingencyTable, independent multinomial
Again, the Bayes factor is different, with the evidence for the alternative dropping to a mere 9:1. As you might expect, the answers would be diffrent again if it were the columns of the contingency table that the experimental design fixed.
Finally, if we turn to hypergeometric sampling in which everything is fixed, we get…
contingencyTableBF
(crosstab,
sampleType =
"hypergeom"
)
#Error in contingencyHypergeometric(as.matrix(data2), a) :
# hypergeometric contingency tables restricted to 2 x 2 tables; see help for contingencyTableBF()
… an error message. Okay, some quick reading through the help files hints that support for larger contingency tables is coming, but it’s not been implemented yet. In the meantime, let’s imagine we have data from the “toy labelling” experiment I described earlier in this section. Specifically, let’s say our data look like this:
toys <-
data.frame
(
stringsAsFactors=
FALSE
,
gender =
c
(
"girl"
,
"boy"
),
pink =
c
(
8
,
2
),
blue =
c
(
2
,
8
)
)
The Bayesian test with hypergeometric sampling gives us this:
contingencyTableBF
(toys,
sampleType =
"hypergeom"
)
#Bayes factor analysis
#--------------
#[1] Non-indep. (a=1) : 8.294321 @plusorminus0%
#
#Against denominator:
# Null, independence, a = 1
#---
#Bayes factor type: BFcontingencyTable, hypergeometric
The Bayes factor of 8:1 provides modest evidence that the labels were being assigned in a way that correlates gender with colour, but it’s not conclusive.
Bayesian
\(t\)
-tests
The second type of statistical inference problem discussed in this book is the comparison between two means, discussed in some detail in the chapter on
\(t\)
-tests (Chapter
13
. If you can remember back that far, you’ll recall that there are several versions of the
\(t\)
-test. The
BayesFactor
package contains a function called
ttestBF()
that is flexible enough to run several different versions of the
\(t\)
-test. I’ll talk a little about Bayesian versions of the independent samples
\(t\)
-tests and the paired samples
\(t\)
-test in this section.
Independent samples
\(t\)
-test
The most common type of
\(t\)
-test is the independent samples
\(t\)
-test, and it arises when you have data that look something like this:
load
(
file.path
(projecthome,
"data"
,
"harpo.Rdata"
))
head
(harpo)
## grade tutor
## 1 65 Anastasia
## 2 72 Bernadette
## 3 66 Bernadette
## 4 74 Anastasia
## 5 73 Anastasia
## 6 71 Bernadette
In this data set, we have two groups of students, those who received lessons from Anastasia and those who took their classes with Bernadette. The question we want to answer is whether there’s any difference in the grades received by these two groups of student. Back in
Chapter@refch
:ttest I suggested you could analyse this kind of data using the
independentSamplesTTest()
function in the
lsr
package. For example, if you want to run a Student’s
\(t\)
-test, you’d use a command like this:
independentSamplesTTest
(
formula =
grade
~
tutor,
data =
harpo,
var.equal =
TRUE
)
##
## Student's independent samples t-test
##
## Outcome variable: grade
## Grouping variable: tutor
##
## Descriptive statistics:
## Anastasia Bernadette
## mean 74.533 69.056
## std dev. 8.999 5.775
##
## Hypotheses:
## null: population means equal for both groups
## alternative: different population means in each group
##
## Test results:
## t-statistic: 2.115
## degrees of freedom: 31
## p-value: 0.043
##
## Other information:
## two-sided 95% confidence interval: [0.197, 10.759]
## estimated effect size (Cohen's d): 0.74
Like most of the functions that I wrote for this book, the
independentSamplesTTest()
is very wordy. It prints out a bunch of descriptive statistics and a reminder of what the null and alternative hypotheses are, before finally getting to the test results. I wrote it that way deliberately, in order to help make things a little clearer for people who are new to statistics.
Again, we obtain a
\(p\)
-value less than 0.05, so we reject the null hypothesis.
What does the Bayesian version of the
\(t\)
-test look like? Using the
ttestBF()
function, we can obtain a Bayesian analog of Student’s independent samples
\(t\)
-test using the following command:
ttestBF
(
formula =
grade
~
tutor,
data =
harpo )
## Bayes factor analysis
## --------------
## [1] Alt., r=0.707 : 1.754927 ±0%
##
## Against denominator:
## Null, mu1-mu2 = 0
## ---
## Bayes factor type: BFindepSample, JZS
Notice that format of this command is pretty standard. As usual we have a
formula
argument in which we specify the outcome variable on the left hand side and the grouping variable on the right. The
data
argument is used to specify the data frame containing the variables. However, notice that there’s no analog of the
var.equal
argument. This is because the
BayesFactor
package does not include an analog of the Welch test, only the Student test.
271
In any case, when you run this command you get this as the output:
So what does all this mean? Just as we saw with the
contingencyTableBF()
function, the output is pretty dense. But, just like last time, there’s not a lot of information here that you actually need to process. Firstly, let’s examine the bottom line. The
BFindepSample
part just tells you that you ran an independent samples
\(t\)
-test, and the
JZS
part is technical information that is a little beyond the scope of this book.
272
Clearly, there’s nothing to worry about in that part. In the line above, the text
Null, mu1-mu2 = 0
is just telling you that the null hypothesis is that there are no differences between means. But you already knew that. So the only part that really matters is this line here:
[1] Alt., r=0.707 : 1.754927 @plusorminus0%
Ignore the
r=0.707
part: it refers to a technical detail that we won’t worry about in this chapter.
273
Instead, you should focus on the part that reads
1.754927
. This is the Bayes factor: the evidence provided by these data are about 1.8:1 in favour of the alternative.
Before moving on, it’s worth highlighting the difference between the orthodox test results and the Bayesian one. According to the orthodox test, we obtained a significant result, though only barely. Nevertheless, many people would happily accept
\(p=.043\)
as reasonably strong evidence for an effect. In contrast, notice that the Bayesian test doesn’t even reach 2:1 odds in favour of an effect, and would be considered very weak evidence at best. In my experience that’s a pretty typical outcome. Bayesian methods usually require more evidence before rejecting the null.
Paired samples
\(t\)
-test
Back in Section
13.5
I discussed the
chico
data frame in which students grades were measured on two tests, and we were interested in finding out whether grades went up from test 1 to test 2. Because every student did both tests, the tool we used to analyse the data was a paired samples
\(t\)
-test. To remind you of what the data look like, here’s the first few cases:
load
(
file.path
(projecthome,
"data"
,
"chico.Rdata"
))
head
(chico)
## id grade_test1 grade_test2
## 1 student1 42.9 44.6
## 2 student2 51.8 54.0
## 3 student3 71.7 72.3
## 4 student4 51.6 53.4
## 5 student5 63.5 63.8
## 6 student6 58.0 59.3
We originally analysed the data using the
pairedSamplesTTest()
function in the
lsr
package, but this time we’ll use the
ttestBF()
function from the
BayesFactor
package to do the same thing. The easiest way to do it with this data set is to use the
x
argument to specify one variable and the
y
argument to specify the other. All we need to do then is specify
paired=TRUE
to tell R that this is a paired samples test. So here’s our command:
ttestBF
(
x =
chico
$
grade_test1,
y =
chico
$
grade_test2,
paired =
TRUE
)
## Bayes factor analysis
## --------------
## [1] Alt., r=0.707 : 5992.05 ±0%
##
## Against denominator:
## Null, mu = 0
## ---
## Bayes factor type: BFoneSample, JZS
At this point, I hope you can read this output without any difficulty. The data provide evidence of about 6000:1 in favour of the alternative. We could probably reject the null with some confidence!
Bayesian regression
Okay, so now we’ve seen Bayesian equivalents to orthodox chi-square tests and
\(t\)
-tests. What’s next? If I were to follow the same progression that I used when developing the orthodox tests you’d expect to see ANOVA next, but I think it’s a little clearer if we start with regression.
A quick refresher
In Chapter
15
I used the
parenthood
data to illustrate the basic ideas behind regression. To remind you of what that data set looks like, here’s the first six observations:
load
(
file.path
(projecthome,
"data"
,
"parenthood.Rdata"
))
head
(parenthood)
## dan.sleep baby.sleep dan.grump day
## 1 7.59 10.18 56 1
## 2 7.91 11.66 60 2
## 3 5.14 7.92 82 3
## 4 7.71 9.61 55 4
## 5 6.68 9.75 67 5
## 6 5.99 5.04 72 6
Back in Chapter
15
I proposed a theory in which my grumpiness (
dan.grump
) on any given day is related to the amount of sleep I got the night before (
dan.sleep
), and possibly to the amount of sleep our baby got (
baby.sleep
), though probably not to the
day
on which we took the measurement. We tested this using a regression model. In order to estimate the regression model we used the
lm()
function, like so:
model <-
lm
(
formula =
dan.grump
~
dan.sleep
+
day
+
baby.sleep,
data =
parenthood
)
The hypothesis tests for each of the terms in the regression model were extracted using the
summary()
function as shown below:
summary
(model)
##
## Call:
## lm(formula = dan.grump ~ dan.sleep + day + baby.sleep, data = parenthood)
##
## Residuals:
## Min 1Q Median 3Q Max
## -10.906 -2.284 -0.295 2.652 11.880
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 126.278707 3.242492 38.945 <2e-16 ***
## dan.sleep -8.969319 0.560007 -16.016 <2e-16 ***
## day -0.004403 0.015262 -0.288 0.774
## baby.sleep 0.015747 0.272955 0.058 0.954
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 4.375 on 96 degrees of freedom
## Multiple R-squared: 0.8163, Adjusted R-squared: 0.8105
## F-statistic: 142.2 on 3 and 96 DF, p-value: < 2.2e-16
When interpreting the results, each row in this table corresponds to one of the possible predictors. The
(Intercept)
term isn’t usually interesting, though it is highly significant. The important thing for our purposes is the fact that
dan.sleep
is significant at
\(p<.001\)
and neither of the other variables are.
The Bayesian version
Okay, so how do we do the same thing using the
BayesFactor
package? The easiest way is to use the
regressionBF()
function instead of
lm()
. As before, we use
formula
to indicate what the full regression model looks like, and the
data
argument to specify the data frame. So the command is:
regressionBF
(
formula =
dan.grump
~
dan.sleep
+
day
+
baby.sleep,
data =
parenthood
)
##
|
| | 0%
|
|========= | 14%
|
|=================== | 29%
|
|============================ | 43%
|
|===================================== | 57%
|
|============================================== | 71%
|
|======================================================== | 86%
|
|=================================================================| 100%
## Bayes factor analysis
## --------------
## [1] dan.sleep : 1.622545e+34 ±0.01%
## [2] day : 0.2724027 ±0%
## [3] baby.sleep : 10018411 ±0%
## [4] dan.sleep + day : 1.016576e+33 ±0%
## [5] dan.sleep + baby.sleep : 9.77022e+32 ±0%
## [6] day + baby.sleep : 2340755 ±0%
## [7] dan.sleep + day + baby.sleep : 7.835625e+31 ±0%
##
## Against denominator:
## Intercept only
## ---
## Bayes factor type: BFlinearModel, JZS
So that’s pretty straightforward: it’s exactly what we’ve been doing throughout the book. The output, however, is a little different from what you get from
lm()
. The format of this is pretty familiar. At the bottom we have some techical rubbish, and at the top we have some information about the Bayes factors. What’s new is the fact that we seem to have
lots
of Bayes factors here. What’s all this about?
The trick to understanding this output is to recognise that if we’re interested in working out which of the 3 predictor variables are related to
dan.grump
, there are actually 8 possible regression models that could be considered. One possibility is the
intercept only model
, in which none of the three variables have an effect. At the other end of the spectrum is the
full model
in which all three variables matter. So what
regressionBF()
does is treat the
intercept only
model as the null hypothesis, and print out the Bayes factors for all other models when compared against that null. For example, if we look at line 4 in the table, we see that the evidence is about
\(10^{33}\)
to 1 in favour of the claim that a model that includes both
dan.sleep
and
day
is better than the intercept only model. Or if we look at line 1, we can see that the odds are about
\(1.6 \times 10^{34}\)
that a model containing the
dan.sleep
variable (but no others) is better than the intercept only model.
Finding the best model
In practice, this isn’t super helpful. In most situations the intercept only model is one that you don’t really care about at all. What I find helpful is to start out by working out which model is the
best
one, and then seeing how well all the alternatives compare to it. Here’s how you do that. In this case, it’s easy enough to see that the best model is actually the one that contains
dan.sleep
only (line 1), because it has the largest Bayes factor. However, if you’ve got a lot of possible models in the output, it’s handy to know that you can use the
head()
function to pick out the best few models. First, we have to go back and save the Bayes factor information to a variable:
models <-
regressionBF
(
formula =
dan.grump
~
dan.sleep
+
day
+
baby.sleep,
data =
parenthood
)
##
|
| | 0%
|
|========= | 14%
|
|=================== | 29%
|
|============================ | 43%
|
|===================================== | 57%
|
|============================================== | 71%
|
|======================================================== | 86%
|
|=================================================================| 100%
Let’s say I want to see the best three models. To do this, I use the
head()
function specifying
n=3
, and here’s what I get as the result:
head
( models,
n =
3
)
## Bayes factor analysis
## --------------
## [1] dan.sleep : 1.622545e+34 ±0.01%
## [2] dan.sleep + day : 1.016576e+33 ±0%
## [3] dan.sleep + baby.sleep : 9.77022e+32 ±0%
##
## Against denominator:
## Intercept only
## ---
## Bayes factor type: BFlinearModel, JZS
This is telling us that the model in line 1 (i.e.,
dan.grump ~ dan.sleep
) is the best one. That’s
almost
what I’m looking for, but it’s still comparing all the models against the intercept only model. That seems silly. What I’d like to know is how big the difference is between the best model and the other good models. For that, there’s this trick:
head
( models
/
max
(models),
n =
3
)
## Bayes factor analysis
## --------------
## [1] dan.sleep : 1 ±0%
## [2] dan.sleep + day : 0.0626532 ±0.01%
## [3] dan.sleep + baby.sleep : 0.0602154 ±0.01%
##
## Against denominator:
## dan.grump ~ dan.sleep
## ---
## Bayes factor type: BFlinearModel, JZS
Notice the bit at the bottom showing that the “denominator” has changed. What that means is that the Bayes factors are now comparing each of those 3 models listed against the
dan.grump ~ dan.sleep
model. Obviously, the Bayes factor in the first line is exactly 1, since that’s just comparing the best model to itself. More to the point, the other two Bayes factors are both less than 1, indicating that they’re all worse than that model. The Bayes factors of 0.06 to 1 imply that the odds for the best model over the second best model are about 16:1. You can work this out by simple arithmetic (i.e.,
\(0.06 / 1 \approx 16\)
), but the other way to do it is to directly compare the models. To see what I mean, here’s the original output:
models
## Bayes factor analysis
## --------------
## [1] dan.sleep : 1.622545e+34 ±0.01%
## [2] day : 0.2724027 ±0%
## [3] baby.sleep : 10018411 ±0%
## [4] dan.sleep + day : 1.016576e+33 ±0%
## [5] dan.sleep + baby.sleep : 9.77022e+32 ±0%
## [6] day + baby.sleep : 2340755 ±0%
## [7] dan.sleep + day + baby.sleep : 7.835625e+31 ±0%
##
## Against denominator:
## Intercept only
## ---
## Bayes factor type: BFlinearModel, JZS
The best model corresponds to row 1 in this table, and the second best model corresponds to row 4. All you have to do to compare these two models is this:
models[
1
]
/
models[
4
]
## Bayes factor analysis
## --------------
## [1] dan.sleep : 15.96088 ±0.01%
##
## Against denominator:
## dan.grump ~ dan.sleep + day
## ---
## Bayes factor type: BFlinearModel, JZS
And there you have it. You’ve found the regression model with the highest Bayes factor (i.e.,
dan.grump ~ dan.sleep
), and you know that the evidence for that model over the next best alternative (i.e.,
dan.grump ~ dan.sleep + day
) is about 16:1.
Bayesian ANOVA
As you can tell, the
BayesFactor
package is pretty flexible, and it can do Bayesian versions of pretty much everything in this book. In fact, it can do a few other neat things that I haven’t covered in the book at all. However, I have to stop somewhere, and so there’s only one other topic I want to cover: Bayesian ANOVA.
A quick refresher
As with the other examples, I think it’s useful to start with a reminder of how I discussed ANOVA earlier in the book. First, let’s remind ourselves of what the data were. The example I used originally is the
clin.trial
data frame, which looks like this
load
(
file.path
(projecthome,
"data"
,
"clinicaltrial.Rdata"
))
head
(clin.trial)
## drug therapy mood.gain
## 1 placebo no.therapy 0.5
## 2 placebo no.therapy 0.3
## 3 placebo no.therapy 0.1
## 4 anxifree no.therapy 0.6
## 5 anxifree no.therapy 0.4
## 6 anxifree no.therapy 0.2
To run our orthodox analysis in earlier chapters we used the
aov()
function to do all the heavy lifting. In Chapter
16
I recommended using the
Anova()
function from the
car
package to produce the ANOVA table, because it uses Type II tests by default. If you’ve forgotten what “Type II tests” are, it might be a good idea to re-read Section
16.10
, because it will become relevant again in a moment. In any case, here’s what our analysis looked like:
library
(car)
model <-
aov
( mood.gain
~
drug
*
therapy,
data =
clin.trial )
Anova
(model)
## Anova Table (Type II tests)
##
## Response: mood.gain
## Sum Sq Df F value Pr(>F)
## drug 3.4533 2 31.7143 1.621e-05 ***
## therapy 0.4672 1 8.5816 0.01262 *
## drug:therapy 0.2711 2 2.4898 0.12460
## Residuals 0.6533 12
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
That’s pretty clearly showing us evidence for a main effect of
drug
at
\(p<.001\)
, an effect of
therapy
at
\(p<.05\)
and no interaction.
The Bayesian version
How do we do the same thing using Bayesian methods? The
BayesFactor
package contains a function called
anovaBF()
that does this for you. It uses a pretty standard
formula
and
data
structure, so the command should look really familiar. Just like we did with regression, it will be useful to save the output to a variable:
models <-
anovaBF
(
formula =
mood.gain
~
drug
*
therapy,
data =
clin.trial
)
##
|
| | 0%
|
|================ | 25%
|
|================================ | 50%
|
|================================================= | 75%
|
|=================================================================| 100%
The output is quite different to the traditional ANOVA, but it’s not too bad once you understand what you’re looking for. Let’s take a look:
models
This looks very similar to the output we obtained from the
regressionBF()
function, and with good reason. Remember what I said back in Section
16.6
: under the hood, ANOVA is no different to regression, and both are just different examples of a linear model. Becasue of this, the
anovaBF()
reports the output in much the same way. For instance, if we want to identify the best model we could use the same commands that we used in the last section. One variant that I find quite useful is this:
models
/
max
(models)
## Bayes factor analysis
## --------------
## [1] drug : 0.3436062 ±1.34%
## [2] therapy : 0.001022285 ±1.34%
## [3] drug + therapy : 1 ±0%
## [4] drug + therapy + drug:therapy : 0.954423 ±1.64%
##
## Against denominator:
## mood.gain ~ drug + therapy
## ---
## Bayes factor type: BFlinearModel, JZS
By “dividing” the
models
output by the best model (i.e.,
max(models)
), what R is doing is using the best model (which in this case is
drugs + therapy
) as the denominator, which gives you a pretty good sense of how close the competitors are. For instance, the model that contains the interaction term is almost as good as the model without the interaction, since the Bayes factor is 0.98. In other words, the data do not clearly indicate whether there is or is not an interaction.
Constructing Bayesian Type II tests
Okay, that’s all well and good, you might be thinking, but what do I report as the alternative to the
\(p\)
-value? In the classical ANOVA table, you get a single
\(p\)
-value for every predictor in the model, so you can talk about the significance of each effect. What’s the Bayesian analog of this?
It’s a good question, but the answer is tricky. Remember what I said in Section
16.10
about ANOVA being complicated. Even in the classical version of ANOVA there are several different “things” that ANOVA might correspond to. Specifically, I discussed how you get different
\(p\)
-values depending on whether you use Type I tests, Type II tests or Type III tests. To work out which Bayes factor is analogous to “the”
\(p\)
-value in a classical ANOVA, you need to work out which version of ANOVA you want an analog for. For the purposes of this section, I’ll assume you want Type II tests, because those are the ones I think are most sensible in general. As I discussed back in Section
16.10
, Type II tests for a two-way ANOVA are reasonably straightforward, but if you have forgotten that section it wouldn’t be a bad idea to read it again before continuing.
Assuming you’ve had a refresher on Type II tests, let’s have a look at how to pull them from the Bayes factor table. Suppose we want to test the main effect of
drug
. The null hypothesis for this test corresponds to a model that includes an effect of
therapy
, but no effect of
drug
. The alternative hypothesis is the model that includes both. In other words, what we want is the Bayes factor corresponding to this comparison:
knitr
::
kable
(tibble
::
tribble
(
~
V1,
~
V2,
"Null model:"
,
"`mood.gain ~ therapy`"
,
"Alternative model:"
,
"`mood.gain ~ therapy + drug`"
),
col.names =
c
(
""
,
""
))
Null model:
mood.gain ~ therapy
Alternative model:
mood.gain ~ therapy + drug
As it happens, we can read the answer to this straight off the table because it corresponds to a comparison between the model in line 2 of the table and the model in line 3: the Bayes factor in this case represents evidence
for
the null of 0.001 to 1. Or, more helpfully, the odds are about 1000 to 1 against the null.
The main effect of
therapy
can be calculated in much the same way. In this case, the null model is the one that contains only an effect of drug, and the alternative is the model that contains both. So the relevant comparison is between lines 2 and 1 in the table. The odds in favour of the null here are only 0.35 to 1. Again, I find it useful to frame things the other way around, so I’d refer to this as evidence of about 3 to 1 in favour of an effect of
therapy
.
Finally, in order to test an interaction effect, the null model here is one that contains both main effects but no interaction. The alternative model adds the interaction. That is:
knitr
::
kable
(tibble
::
tribble
(
~
V1,
~
V2,
"Null model:"
,
"`mood.gain ~ drug + therapy`"
,
"Alternative model:"
,
"`mood.gain ~ drug + therapy + drug:therapy`"
),
col.names =
c
(
""
,
""
))
Null model:
mood.gain ~ drug + therapy
Alternative model:
mood.gain ~ drug + therapy + drug:therapy
If we look those two models up in the table, we see that this comparison is between the models on lines 3 and 4 of the table. The odds of 0.98 to 1 imply that these two models are fairly evenly matched.
You might be thinking that this is all pretty laborious, and I’ll concede that’s true. At some stage I might consider adding a function to the
lsr
package that would automate this process and construct something like a “Bayesian Type II ANOVA table” from the output of the
anovaBF()
function. However, I haven’t had time to do this yet, nor have I made up my mind about whether it’s really a good idea to do this. In the meantime, I thought I should show you the trick for how I do this in practice. The command that I use when I want to grab the right Bayes factors for a Type II ANOVA is this one:
max
(models)
/
models
## denominator
## numerator drug therapy drug + therapy
## drug + therapy 2.910308 978.2007 1
## denominator
## numerator drug + therapy + drug:therapy
## drug + therapy 1.047753
The output isn’t quite so pretty as the last one, but the nice thing is that you can read off everything you need. The best model is
drug + therapy
, so all the other models are being compared to that. What’s the Bayes factor
for
the main effect of
drug
? The relevant null hypothesis is the one that contains only
therapy
, and the Bayes factor in question is 954:1. The main effect of
therapy
is weaker, and the evidence here is only 2.8:1. Finally, the evidence
against
an interaction is very weak, at 1.01:1.
Reading the results off this table is sort of counterintuitive, because you have to read off the answers from the “wrong” part of the table. For instance, the evidence for an effect of
drug
can be read from the column labelled
therapy
, which is pretty damned weird. To be fair to the authors of the package, I don’t think they ever intended for the
anovaBF()
function to be used this way. My understanding
274
is that their view is simply that you should find the best model and report that model: there’s no inherent reason why a Bayesian ANOVA should try to follow the exact same design as an orthodox ANOVA.
275
In any case, if you know what you’re looking for, you can look at this table and then report the results of the Bayesian analysis in a way that is pretty closely analogous to how you’d report a regular Type II ANOVA. As I mentioned earlier, there’s still no convention on how to do that, but I usually go for something like this:
A Bayesian Type II ANOVA found evidence for main effects of drug (Bayes factor: 954:1) and therapy (Bayes factor: 3:1), but no clear evidence for or against an interaction (Bayes factor: 1:1).
Summary
The first half of this chapter was focused primarily on the theoretical underpinnings of Bayesian statistics. I introduced the mathematics for how Bayesian inference works (Section
17.1
), and gave a very basic overview of how Bayesian hypothesis testing is typically done (Section
17.2
). Finally, I devoted some space to talking about why I think Bayesian methods are worth using (Section
17.3
.
The second half of the chapter was a lot more practical, and focused on tools provided by the
BayesFactor
package. Specifically, I talked about using the
contingencyTableBF()
function to do Bayesian analogs of chi-square tests (Section
17.6
, the
ttestBF()
function to do Bayesian
\(t\)
-tests, (Section
17.7
), the
regressionBF()
function to do Bayesian regressions, and finally the
anovaBF()
function for Bayesian ANOVA.
If you’re interested in learning more about the Bayesian approach, there are many good books you could look into. John Kruschke’s book
Doing Bayesian Data Analysis
is a pretty good place to start
(Kruschke
2011
)
, and is a nice mix of theory and practice. His approach is a little different to the “Bayes factor” approach that I’ve discussed here, so you won’t be covering the same ground. If you’re a cognitive psychologist, you might want to check out Michael Lee and E.J. Wagenmakers’ book
Bayesian Cognitive Modeling
(Lee and Wagenmakers
2014
)
. I picked these two because I think they’re especially useful for people in my discipline, but there’s a lot of good books out there, so look around!
References
Jeffreys, Harold. 1961.
The Theory of Probability
. 3rd ed. Oxford.
Kass, Robert E., and Adrian E. Raftery. 1995. “Bayes Factors.”
Journal of the American Statistical Association
90: 773–95.
Fisher, R. 1925.
Statistical Methods for Research Workers
. Edinburgh, UK: Oliver; Boyd.
Johnson, Valen E. 2013. “Revised Standards for Statistical Evidence.”
Proceedings of the National Academy of Sciences
, no. 48: 19313–7.
Gunel, Erdogan, and James Dickey. 1974. “Bayes Factors for Independence in Contingency Tables.”
Biometrika
, 545–57.
Kruschke, J. K. 2011.
Doing Bayesian Data Analysis: A Tutorial with R and BUGS
. Burlington, MA: Academic Press.
Lee, Michael D, and Eric-Jan Wagenmakers. 2014.
Bayesian Cognitive Modeling: A Practical Course
. Cambridge University Press.
http://en.wikiquote.org/wiki/David_Hume
↩
http://en.wikipedia.org/wiki/Climate_of_Adelaide
↩
It’s a leap of faith, I know, but let’s run with it okay?
↩
Um. I hate to bring this up, but some statisticians would object to me using the word “likelihood” here. The problem is that the word “likelihood” has a very specific meaning in frequentist statistics, and it’s not quite the same as what it means in Bayesian statistics. As far as I can tell, Bayesians didn’t originally have any agreed upon name for the likelihood, and so it became common practice for people to use the frequentist terminology. This wouldn’t have been a problem, except for the fact that the way that Bayesians use the word turns out to be quite different to the way frequentists do. This isn’t the place for yet another lengthy history lesson, but to put it crudely: when a Bayesian says “
a
likelihood function” they’re usually referring one of the
rows
of the table. When a frequentist says the same thing, they’re referring to the same table, but to them “
a
likelihood function” almost always refers to one of the
columns
. This distinction matters in some contexts, but it’s not important for our purposes.
↩
If we were being a bit more sophisticated, we could extend the example to accommodate the possibility that I’m lying about the umbrella. But let’s keep things simple, shall we?
↩
You might notice that this equation is actually a restatement of the same basic rule I listed at the start of the last section. If you multiply both sides of the equation by
\(P(d)\)
, then you get
\(P(d) P(h| d) = P(d,h)\)
, which is the rule for how joint probabilities are calculated. So I’m not actually introducing any “new” rules here, I’m just using the same rule in a different way.
↩
Obviously, this is a highly simplified story. All the complexity of real life Bayesian hypothesis testing comes down to how you calculate the likelihood
\(P(d|h)\)
when the hypothesis
\(h\)
is a complex and vague thing. I’m not going to talk about those complexities in this book, but I do want to highlight that although this simple story is true as far as it goes, real life is messier than I’m able to cover in an introductory stats textbook.
↩
http://www.imdb.com/title/tt0093779/quotes
. I should note in passing that I’m not the first person to use this quote to complain about frequentist methods. Rich Morey and colleagues had the idea first. I’m shamelessly stealing it because it’s such an awesome pull quote to use in this context and I refuse to miss any opportunity to quote
The Princess Bride
.
↩
http://about.abc.net.au/reports-publications/appreciation-survey-summary-report-2013/
↩
http://knowyourmeme.com/memes/the-cake-is-a-lie
↩
In the interests of being completely honest, I should acknowledge that not all orthodox statistical tests that rely on this silly assumption. There are a number of
sequential analysis
tools that are sometimes used in clinical trials and the like. These methods are built on the assumption that data are analysed as they arrive, and these tests aren’t horribly broken in the way I’m complaining about here. However, sequential analysis methods are constructed in a very different fashion to the “standard” version of null hypothesis testing. They don’t make it into any introductory textbooks, and they’re not very widely used in the psychological literature. The concern I’m raising here is valid for every single orthodox test I’ve presented so far, and for almost every test I’ve seen reported in the papers I read.
↩
A related problem:
http://xkcd.com/1478/
↩
Some readers might wonder why I picked 3:1 rather than 5:1, given that
Johnson (
2013
)
suggests that
\(p=.05\)
lies somewhere in that range. I did so in order to be charitable to the
\(p\)
-value. If I’d chosen a 5:1 Bayes factor instead, the results would look even better for the Bayesian approach.
↩
http://www.quotationspage.com/quotes/Ambrosius_Macrobius/
↩
Okay, I just
know
that some knowledgeable frequentists will read this and start complaining about this section. Look, I’m not dumb. I absolutely know that if you adopt a sequential analysis perspective you can avoid these errors within the orthodox framework. I also know that you can explictly design studies with interim analyses in mind. So yes, in one sense I’m attacking a “straw man” version of orthodox methods. However, the straw man that I’m attacking is the one that
is used by almost every single practitioner
. If it ever reaches the point where sequential methods become the norm among experimental psychologists and I’m no longer forced to read 20 extremely dubious ANOVAs a day, I promise I’ll rewrite this section and dial down the vitriol. But until that day arrives, I stand by my claim that
default
Bayes factor methods are much more robust in the face of data analysis practices as they exist in the real world.
Default
orthodox methods suck, and we all know it.
↩
If you’re desperate to know, you can find all the gory details in Gunel and Dickey (1974). However, that’s a pretty technical paper. The help documentation to the
contingencyTableBF()
gives this explanation: “the argument
priorConcentration
indexes the expected deviation from the null hypothesis under the alternative, and corresponds to Gunel and Dickey’s (1974)
\(a\)
parameter.” As I write this I’m about halfway through the Gunel and Dickey paper, and I agree that setting
\(a=1\)
is a pretty sensible default choice, since it corresponds to an assumption that you have very little
a priori
knowledge about the contingency table.
↩
In some of the later examples, you’ll see that this number is not always 0%. This is because the
BayesFactor
package often has to run some simulations to compute approximate Bayes factors. So the answers you get won’t always be identical when you run the command a second time. That’s why the output of these functions tells you what the margin for error is.
↩
Apparently this omission is deliberate. I have this vague recollection that I spoke to Jeff Rouder about this once, and his opinion was that when homogeneity of variance is violated the results of a
\(t\)
-test are uninterpretable. I can see the argument for this, but I’ve never really held a strong opinion myself. (Jeff, if you never said that, I’m sorry)
↩
Just in case you’re interested: the “JZS” part of the output relates to how the Bayesian test expresses the prior uncertainty about the variance
\(\sigma^2\)
, and it’s short for the names of three people: “Jeffreys Zellner Siow”. See
Rouder et al. (
2009
)
for details.
↩
Again, in case you care … the null hypothesis here specifies an effect size of 0, since the two means are identical. The alternative hypothesis states that there
is
an effect, but it doesn’t specify exactly how big the effect will be. The
\(r\)
value here relates to how big the effect is expected to be according to the alternative. You can type
?ttestBF
to get more details.
↩
Again, guys, sorry if I’ve misread you.
↩
I don’t even disagree with them: it’s
not
at all obvious why a Bayesian ANOVA should reproduce (say) the same set of model comparisons that the Type II testing strategy uses. It’s precisely because of the fact that I haven’t really come to any strong conclusions that I haven’t added anything to the
lsr
package to make Bayesian Type II tests easier to produce.
↩ |
| Markdown | - [Learning Statistics with R](https://learningstatisticswithr.com/book/)
- [Overview](https://learningstatisticswithr.com/book/index.html)
- [Licensing](https://learningstatisticswithr.com/book/licensing.html)
- [Dedication](https://learningstatisticswithr.com/book/dedication.html)
- [Preface](https://learningstatisticswithr.com/book/preface.html)
- [**0\.1** Preface to Version 0.6.1](https://learningstatisticswithr.com/book/preface.html#preface-to-version-0.6.1)
- [**0\.2** Preface to Version 0.6](https://learningstatisticswithr.com/book/preface.html#preface-to-version-0.6)
- [**0\.3** Preface to Version 0.5](https://learningstatisticswithr.com/book/preface.html#preface-to-version-0.5)
- [**0\.4** Preface to Version 0.4](https://learningstatisticswithr.com/book/preface.html#preface-to-version-0.4)
- [**0\.5** Preface to Version 0.3](https://learningstatisticswithr.com/book/preface.html#preface-to-version-0.3)
- [Part I. Background](https://learningstatisticswithr.com/book/part-i-background.html)
- [**1** Why do we learn statistics?](https://learningstatisticswithr.com/book/why-do-we-learn-statistics.html)
- [**1\.1** On the psychology of statistics](https://learningstatisticswithr.com/book/why-do-we-learn-statistics.html#whywhywhy)
- [**1\.1.1** The curse of belief bias](https://learningstatisticswithr.com/book/why-do-we-learn-statistics.html#the-curse-of-belief-bias)
- [**1\.2** The cautionary tale of Simpson’s paradox](https://learningstatisticswithr.com/book/why-do-we-learn-statistics.html#the-cautionary-tale-of-simpsons-paradox)
- [**1\.3** Statistics in psychology](https://learningstatisticswithr.com/book/why-do-we-learn-statistics.html#statistics-in-psychology)
- [**1\.4** Statistics in everyday life](https://learningstatisticswithr.com/book/why-do-we-learn-statistics.html#statistics-in-everyday-life)
- [**1\.5** There’s more to research methods than statistics](https://learningstatisticswithr.com/book/why-do-we-learn-statistics.html#theres-more-to-research-methods-than-statistics)
- [**2** A brief introduction to research design](https://learningstatisticswithr.com/book/studydesign.html)
- [**2\.1** Introduction to psychological measurement](https://learningstatisticswithr.com/book/studydesign.html#measurement)
- [**2\.1.1** Some thoughts about psychological measurement](https://learningstatisticswithr.com/book/studydesign.html#some-thoughts-about-psychological-measurement)
- [**2\.1.2** Operationalisation: defining your measurement](https://learningstatisticswithr.com/book/studydesign.html#operationalisation-defining-your-measurement)
- [**2\.2** Scales of measurement](https://learningstatisticswithr.com/book/studydesign.html#scales)
- [**2\.2.1** Nominal scale](https://learningstatisticswithr.com/book/studydesign.html#nominal-scale)
- [**2\.2.2** Ordinal scale](https://learningstatisticswithr.com/book/studydesign.html#ordinal-scale)
- [**2\.2.3** Interval scale](https://learningstatisticswithr.com/book/studydesign.html#interval-scale)
- [**2\.2.4** Ratio scale](https://learningstatisticswithr.com/book/studydesign.html#ratio-scale)
- [**2\.2.5** Continuous versus discrete variables](https://learningstatisticswithr.com/book/studydesign.html#continuousdiscrete)
- [**2\.2.6** Some complexities](https://learningstatisticswithr.com/book/studydesign.html#some-complexities)
- [**2\.3** Assessing the reliability of a measurement](https://learningstatisticswithr.com/book/studydesign.html#reliability)
- [**2\.4** The “role” of variables: predictors and outcomes](https://learningstatisticswithr.com/book/studydesign.html#ivdv)
- [**2\.5** Experimental and non-experimental research](https://learningstatisticswithr.com/book/studydesign.html#researchdesigns)
- [**2\.5.1** Experimental research](https://learningstatisticswithr.com/book/studydesign.html#experimental-research)
- [**2\.5.2** Non-experimental research](https://learningstatisticswithr.com/book/studydesign.html#non-experimental-research)
- [**2\.6** Assessing the validity of a study](https://learningstatisticswithr.com/book/studydesign.html#validity)
- [**2\.6.1** Internal validity](https://learningstatisticswithr.com/book/studydesign.html#internal-validity)
- [**2\.6.2** External validity](https://learningstatisticswithr.com/book/studydesign.html#external-validity)
- [**2\.6.3** Construct validity](https://learningstatisticswithr.com/book/studydesign.html#construct-validity)
- [**2\.6.4** Face validity](https://learningstatisticswithr.com/book/studydesign.html#face-validity)
- [**2\.6.5** Ecological validity](https://learningstatisticswithr.com/book/studydesign.html#ecological-validity)
- [**2\.7** Confounds, artifacts and other threats to validity](https://learningstatisticswithr.com/book/studydesign.html#confounds-artifacts-and-other-threats-to-validity)
- [**2\.7.1** History effects](https://learningstatisticswithr.com/book/studydesign.html#history-effects)
- [**2\.7.2** Maturation effects](https://learningstatisticswithr.com/book/studydesign.html#maturation-effects)
- [**2\.7.3** Repeated testing effects](https://learningstatisticswithr.com/book/studydesign.html#repeated-testing-effects)
- [**2\.7.4** Selection bias](https://learningstatisticswithr.com/book/studydesign.html#selection-bias)
- [**2\.7.5** Differential attrition](https://learningstatisticswithr.com/book/studydesign.html#differentialattrition)
- [**2\.7.6** Non-response bias](https://learningstatisticswithr.com/book/studydesign.html#non-response-bias)
- [**2\.7.7** Regression to the mean](https://learningstatisticswithr.com/book/studydesign.html#regression-to-the-mean)
- [**2\.7.8** Experimenter bias](https://learningstatisticswithr.com/book/studydesign.html#experimenter-bias)
- [**2\.7.9** Demand effects and reactivity](https://learningstatisticswithr.com/book/studydesign.html#demand-effects-and-reactivity)
- [**2\.7.10** Placebo effects](https://learningstatisticswithr.com/book/studydesign.html#placebo-effects)
- [**2\.7.11** Situation, measurement and subpopulation effects](https://learningstatisticswithr.com/book/studydesign.html#situation-measurement-and-subpopulation-effects)
- [**2\.7.12** Fraud, deception and self-deception](https://learningstatisticswithr.com/book/studydesign.html#fraud-deception-and-self-deception)
- [**2\.8** Summary](https://learningstatisticswithr.com/book/studydesign.html#summary)
- [Part II. An introduction to R](https://learningstatisticswithr.com/book/part-ii-an-introduction-to-r.html)
- [**3** Getting started with R](https://learningstatisticswithr.com/book/introR.html)
- [**3\.1** Installing R](https://learningstatisticswithr.com/book/introR.html#gettingR)
- [**3\.1.1** Installing R on a Windows computer](https://learningstatisticswithr.com/book/introR.html#installing-r-on-a-windows-computer)
- [**3\.1.2** Installing R on a Mac](https://learningstatisticswithr.com/book/introR.html#installing-r-on-a-mac)
- [**3\.1.3** Installing R on a Linux computer](https://learningstatisticswithr.com/book/introR.html#installing-r-on-a-linux-computer)
- [**3\.1.4** Downloading and installing RStudio](https://learningstatisticswithr.com/book/introR.html#installingrstudio)
- [**3\.1.5** Starting up R](https://learningstatisticswithr.com/book/introR.html#startingR)
- [**3\.2** Typing commands at the R console](https://learningstatisticswithr.com/book/introR.html#firstcommand)
- [**3\.2.1** An important digression about formatting](https://learningstatisticswithr.com/book/introR.html#an-important-digression-about-formatting)
- [**3\.2.2** Be very careful to avoid typos](https://learningstatisticswithr.com/book/introR.html#be-very-careful-to-avoid-typos)
- [**3\.2.3** R is (a bit) flexible with spacing](https://learningstatisticswithr.com/book/introR.html#r-is-a-bit-flexible-with-spacing)
- [**3\.2.4** R can sometimes tell that you’re not finished yet (but not often)](https://learningstatisticswithr.com/book/introR.html#r-can-sometimes-tell-that-youre-not-finished-yet-but-not-often)
- [**3\.3** Doing simple calculations with R](https://learningstatisticswithr.com/book/introR.html#arithmetic)
- [**3\.3.1** Adding, subtracting, multiplying and dividing](https://learningstatisticswithr.com/book/introR.html#adding-subtracting-multiplying-and-dividing)
- [**3\.3.2** Taking powers](https://learningstatisticswithr.com/book/introR.html#taking-powers)
- [**3\.3.3** Doing calculations in the right order](https://learningstatisticswithr.com/book/introR.html#bedmas)
- [**3\.4** Storing a number as a variable](https://learningstatisticswithr.com/book/introR.html#assign)
- [**3\.4.1** Variable assignment using `<-` and `->`](https://learningstatisticswithr.com/book/introR.html#variable-assignment-using---and--)
- [**3\.4.2** Doing calculations using variables](https://learningstatisticswithr.com/book/introR.html#doing-calculations-using-variables)
- [**3\.4.3** Rules and conventions for naming variables](https://learningstatisticswithr.com/book/introR.html#rules-and-conventions-for-naming-variables)
- [**3\.5** Using functions to do calculations](https://learningstatisticswithr.com/book/introR.html#usingfunctions)
- [**3\.5.1** Function arguments, their names and their defaults](https://learningstatisticswithr.com/book/introR.html#functionarguments)
- [**3\.6** Letting RStudio help you with your commands](https://learningstatisticswithr.com/book/introR.html#RStudio1)
- [**3\.6.1** Autocomplete using “tab”](https://learningstatisticswithr.com/book/introR.html#autocomplete-using-tab)
- [**3\.6.2** Browsing your command history](https://learningstatisticswithr.com/book/introR.html#browsing-your-command-history)
- [**3\.7** Storing many numbers as a vector](https://learningstatisticswithr.com/book/introR.html#vectors)
- [**3\.7.1** Creating a vector](https://learningstatisticswithr.com/book/introR.html#creating-a-vector)
- [**3\.7.2** A handy digression](https://learningstatisticswithr.com/book/introR.html#a-handy-digression)
- [**3\.7.3** Getting information out of vectors](https://learningstatisticswithr.com/book/introR.html#vectorsubset)
- [**3\.7.4** Altering the elements of a vector](https://learningstatisticswithr.com/book/introR.html#altering-the-elements-of-a-vector)
- [**3\.7.5** Useful things to know about vectors](https://learningstatisticswithr.com/book/introR.html#veclength)
- [**3\.8** Storing text data](https://learningstatisticswithr.com/book/introR.html#text)
- [**3\.8.1** Working with text](https://learningstatisticswithr.com/book/introR.html#simpletext)
- [**3\.9** Storing “true or false” data](https://learningstatisticswithr.com/book/introR.html#logicals)
- [**3\.9.1** Assessing mathematical truths](https://learningstatisticswithr.com/book/introR.html#assessing-mathematical-truths)
- [**3\.9.2** Logical operations](https://learningstatisticswithr.com/book/introR.html#logical-operations)
- [**3\.9.3** Storing and using logical data](https://learningstatisticswithr.com/book/introR.html#storing-and-using-logical-data)
- [**3\.9.4** Vectors of logicals](https://learningstatisticswithr.com/book/introR.html#vectors-of-logicals)
- [**3\.9.5** Applying logical operation to text](https://learningstatisticswithr.com/book/introR.html#logictext)
- [**3\.10** Indexing vectors](https://learningstatisticswithr.com/book/introR.html#indexing)
- [**3\.10.1** Extracting multiple elements](https://learningstatisticswithr.com/book/introR.html#extracting-multiple-elements)
- [**3\.10.2** Logical indexing](https://learningstatisticswithr.com/book/introR.html#logical-indexing)
- [**3\.11** Quitting R](https://learningstatisticswithr.com/book/introR.html#quitting-r)
- [**3\.12** Summary](https://learningstatisticswithr.com/book/introR.html#summary-1)
- [**4** Additional R concepts](https://learningstatisticswithr.com/book/mechanics.html)
- [**4\.1** Using comments](https://learningstatisticswithr.com/book/mechanics.html#comments)
- [**4\.2** Installing and loading packages](https://learningstatisticswithr.com/book/mechanics.html#packageinstall)
- [**4\.2.1** The package panel in RStudio](https://learningstatisticswithr.com/book/mechanics.html#the-package-panel-in-rstudio)
- [**4\.2.2** Loading a package](https://learningstatisticswithr.com/book/mechanics.html#packageload)
- [**4\.2.3** Unloading a package](https://learningstatisticswithr.com/book/mechanics.html#packageunload)
- [**4\.2.4** A few extra comments](https://learningstatisticswithr.com/book/mechanics.html#a-few-extra-comments)
- [**4\.2.5** Downloading new packages](https://learningstatisticswithr.com/book/mechanics.html#downloading-new-packages)
- [**4\.2.6** Updating R and R packages](https://learningstatisticswithr.com/book/mechanics.html#updating-r-and-r-packages)
- [**4\.2.7** What packages does this book use?](https://learningstatisticswithr.com/book/mechanics.html#what-packages-does-this-book-use)
- [**4\.3** Managing the workspace](https://learningstatisticswithr.com/book/mechanics.html#workspace)
- [**4\.3.1** Listing the contents of the workspace](https://learningstatisticswithr.com/book/mechanics.html#listing-the-contents-of-the-workspace)
- [**4\.3.2** Removing variables from the workspace](https://learningstatisticswithr.com/book/mechanics.html#removing-variables-from-the-workspace)
- [**4\.4** Navigating the file system](https://learningstatisticswithr.com/book/mechanics.html#navigation)
- [**4\.4.1** The file system itself](https://learningstatisticswithr.com/book/mechanics.html#filesystem)
- [**4\.4.2** Navigating the file system using the R console](https://learningstatisticswithr.com/book/mechanics.html#navigationR)
- [**4\.4.3** Why do the Windows paths use the wrong slash?](https://learningstatisticswithr.com/book/mechanics.html#why-do-the-windows-paths-use-the-wrong-slash)
- [**4\.4.4** Navigating the file system using the RStudio file panel](https://learningstatisticswithr.com/book/mechanics.html#nav3)
- [**4\.5** Loading and saving data](https://learningstatisticswithr.com/book/mechanics.html#load)
- [**4\.5.1** Loading workspace files using R](https://learningstatisticswithr.com/book/mechanics.html#loading-workspace-files-using-r)
- [**4\.5.2** Loading workspace files using RStudio](https://learningstatisticswithr.com/book/mechanics.html#loading-workspace-files-using-rstudio)
- [**4\.5.3** Importing data from CSV files using `loadingcsv`](https://learningstatisticswithr.com/book/mechanics.html#importing-data-from-csv-files-using-loadingcsv)
- [**4\.5.4** Importing data from CSV files using RStudio](https://learningstatisticswithr.com/book/mechanics.html#importing-data-from-csv-files-using-rstudio)
- [**4\.5.5** Saving a workspace file using `save`](https://learningstatisticswithr.com/book/mechanics.html#saving-a-workspace-file-using-save)
- [**4\.5.6** Saving a workspace file using RStudio](https://learningstatisticswithr.com/book/mechanics.html#save1)
- [**4\.5.7** Other things you might want to save](https://learningstatisticswithr.com/book/mechanics.html#other-things-you-might-want-to-save)
- [**4\.6** Useful things to know about variables](https://learningstatisticswithr.com/book/mechanics.html#useful)
- [**4\.6.1** Special values](https://learningstatisticswithr.com/book/mechanics.html#specials)
- [**4\.6.2** Assigning names to vector elements](https://learningstatisticswithr.com/book/mechanics.html#names)
- [**4\.6.3** Variable classes](https://learningstatisticswithr.com/book/mechanics.html#variable-classes)
- [**4\.7** Factors](https://learningstatisticswithr.com/book/mechanics.html#factors)
- [**4\.7.1** Introducing factors](https://learningstatisticswithr.com/book/mechanics.html#introducing-factors)
- [**4\.7.2** Labelling the factor levels](https://learningstatisticswithr.com/book/mechanics.html#labelling-the-factor-levels)
- [**4\.7.3** Moving on…](https://learningstatisticswithr.com/book/mechanics.html#moving-on)
- [**4\.8** Data frames](https://learningstatisticswithr.com/book/mechanics.html#dataframes)
- [**4\.8.1** Introducing data frames](https://learningstatisticswithr.com/book/mechanics.html#introducing-data-frames)
- [**4\.8.2** Pulling out the contents of the data frame using `$`](https://learningstatisticswithr.com/book/mechanics.html#pulling-out-the-contents-of-the-data-frame-using)
- [**4\.8.3** Getting information about a data frame](https://learningstatisticswithr.com/book/mechanics.html#getting-information-about-a-data-frame)
- [**4\.8.4** Looking for more on data frames?](https://learningstatisticswithr.com/book/mechanics.html#looking-for-more-on-data-frames)
- [**4\.9** Lists](https://learningstatisticswithr.com/book/mechanics.html#lists)
- [**4\.10** Formulas](https://learningstatisticswithr.com/book/mechanics.html#formulas)
- [**4\.11** Generic functions](https://learningstatisticswithr.com/book/mechanics.html#generics)
- [**4\.12** Getting help](https://learningstatisticswithr.com/book/mechanics.html#help)
- [**4\.12.1** How to read the help documentation](https://learningstatisticswithr.com/book/mechanics.html#how-to-read-the-help-documentation)
- [**4\.12.2** Other resources](https://learningstatisticswithr.com/book/mechanics.html#other-resources)
- [**4\.13** Summary](https://learningstatisticswithr.com/book/mechanics.html#summary-2)
- [Part III. Working with data](https://learningstatisticswithr.com/book/part-iii-working-with-data.html)
- [**5** Descriptive statistics](https://learningstatisticswithr.com/book/descriptives.html)
- [**5\.1** Measures of central tendency](https://learningstatisticswithr.com/book/descriptives.html#centraltendency)
- [**5\.1.1** The mean](https://learningstatisticswithr.com/book/descriptives.html#mean)
- [**5\.1.2** Calculating the mean in R](https://learningstatisticswithr.com/book/descriptives.html#calculating-the-mean-in-r)
- [**5\.1.3** The median](https://learningstatisticswithr.com/book/descriptives.html#median)
- [**5\.1.4** Mean or median? What’s the difference?](https://learningstatisticswithr.com/book/descriptives.html#mean-or-median-whats-the-difference)
- [**5\.1.5** A real life example](https://learningstatisticswithr.com/book/descriptives.html#housingpriceexample)
- [**5\.1.6** Trimmed mean](https://learningstatisticswithr.com/book/descriptives.html#trimmedmean)
- [**5\.1.7** Mode](https://learningstatisticswithr.com/book/descriptives.html#mode)
- [**5\.2** Measures of variability](https://learningstatisticswithr.com/book/descriptives.html#var)
- [**5\.2.1** Range](https://learningstatisticswithr.com/book/descriptives.html#range)
- [**5\.2.2** Interquartile range](https://learningstatisticswithr.com/book/descriptives.html#interquartile-range)
- [**5\.2.3** Mean absolute deviation](https://learningstatisticswithr.com/book/descriptives.html#aad)
- [**5\.2.4** Variance](https://learningstatisticswithr.com/book/descriptives.html#variance)
- [**5\.2.5** Standard deviation](https://learningstatisticswithr.com/book/descriptives.html#sd)
- [**5\.2.6** Median absolute deviation](https://learningstatisticswithr.com/book/descriptives.html#mad)
- [**5\.2.7** Which measure to use?](https://learningstatisticswithr.com/book/descriptives.html#which-measure-to-use)
- [**5\.3** Skew and kurtosis](https://learningstatisticswithr.com/book/descriptives.html#skewandkurtosis)
- [**5\.4** Getting an overall summary of a variable](https://learningstatisticswithr.com/book/studydesign.html#summary)
- [**5\.4.1** “Summarising” a variable](https://learningstatisticswithr.com/book/descriptives.html#summarising-a-variable)
- [**5\.4.2** “Summarising” a data frame](https://learningstatisticswithr.com/book/descriptives.html#summarising-a-data-frame)
- [**5\.4.3** “Describing” a data frame](https://learningstatisticswithr.com/book/descriptives.html#describing-a-data-frame)
- [**5\.5** Descriptive statistics separately for each group](https://learningstatisticswithr.com/book/descriptives.html#groupdescriptives)
- [**5\.6** Standard scores](https://learningstatisticswithr.com/book/descriptives.html#zscore)
- [**5\.7** Correlations](https://learningstatisticswithr.com/book/descriptives.html#correl)
- [**5\.7.1** The data](https://learningstatisticswithr.com/book/descriptives.html#the-data)
- [**5\.7.2** The strength and direction of a relationship](https://learningstatisticswithr.com/book/descriptives.html#the-strength-and-direction-of-a-relationship)
- [**5\.7.3** The correlation coefficient](https://learningstatisticswithr.com/book/descriptives.html#the-correlation-coefficient)
- [**5\.7.4** Calculating correlations in R](https://learningstatisticswithr.com/book/descriptives.html#calculating-correlations-in-r)
- [**5\.7.5** Interpreting a correlation](https://learningstatisticswithr.com/book/descriptives.html#interpretingcorrelations)
- [**5\.7.6** Spearman’s rank correlations](https://learningstatisticswithr.com/book/descriptives.html#spearmans-rank-correlations)
- [**5\.7.7** The `correlate()` function](https://learningstatisticswithr.com/book/descriptives.html#the-correlate-function)
- [**5\.8** Handling missing values](https://learningstatisticswithr.com/book/descriptives.html#missing)
- [**5\.8.1** The single variable case](https://learningstatisticswithr.com/book/descriptives.html#the-single-variable-case)
- [**5\.8.2** Missing values in pairwise calculations](https://learningstatisticswithr.com/book/descriptives.html#missing-values-in-pairwise-calculations)
- [**5\.9** Summary](https://learningstatisticswithr.com/book/descriptives.html#summary-3)
- [**5\.10** Epilogue: Good descriptive statistics are descriptive\!](https://learningstatisticswithr.com/book/descriptives.html#epilogue-good-descriptive-statistics-are-descriptive)
- [**6** Drawing graphs](https://learningstatisticswithr.com/book/graphics.html)
- [**6\.1** An overview of R graphics](https://learningstatisticswithr.com/book/graphics.html#rgraphics)
- [**6\.2** An introduction to plotting](https://learningstatisticswithr.com/book/graphics.html#introplotting)
- [**6\.2.1** A tedious digression](https://learningstatisticswithr.com/book/graphics.html#a-tedious-digression)
- [**6\.2.2** Customising the title and the axis labels](https://learningstatisticswithr.com/book/graphics.html#figtitles)
- [**6\.2.3** Changing the plot type](https://learningstatisticswithr.com/book/graphics.html#changing-the-plot-type)
- [**6\.2.4** Changing other features of the plot](https://learningstatisticswithr.com/book/graphics.html#changing-other-features-of-the-plot)
- [**6\.2.5** Changing the appearance of the axes](https://learningstatisticswithr.com/book/graphics.html#changing-the-appearance-of-the-axes)
- [**6\.2.6** Don’t panic](https://learningstatisticswithr.com/book/graphics.html#dont-panic)
- [**6\.3** Histograms](https://learningstatisticswithr.com/book/graphics.html#hist)
- [**6\.3.1** Visual style of your histogram](https://learningstatisticswithr.com/book/graphics.html#visual-style-of-your-histogram)
- [**6\.4** Stem and leaf plots](https://learningstatisticswithr.com/book/graphics.html#stem)
- [**6\.5** Boxplots](https://learningstatisticswithr.com/book/graphics.html#boxplots)
- [**6\.5.1** Visual style of your boxplot](https://learningstatisticswithr.com/book/graphics.html#visual-style-of-your-boxplot)
- [**6\.5.2** Using box plots to detect outliers](https://learningstatisticswithr.com/book/graphics.html#boxplotoutliers)
- [**6\.5.3** Drawing multiple boxplots](https://learningstatisticswithr.com/book/graphics.html#multipleboxplots)
- [**6\.6** Scatterplots](https://learningstatisticswithr.com/book/graphics.html#scatterplots)
- [**6\.6.1** More elaborate options](https://learningstatisticswithr.com/book/graphics.html#more-elaborate-options)
- [**6\.7** Bar graphs](https://learningstatisticswithr.com/book/graphics.html#bargraph)
- [**6\.7.1** Changing global settings using par()](https://learningstatisticswithr.com/book/graphics.html#par)
- [**6\.8** Saving image files using R and Rstudio](https://learningstatisticswithr.com/book/graphics.html#saveimage)
- [**6\.8.1** The ugly details (advanced)](https://learningstatisticswithr.com/book/graphics.html#the-ugly-details-advanced)
- [**6\.9** Summary](https://learningstatisticswithr.com/book/graphics.html#summary-4)
- [**7** Pragmatic matters](https://learningstatisticswithr.com/book/datahandling.html)
- [**7\.1** Tabulating and cross-tabulating data](https://learningstatisticswithr.com/book/datahandling.html#freqtables)
- [**7\.1.1** Creating tables from vectors](https://learningstatisticswithr.com/book/datahandling.html#creating-tables-from-vectors)
- [**7\.1.2** Creating tables from data frames](https://learningstatisticswithr.com/book/datahandling.html#creating-tables-from-data-frames)
- [**7\.1.3** Converting a table of counts to a table of proportions](https://learningstatisticswithr.com/book/datahandling.html#converting-a-table-of-counts-to-a-table-of-proportions)
- [**7\.1.4** Low level tabulation](https://learningstatisticswithr.com/book/datahandling.html#low-level-tabulation)
- [**7\.2** Transforming and recoding a variable](https://learningstatisticswithr.com/book/datahandling.html#transform)
- [**7\.2.1** Creating a transformed variable](https://learningstatisticswithr.com/book/datahandling.html#creating-a-transformed-variable)
- [**7\.2.2** Cutting a numeric variable into categories](https://learningstatisticswithr.com/book/datahandling.html#cutting-a-numeric-variable-into-categories)
- [**7\.3** A few more mathematical functions and operations](https://learningstatisticswithr.com/book/datahandling.html#mathfunc)
- [**7\.3.1** Rounding a number](https://learningstatisticswithr.com/book/datahandling.html#rounding-a-number)
- [**7\.3.2** Modulus and integer division](https://learningstatisticswithr.com/book/datahandling.html#modulus-and-integer-division)
- [**7\.3.3** Logarithms and exponentials](https://learningstatisticswithr.com/book/datahandling.html#logarithms-and-exponentials)
- [**7\.4** Extracting a subset of a vector](https://learningstatisticswithr.com/book/datahandling.html#subset)
- [**7\.4.1** Refresher](https://learningstatisticswithr.com/book/datahandling.html#refresher)
- [**7\.4.2** Using `%in%` to match multiple cases](https://learningstatisticswithr.com/book/datahandling.html#using-in-to-match-multiple-cases)
- [**7\.4.3** Using negative indices to drop elements](https://learningstatisticswithr.com/book/datahandling.html#using-negative-indices-to-drop-elements)
- [**7\.4.4** Splitting a vector by group](https://learningstatisticswithr.com/book/datahandling.html#splitting-a-vector-by-group)
- [**7\.5** Extracting a subset of a data frame](https://learningstatisticswithr.com/book/datahandling.html#subsetdataframe)
- [**7\.5.1** Using the `subset()` function](https://learningstatisticswithr.com/book/datahandling.html#using-the-subset-function)
- [**7\.5.2** Using square brackets: I. Rows and columns](https://learningstatisticswithr.com/book/datahandling.html#using-square-brackets-i.-rows-and-columns)
- [**7\.5.3** Using square brackets: II. Some elaborations](https://learningstatisticswithr.com/book/datahandling.html#using-square-brackets-ii.-some-elaborations)
- [**7\.5.4** Using square brackets: III. Understanding “dropping”](https://learningstatisticswithr.com/book/datahandling.html#dropping)
- [**7\.5.5** Using square brackets: IV. Columns only](https://learningstatisticswithr.com/book/datahandling.html#using-square-brackets-iv.-columns-only)
- [**7\.6** Sorting, flipping and merging data](https://learningstatisticswithr.com/book/datahandling.html#sort)
- [**7\.6.1** Sorting a numeric or character vector](https://learningstatisticswithr.com/book/datahandling.html#sorting-a-numeric-or-character-vector)
- [**7\.6.2** Sorting a factor](https://learningstatisticswithr.com/book/datahandling.html#sorting-a-factor)
- [**7\.6.3** Sorting a data frame](https://learningstatisticswithr.com/book/datahandling.html#sortframe)
- [**7\.6.4** Binding vectors together](https://learningstatisticswithr.com/book/datahandling.html#binding-vectors-together)
- [**7\.6.5** Binding multiple copies of the same vector together](https://learningstatisticswithr.com/book/datahandling.html#binding-multiple-copies-of-the-same-vector-together)
- [**7\.6.6** Transposing a matrix or data frame](https://learningstatisticswithr.com/book/datahandling.html#transposing-a-matrix-or-data-frame)
- [**7\.7** Reshaping a data frame](https://learningstatisticswithr.com/book/datahandling.html#reshape)
- [**7\.7.1** Long form and wide form data](https://learningstatisticswithr.com/book/datahandling.html#long-form-and-wide-form-data)
- [**7\.7.2** Reshaping data using `wideToLong()`](https://learningstatisticswithr.com/book/datahandling.html#reshaping-data-using-widetolong)
- [**7\.7.3** Reshaping data using `longToWide()`](https://learningstatisticswithr.com/book/datahandling.html#reshaping-data-using-longtowide)
- [**7\.7.4** Reshaping with multiple within-subject factors](https://learningstatisticswithr.com/book/datahandling.html#reshaping-with-multiple-within-subject-factors)
- [**7\.7.5** What other options are there?](https://learningstatisticswithr.com/book/datahandling.html#what-other-options-are-there)
- [**7\.8** Working with text](https://learningstatisticswithr.com/book/datahandling.html#textprocessing)
- [**7\.8.1** Shortening a string](https://learningstatisticswithr.com/book/datahandling.html#shortening-a-string)
- [**7\.8.2** Pasting strings together](https://learningstatisticswithr.com/book/datahandling.html#pasting-strings-together)
- [**7\.8.3** Splitting strings](https://learningstatisticswithr.com/book/datahandling.html#splitting-strings)
- [**7\.8.4** Making simple conversions](https://learningstatisticswithr.com/book/datahandling.html#making-simple-conversions)
- [**7\.8.5** Applying logical operations to text](https://learningstatisticswithr.com/book/datahandling.html#logictext2)
- [**7\.8.6** Concatenating and printing with `cat()`](https://learningstatisticswithr.com/book/datahandling.html#concatenating-and-printing-with-cat)
- [**7\.8.7** Using escape characters in text](https://learningstatisticswithr.com/book/datahandling.html#escapechars)
- [**7\.8.8** Matching and substituting text](https://learningstatisticswithr.com/book/datahandling.html#matching-and-substituting-text)
- [**7\.8.9** Regular expressions (not really)](https://learningstatisticswithr.com/book/datahandling.html#regex)
- [**7\.9** Reading unusual data files](https://learningstatisticswithr.com/book/datahandling.html#importing)
- [**7\.9.1** Loading data from text files](https://learningstatisticswithr.com/book/datahandling.html#loading-data-from-text-files)
- [**7\.9.2** Loading data from SPSS (and other statistics packages)](https://learningstatisticswithr.com/book/datahandling.html#loading-data-from-spss-and-other-statistics-packages)
- [**7\.9.3** Loading Excel files](https://learningstatisticswithr.com/book/datahandling.html#loading-excel-files)
- [**7\.9.4** Loading Matlab (& Octave) files](https://learningstatisticswithr.com/book/datahandling.html#loading-matlab-octave-files)
- [**7\.9.5** Saving other kinds of data](https://learningstatisticswithr.com/book/datahandling.html#saving-other-kinds-of-data)
- [**7\.9.6** Are we done yet?](https://learningstatisticswithr.com/book/datahandling.html#are-we-done-yet)
- [**7\.10** Coercing data from one class to another](https://learningstatisticswithr.com/book/datahandling.html#coercion)
- [**7\.11** Other useful data structures](https://learningstatisticswithr.com/book/datahandling.html#datastructures)
- [**7\.11.1** Matrices](https://learningstatisticswithr.com/book/datahandling.html#matrix)
- [**7\.11.2** Ordered factors](https://learningstatisticswithr.com/book/datahandling.html#orderedfactors)
- [**7\.11.3** Dates and times](https://learningstatisticswithr.com/book/datahandling.html#dates)
- [**7\.12** Miscellaneous topics](https://learningstatisticswithr.com/book/datahandling.html#miscdatahandling)
- [**7\.12.1** The problems with floating point arithmetic](https://learningstatisticswithr.com/book/datahandling.html#the-problems-with-floating-point-arithmetic)
- [**7\.12.2** The recycling rule](https://learningstatisticswithr.com/book/datahandling.html#recycling)
- [**7\.12.3** An introduction to environments](https://learningstatisticswithr.com/book/datahandling.html#environments)
- [**7\.12.4** Attaching a data frame](https://learningstatisticswithr.com/book/datahandling.html#attaching-a-data-frame)
- [**7\.13** Summary](https://learningstatisticswithr.com/book/datahandling.html#summary-5)
- [**8** Basic programming](https://learningstatisticswithr.com/book/scripting.html)
- [**8\.1** Scripts](https://learningstatisticswithr.com/book/scripting.html#scripts)
- [**8\.1.1** Why use scripts?](https://learningstatisticswithr.com/book/scripting.html#why-use-scripts)
- [**8\.1.2** Our first script](https://learningstatisticswithr.com/book/scripting.html#our-first-script)
- [**8\.1.3** Using Rstudio to write scripts](https://learningstatisticswithr.com/book/scripting.html#using-rstudio-to-write-scripts)
- [**8\.1.4** Commenting your script](https://learningstatisticswithr.com/book/scripting.html#commenting-your-script)
- [**8\.1.5** Differences between scripts and the command line](https://learningstatisticswithr.com/book/scripting.html#differences-between-scripts-and-the-command-line)
- [**8\.1.6** Done\!](https://learningstatisticswithr.com/book/scripting.html#done)
- [**8\.2** Loops](https://learningstatisticswithr.com/book/scripting.html#loops)
- [**8\.2.1** The `while` loop](https://learningstatisticswithr.com/book/scripting.html#the-while-loop)
- [**8\.2.2** The `for` loop](https://learningstatisticswithr.com/book/scripting.html#for)
- [**8\.2.3** A more realistic example of a loop](https://learningstatisticswithr.com/book/scripting.html#a-more-realistic-example-of-a-loop)
- [**8\.3** Conditional statements](https://learningstatisticswithr.com/book/scripting.html#if)
- [**8\.4** Writing functions](https://learningstatisticswithr.com/book/scripting.html#functions)
- [**8\.4.1** Function arguments revisited](https://learningstatisticswithr.com/book/scripting.html#dotsargument)
- [**8\.4.2** There’s more to functions than this](https://learningstatisticswithr.com/book/scripting.html#theres-more-to-functions-than-this)
- [**8\.5** Implicit loops](https://learningstatisticswithr.com/book/scripting.html#vectorised)
- [**8\.6** Summary](https://learningstatisticswithr.com/book/scripting.html#summary-6)
- [Part IV. Statistical theory](https://learningstatisticswithr.com/book/part-iv-statistical-theory.html)
- [On the limits of logical reasoning](https://learningstatisticswithr.com/book/part-iv-statistical-theory.html#on-the-limits-of-logical-reasoning)
- [Learning without making assumptions is a myth](https://learningstatisticswithr.com/book/part-iv-statistical-theory.html#learning-without-making-assumptions-is-a-myth)
- [**9** Introduction to probability](https://learningstatisticswithr.com/book/probability.html)
- [**9\.1** How are probability and statistics different?](https://learningstatisticswithr.com/book/probability.html#probstats)
- [**9\.2** What does probability mean?](https://learningstatisticswithr.com/book/probability.html#probmeaning)
- [**9\.2.1** The frequentist view](https://learningstatisticswithr.com/book/probability.html#the-frequentist-view)
- [**9\.2.2** The Bayesian view](https://learningstatisticswithr.com/book/probability.html#the-bayesian-view)
- [**9\.2.3** What’s the difference? And who is right?](https://learningstatisticswithr.com/book/probability.html#whats-the-difference-and-who-is-right)
- [**9\.3** Basic probability theory](https://learningstatisticswithr.com/book/probability.html#basicprobability)
- [**9\.3.1** Introducing probability distributions](https://learningstatisticswithr.com/book/probability.html#introducing-probability-distributions)
- [**9\.4** The binomial distribution](https://learningstatisticswithr.com/book/probability.html#binomial)
- [**9\.4.1** Introducing the binomial](https://learningstatisticswithr.com/book/probability.html#introducing-the-binomial)
- [**9\.4.2** Working with the binomial distribution in R](https://learningstatisticswithr.com/book/probability.html#working-with-the-binomial-distribution-in-r)
- [**9\.5** The normal distribution](https://learningstatisticswithr.com/book/probability.html#normal)
- [**9\.5.1** Probability density](https://learningstatisticswithr.com/book/probability.html#density)
- [**9\.6** Other useful distributions](https://learningstatisticswithr.com/book/probability.html#otherdists)
- [**9\.7** Summary](https://learningstatisticswithr.com/book/probability.html#summary-7)
- [**10** Estimating unknown quantities from a sample](https://learningstatisticswithr.com/book/estimation.html)
- [**10\.1** Samples, populations and sampling](https://learningstatisticswithr.com/book/estimation.html#srs)
- [**10\.1.1** Defining a population](https://learningstatisticswithr.com/book/estimation.html#pop)
- [**10\.1.2** Simple random samples](https://learningstatisticswithr.com/book/estimation.html#simple-random-samples)
- [**10\.1.3** Most samples are not simple random samples](https://learningstatisticswithr.com/book/estimation.html#most-samples-are-not-simple-random-samples)
- [**10\.1.4** How much does it matter if you don’t have a simple random sample?](https://learningstatisticswithr.com/book/estimation.html#how-much-does-it-matter-if-you-dont-have-a-simple-random-sample)
- [**10\.1.5** Population parameters and sample statistics](https://learningstatisticswithr.com/book/estimation.html#population-parameters-and-sample-statistics)
- [**10\.2** The law of large numbers](https://learningstatisticswithr.com/book/estimation.html#lawlargenumbers)
- [**10\.3** Sampling distributions and the central limit theorem](https://learningstatisticswithr.com/book/estimation.html#samplesandclt)
- [**10\.3.1** Sampling distribution of the mean](https://learningstatisticswithr.com/book/estimation.html#samplingdists)
- [**10\.3.2** Sampling distributions exist for any sample statistic\!](https://learningstatisticswithr.com/book/estimation.html#sampling-distributions-exist-for-any-sample-statistic)
- [**10\.3.3** The central limit theorem](https://learningstatisticswithr.com/book/estimation.html#clt)
- [**10\.4** Estimating population parameters](https://learningstatisticswithr.com/book/estimation.html#pointestimates)
- [**10\.4.1** Estimating the population mean](https://learningstatisticswithr.com/book/estimation.html#estimating-the-population-mean)
- [**10\.4.2** Estimating the population standard deviation](https://learningstatisticswithr.com/book/estimation.html#estimating-the-population-standard-deviation)
- [**10\.5** Estimating a confidence interval](https://learningstatisticswithr.com/book/estimation.html#ci)
- [**10\.5.1** A slight mistake in the formula](https://learningstatisticswithr.com/book/estimation.html#a-slight-mistake-in-the-formula)
- [**10\.5.2** Interpreting a confidence interval](https://learningstatisticswithr.com/book/estimation.html#interpreting-a-confidence-interval)
- [**10\.5.3** Calculating confidence intervals in R](https://learningstatisticswithr.com/book/estimation.html#calculating-confidence-intervals-in-r)
- [**10\.5.4** Plotting confidence intervals in R](https://learningstatisticswithr.com/book/estimation.html#ciplots)
- [**10\.6** Summary](https://learningstatisticswithr.com/book/estimation.html#summary-8)
- [**11** Hypothesis testing](https://learningstatisticswithr.com/book/hypothesistesting.html)
- [**11\.1** A menagerie of hypotheses](https://learningstatisticswithr.com/book/hypothesistesting.html#hypotheses)
- [**11\.1.1** Research hypotheses versus statistical hypotheses](https://learningstatisticswithr.com/book/hypothesistesting.html#research-hypotheses-versus-statistical-hypotheses)
- [**11\.1.2** Null hypotheses and alternative hypotheses](https://learningstatisticswithr.com/book/hypothesistesting.html#null-hypotheses-and-alternative-hypotheses)
- [**11\.2** Two types of errors](https://learningstatisticswithr.com/book/hypothesistesting.html#errortypes)
- [**11\.3** Test statistics and sampling distributions](https://learningstatisticswithr.com/book/hypothesistesting.html#teststatistics)
- [**11\.4** Making decisions](https://learningstatisticswithr.com/book/hypothesistesting.html#decisionmaking)
- [**11\.4.1** Critical regions and critical values](https://learningstatisticswithr.com/book/hypothesistesting.html#critical-regions-and-critical-values)
- [**11\.4.2** A note on statistical “significance”](https://learningstatisticswithr.com/book/hypothesistesting.html#a-note-on-statistical-significance)
- [**11\.4.3** The difference between one sided and two sided tests](https://learningstatisticswithr.com/book/hypothesistesting.html#onesidedtests)
- [**11\.5** The \\(p\\) value of a test](https://learningstatisticswithr.com/book/hypothesistesting.html#pvalue)
- [**11\.5.1** A softer view of decision making](https://learningstatisticswithr.com/book/hypothesistesting.html#a-softer-view-of-decision-making)
- [**11\.5.2** The probability of extreme data](https://learningstatisticswithr.com/book/hypothesistesting.html#the-probability-of-extreme-data)
- [**11\.5.3** A common mistake](https://learningstatisticswithr.com/book/hypothesistesting.html#a-common-mistake)
- [**11\.6** Reporting the results of a hypothesis test](https://learningstatisticswithr.com/book/hypothesistesting.html#writeup)
- [**11\.6.1** The issue](https://learningstatisticswithr.com/book/hypothesistesting.html#the-issue)
- [**11\.6.2** Two proposed solutions](https://learningstatisticswithr.com/book/hypothesistesting.html#two-proposed-solutions)
- [**11\.7** Running the hypothesis test in practice](https://learningstatisticswithr.com/book/hypothesistesting.html#running-the-hypothesis-test-in-practice)
- [**11\.8** Effect size, sample size and power](https://learningstatisticswithr.com/book/hypothesistesting.html#effectsize)
- [**11\.8.1** The power function](https://learningstatisticswithr.com/book/hypothesistesting.html#the-power-function)
- [**11\.8.2** Effect size](https://learningstatisticswithr.com/book/hypothesistesting.html#effect-size)
- [**11\.8.3** Increasing the power of your study](https://learningstatisticswithr.com/book/hypothesistesting.html#increasing-the-power-of-your-study)
- [**11\.9** Some issues to consider](https://learningstatisticswithr.com/book/hypothesistesting.html#nhstmess)
- [**11\.9.1** Neyman versus Fisher](https://learningstatisticswithr.com/book/hypothesistesting.html#neyman-versus-fisher)
- [**11\.9.2** Bayesians versus frequentists](https://learningstatisticswithr.com/book/hypothesistesting.html#bayesians-versus-frequentists)
- [**11\.9.3** Traps](https://learningstatisticswithr.com/book/hypothesistesting.html#traps)
- [**11\.10** Summary](https://learningstatisticswithr.com/book/hypothesistesting.html#summary-9)
- [Part V. Statistical tools](https://learningstatisticswithr.com/book/part-v-statistical-tools.html)
- [**12** Categorical data analysis](https://learningstatisticswithr.com/book/chisquare.html)
- [**12\.1** The \\(\\chi^2\\) goodness-of-fit test](https://learningstatisticswithr.com/book/chisquare.html#goftest)
- [**12\.1.1** The cards data](https://learningstatisticswithr.com/book/chisquare.html#the-cards-data)
- [**12\.1.2** The null hypothesis and the alternative hypothesis](https://learningstatisticswithr.com/book/chisquare.html#the-null-hypothesis-and-the-alternative-hypothesis)
- [**12\.1.3** The “goodness of fit” test statistic](https://learningstatisticswithr.com/book/chisquare.html#the-goodness-of-fit-test-statistic)
- [**12\.1.4** The sampling distribution of the GOF statistic (advanced)](https://learningstatisticswithr.com/book/chisquare.html#the-sampling-distribution-of-the-gof-statistic-advanced)
- [**12\.1.5** Degrees of freedom](https://learningstatisticswithr.com/book/chisquare.html#degrees-of-freedom)
- [**12\.1.6** Testing the null hypothesis](https://learningstatisticswithr.com/book/chisquare.html#testing-the-null-hypothesis)
- [**12\.1.7** Doing the test in R](https://learningstatisticswithr.com/book/chisquare.html#gofTestInR)
- [**12\.1.8** Specifying a different null hypothesis](https://learningstatisticswithr.com/book/chisquare.html#specifying-a-different-null-hypothesis)
- [**12\.1.9** How to report the results of the test](https://learningstatisticswithr.com/book/chisquare.html#chisqreport)
- [**12\.1.10** A comment on statistical notation (advanced)](https://learningstatisticswithr.com/book/chisquare.html#a-comment-on-statistical-notation-advanced)
- [**12\.2** The \\(\\chi^2\\) test of independence (or association)](https://learningstatisticswithr.com/book/chisquare.html#chisqindependence)
- [**12\.2.1** Constructing our hypothesis test](https://learningstatisticswithr.com/book/chisquare.html#constructing-our-hypothesis-test)
- [**12\.2.2** Doing the test in R](https://learningstatisticswithr.com/book/chisquare.html#AssocTestInR)
- [**12\.2.3** Postscript](https://learningstatisticswithr.com/book/chisquare.html#postscript)
- [**12\.3** The continuity correction](https://learningstatisticswithr.com/book/chisquare.html#yates)
- [**12\.4** Effect size](https://learningstatisticswithr.com/book/chisquare.html#chisqeffectsize)
- [**12\.5** Assumptions of the test(s)](https://learningstatisticswithr.com/book/chisquare.html#chisqassumptions)
- [**12\.6** The most typical way to do chi-square tests in R](https://learningstatisticswithr.com/book/chisquare.html#chisq.test)
- [**12\.7** The Fisher exact test](https://learningstatisticswithr.com/book/chisquare.html#fisherexacttest)
- [**12\.8** The McNemar test](https://learningstatisticswithr.com/book/chisquare.html#mcnemar)
- [**12\.8.1** Doing the McNemar test in R](https://learningstatisticswithr.com/book/chisquare.html#doing-the-mcnemar-test-in-r)
- [**12\.9** What’s the difference between McNemar and independence?](https://learningstatisticswithr.com/book/chisquare.html#whats-the-difference-between-mcnemar-and-independence)
- [**12\.10** Summary](https://learningstatisticswithr.com/book/chisquare.html#summary-10)
- [**13** Comparing two means](https://learningstatisticswithr.com/book/ttest.html)
- [**13\.1** The one-sample \\(z\\)\-test](https://learningstatisticswithr.com/book/ttest.html#the-one-sample-z-test)
- [**13\.1.1** The inference problem that the test addresses](https://learningstatisticswithr.com/book/ttest.html#the-inference-problem-that-the-test-addresses)
- [**13\.1.2** Constructing the hypothesis test](https://learningstatisticswithr.com/book/ttest.html#constructing-the-hypothesis-test)
- [**13\.1.3** A worked example using R](https://learningstatisticswithr.com/book/ttest.html#a-worked-example-using-r)
- [**13\.1.4** Assumptions of the \\(z\\)\-test](https://learningstatisticswithr.com/book/ttest.html#zassumptions)
- [**13\.2** The one-sample \\(t\\)\-test](https://learningstatisticswithr.com/book/ttest.html#onesamplettest)
- [**13\.2.1** Introducing the \\(t\\)\-test](https://learningstatisticswithr.com/book/ttest.html#introducing-the-t-test)
- [**13\.2.2** Doing the test in R](https://learningstatisticswithr.com/book/ttest.html#doing-the-test-in-r)
- [**13\.2.3** Assumptions of the one sample \\(t\\)\-test](https://learningstatisticswithr.com/book/ttest.html#ttestoneassumptions)
- [**13\.3** The independent samples \\(t\\)\-test (Student test)](https://learningstatisticswithr.com/book/ttest.html#studentttest)
- [**13\.3.1** The data](https://learningstatisticswithr.com/book/ttest.html#the-data-1)
- [**13\.3.2** Introducing the test](https://learningstatisticswithr.com/book/ttest.html#introducing-the-test)
- [**13\.3.3** A “pooled estimate” of the standard deviation](https://learningstatisticswithr.com/book/ttest.html#a-pooled-estimate-of-the-standard-deviation)
- [**13\.3.4** The same pooled estimate, described differently](https://learningstatisticswithr.com/book/ttest.html#the-same-pooled-estimate-described-differently)
- [**13\.3.5** Completing the test](https://learningstatisticswithr.com/book/ttest.html#completing-the-test)
- [**13\.3.6** Doing the test in R](https://learningstatisticswithr.com/book/ttest.html#doing-the-test-in-r-1)
- [**13\.3.7** Positive and negative \\(t\\) values](https://learningstatisticswithr.com/book/ttest.html#positive-and-negative-t-values)
- [**13\.3.8** Assumptions of the test](https://learningstatisticswithr.com/book/ttest.html#studentassumptions)
- [**13\.4** The independent samples \\(t\\)\-test (Welch test)](https://learningstatisticswithr.com/book/ttest.html#welchttest)
- [**13\.4.1** Doing the test in R](https://learningstatisticswithr.com/book/ttest.html#doing-the-test-in-r-2)
- [**13\.4.2** Assumptions of the test](https://learningstatisticswithr.com/book/ttest.html#assumptions-of-the-test)
- [**13\.5** The paired-samples \\(t\\)\-test](https://learningstatisticswithr.com/book/ttest.html#pairedsamplesttest)
- [**13\.5.1** The data](https://learningstatisticswithr.com/book/ttest.html#the-data-2)
- [**13\.5.2** What is the paired samples \\(t\\)\-test?](https://learningstatisticswithr.com/book/ttest.html#what-is-the-paired-samples-t-test)
- [**13\.5.3** Doing the test in R, part 1](https://learningstatisticswithr.com/book/ttest.html#doing-the-test-in-r-part-1)
- [**13\.5.4** Doing the test in R, part 2](https://learningstatisticswithr.com/book/ttest.html#doing-the-test-in-r-part-2)
- [**13\.6** One sided tests](https://learningstatisticswithr.com/book/ttest.html#one-sided-tests)
- [**13\.7** Using the t.test() function](https://learningstatisticswithr.com/book/ttest.html#ttestfunction)
- [**13\.8** Effect size](https://learningstatisticswithr.com/book/ttest.html#cohensd)
- [**13\.8.1** Cohen’s \\(d\\) from one sample](https://learningstatisticswithr.com/book/ttest.html#cohens-d-from-one-sample)
- [**13\.8.2** Cohen’s \\(d\\) from a Student \\(t\\) test](https://learningstatisticswithr.com/book/ttest.html#cohens-d-from-a-student-t-test)
- [**13\.8.3** Cohen’s \\(d\\) from a Welch test](https://learningstatisticswithr.com/book/ttest.html#cohens-d-from-a-welch-test)
- [**13\.8.4** Cohen’s \\(d\\) from a paired-samples test](https://learningstatisticswithr.com/book/ttest.html#cohens-d-from-a-paired-samples-test)
- [**13\.9** Checking the normality of a sample](https://learningstatisticswithr.com/book/ttest.html#shapiro)
- [**13\.9.1** QQ plots](https://learningstatisticswithr.com/book/ttest.html#qq-plots)
- [**13\.9.2** Shapiro-Wilk tests](https://learningstatisticswithr.com/book/ttest.html#shapiro-wilk-tests)
- [**13\.10** Testing non-normal data with Wilcoxon tests](https://learningstatisticswithr.com/book/ttest.html#wilcox)
- [**13\.10.1** Two sample Wilcoxon test](https://learningstatisticswithr.com/book/ttest.html#two-sample-wilcoxon-test)
- [**13\.10.2** One sample Wilcoxon test](https://learningstatisticswithr.com/book/ttest.html#one-sample-wilcoxon-test)
- [**13\.11** Summary](https://learningstatisticswithr.com/book/ttest.html#summary-11)
- [**14** Comparing several means (one-way ANOVA)](https://learningstatisticswithr.com/book/anova.html)
- [**14\.1** An illustrative data set](https://learningstatisticswithr.com/book/anova.html#anxifree)
- [**14\.2** How ANOVA works](https://learningstatisticswithr.com/book/anova.html#anovaintro)
- [**14\.2.1** Two formulas for the variance of \\(Y\\)](https://learningstatisticswithr.com/book/anova.html#two-formulas-for-the-variance-of-y)
- [**14\.2.2** From variances to sums of squares](https://learningstatisticswithr.com/book/anova.html#from-variances-to-sums-of-squares)
- [**14\.2.3** From sums of squares to the \\(F\\)\-test](https://learningstatisticswithr.com/book/anova.html#from-sums-of-squares-to-the-f-test)
- [**14\.2.4** The model for the data and the meaning of \\(F\\) (advanced)](https://learningstatisticswithr.com/book/anova.html#anovamodel)
- [**14\.2.5** A worked example](https://learningstatisticswithr.com/book/anova.html#anovacalc)
- [**14\.3** Running an ANOVA in R](https://learningstatisticswithr.com/book/anova.html#introduceaov)
- [**14\.3.1** Using the `aov()` function to specify your ANOVA](https://learningstatisticswithr.com/book/anova.html#using-the-aov-function-to-specify-your-anova)
- [**14\.3.2** Understanding what the `aov()` function produces](https://learningstatisticswithr.com/book/anova.html#aovobjects)
- [**14\.3.3** Running the hypothesis tests for the ANOVA](https://learningstatisticswithr.com/book/anova.html#running-the-hypothesis-tests-for-the-anova)
- [**14\.4** Effect size](https://learningstatisticswithr.com/book/anova.html#etasquared)
- [**14\.5** Multiple comparisons and post hoc tests](https://learningstatisticswithr.com/book/anova.html#posthoc)
- [**14\.5.1** Running “pairwise” \\(t\\)\-tests](https://learningstatisticswithr.com/book/anova.html#running-pairwise-t-tests)
- [**14\.5.2** Corrections for multiple testing](https://learningstatisticswithr.com/book/anova.html#corrections-for-multiple-testing)
- [**14\.5.3** Bonferroni corrections](https://learningstatisticswithr.com/book/anova.html#bonferroni-corrections)
- [**14\.5.4** Holm corrections](https://learningstatisticswithr.com/book/anova.html#holm-corrections)
- [**14\.5.5** Writing up the post hoc test](https://learningstatisticswithr.com/book/anova.html#writing-up-the-post-hoc-test)
- [**14\.6** Assumptions of one-way ANOVA](https://learningstatisticswithr.com/book/anova.html#anovaassumptions)
- [**14\.6.1** How robust is ANOVA?](https://learningstatisticswithr.com/book/anova.html#how-robust-is-anova)
- [**14\.7** Checking the homogeneity of variance assumption](https://learningstatisticswithr.com/book/anova.html#levene)
- [**14\.7.1** Running the Levene’s test in R](https://learningstatisticswithr.com/book/anova.html#running-the-levenes-test-in-r)
- [**14\.7.2** Additional comments](https://learningstatisticswithr.com/book/anova.html#additional-comments)
- [**14\.8** Removing the homogeneity of variance assumption](https://learningstatisticswithr.com/book/anova.html#welchoneway)
- [**14\.9** Checking the normality assumption](https://learningstatisticswithr.com/book/anova.html#anovanormality)
- [**14\.10** Removing the normality assumption](https://learningstatisticswithr.com/book/anova.html#kruskalwallis)
- [**14\.10.1** The logic behind the Kruskal-Wallis test](https://learningstatisticswithr.com/book/anova.html#the-logic-behind-the-kruskal-wallis-test)
- [**14\.10.2** Additional details](https://learningstatisticswithr.com/book/anova.html#additional-details)
- [**14\.10.3** How to run the Kruskal-Wallis test in R](https://learningstatisticswithr.com/book/anova.html#how-to-run-the-kruskal-wallis-test-in-r)
- [**14\.11** On the relationship between ANOVA and the Student \\(t\\) test](https://learningstatisticswithr.com/book/anova.html#anovaandt)
- [**14\.12** Summary](https://learningstatisticswithr.com/book/anova.html#summary-12)
- [**15** Linear regression](https://learningstatisticswithr.com/book/regression.html)
- [**15\.1** What is a linear regression model?](https://learningstatisticswithr.com/book/regression.html#introregression)
- [**15\.2** Estimating a linear regression model](https://learningstatisticswithr.com/book/regression.html#regressionestimation)
- [**15\.2.1** Using the `lm()` function](https://learningstatisticswithr.com/book/regression.html#lm)
- [**15\.2.2** Interpreting the estimated model](https://learningstatisticswithr.com/book/regression.html#interpreting-the-estimated-model)
- [**15\.3** Multiple linear regression](https://learningstatisticswithr.com/book/regression.html#multipleregression)
- [**15\.3.1** Doing it in R](https://learningstatisticswithr.com/book/regression.html#doing-it-in-r)
- [**15\.3.2** Formula for the general case](https://learningstatisticswithr.com/book/regression.html#formula-for-the-general-case)
- [**15\.4** Quantifying the fit of the regression model](https://learningstatisticswithr.com/book/regression.html#r2)
- [**15\.4.1** The \\(R^2\\) value](https://learningstatisticswithr.com/book/regression.html#the-r2-value)
- [**15\.4.2** The relationship between regression and correlation](https://learningstatisticswithr.com/book/regression.html#the-relationship-between-regression-and-correlation)
- [**15\.4.3** The adjusted \\(R^2\\) value](https://learningstatisticswithr.com/book/regression.html#the-adjusted-r2-value)
- [**15\.5** Hypothesis tests for regression models](https://learningstatisticswithr.com/book/regression.html#regressiontests)
- [**15\.5.1** Testing the model as a whole](https://learningstatisticswithr.com/book/regression.html#testing-the-model-as-a-whole)
- [**15\.5.2** Tests for individual coefficients](https://learningstatisticswithr.com/book/regression.html#tests-for-individual-coefficients)
- [**15\.5.3** Running the hypothesis tests in R](https://learningstatisticswithr.com/book/regression.html#regressionsummary)
- [**15\.6** Testing the significance of a correlation](https://learningstatisticswithr.com/book/regression.html#corrhyp)
- [**15\.6.1** Hypothesis tests for a single correlation](https://learningstatisticswithr.com/book/regression.html#hypothesis-tests-for-a-single-correlation)
- [**15\.6.2** Hypothesis tests for all pairwise correlations](https://learningstatisticswithr.com/book/regression.html#corrhyp2)
- [**15\.7** Regarding regression coefficients](https://learningstatisticswithr.com/book/regression.html#regressioncoefs)
- [**15\.7.1** Confidence intervals for the coefficients](https://learningstatisticswithr.com/book/regression.html#confidence-intervals-for-the-coefficients)
- [**15\.7.2** Calculating standardised regression coefficients](https://learningstatisticswithr.com/book/regression.html#calculating-standardised-regression-coefficients)
- [**15\.8** Assumptions of regression](https://learningstatisticswithr.com/book/regression.html#regressionassumptions)
- [**15\.9** Model checking](https://learningstatisticswithr.com/book/regression.html#regressiondiagnostics)
- [**15\.9.1** Three kinds of residuals](https://learningstatisticswithr.com/book/regression.html#three-kinds-of-residuals)
- [**15\.9.2** Three kinds of anomalous data](https://learningstatisticswithr.com/book/regression.html#regressionoutliers)
- [**15\.9.3** Checking the normality of the residuals](https://learningstatisticswithr.com/book/regression.html#regressionnormality)
- [**15\.9.4** Checking the linearity of the relationship](https://learningstatisticswithr.com/book/regression.html#regressionlinearity)
- [**15\.9.5** Checking the homogeneity of variance](https://learningstatisticswithr.com/book/regression.html#regressionhomogeneity)
- [**15\.9.6** Checking for collinearity](https://learningstatisticswithr.com/book/regression.html#regressioncollinearity)
- [**15\.10** Model selection](https://learningstatisticswithr.com/book/regression.html#modelselreg)
- [**15\.10.1** Backward elimination](https://learningstatisticswithr.com/book/regression.html#backward-elimination)
- [**15\.10.2** Forward selection](https://learningstatisticswithr.com/book/regression.html#forward-selection)
- [**15\.10.3** A caveat](https://learningstatisticswithr.com/book/regression.html#a-caveat)
- [**15\.10.4** Comparing two regression models](https://learningstatisticswithr.com/book/regression.html#comparing-two-regression-models)
- [**15\.11** Summary](https://learningstatisticswithr.com/book/regression.html#summary-13)
- [**16** Factorial ANOVA](https://learningstatisticswithr.com/book/anova2.html)
- [**16\.1** Factorial ANOVA 1: balanced designs, no interactions](https://learningstatisticswithr.com/book/anova2.html#factorialanovasimple)
- [**16\.1.1** What hypotheses are we testing?](https://learningstatisticswithr.com/book/anova2.html#factanovahyp)
- [**16\.1.2** Running the analysis in R](https://learningstatisticswithr.com/book/anova2.html#running-the-analysis-in-r)
- [**16\.1.3** How are the sum of squares calculated?](https://learningstatisticswithr.com/book/anova2.html#how-are-the-sum-of-squares-calculated)
- [**16\.1.4** What are our degrees of freedom?](https://learningstatisticswithr.com/book/anova2.html#what-are-our-degrees-of-freedom)
- [**16\.1.5** Factorial ANOVA versus one-way ANOVAs](https://learningstatisticswithr.com/book/anova2.html#factorial-anova-versus-one-way-anovas)
- [**16\.1.6** What kinds of outcomes does this analysis capture?](https://learningstatisticswithr.com/book/anova2.html#what-kinds-of-outcomes-does-this-analysis-capture)
- [**16\.2** Factorial ANOVA 2: balanced designs, interactions allowed](https://learningstatisticswithr.com/book/anova2.html#interactions)
- [**16\.2.1** What exactly *is an interaction effect?*](https://learningstatisticswithr.com/book/anova2.html#what-exactly-is-an-interaction-effect)
- [**16\.2.2** Calculating sums of squares for the interaction](https://learningstatisticswithr.com/book/anova2.html#calculating-sums-of-squares-for-the-interaction)
- [**16\.2.3** Degrees of freedom for the interaction](https://learningstatisticswithr.com/book/anova2.html#degrees-of-freedom-for-the-interaction)
- [**16\.2.4** Running the ANOVA in R](https://learningstatisticswithr.com/book/anova2.html#running-the-anova-in-r)
- [**16\.2.5** Interpreting the results](https://learningstatisticswithr.com/book/anova2.html#interpreting-the-results)
- [**16\.3** Effect size, estimated means, and confidence intervals](https://learningstatisticswithr.com/book/anova2.html#effectsizefactorialanova)
- [**16\.3.1** Effect sizes](https://learningstatisticswithr.com/book/anova2.html#effect-sizes)
- [**16\.3.2** Estimated group means](https://learningstatisticswithr.com/book/anova2.html#estimated-group-means)
- [**16\.4** Assumption checking](https://learningstatisticswithr.com/book/anova2.html#factorialanovaassumptions)
- [**16\.4.1** Levene test for homogeneity of variance](https://learningstatisticswithr.com/book/anova2.html#levene-test-for-homogeneity-of-variance)
- [**16\.4.2** Normality of residuals](https://learningstatisticswithr.com/book/anova2.html#normality-of-residuals)
- [**16\.5** The \\(F\\) test as a model comparison](https://learningstatisticswithr.com/book/anova2.html#omnibusF)
- [**16\.5.1** The \\(F\\) test comparing two models](https://learningstatisticswithr.com/book/anova2.html#the-f-test-comparing-two-models)
- [**16\.6** ANOVA as a linear model](https://learningstatisticswithr.com/book/anova2.html#anovalm)
- [**16\.6.1** Some data](https://learningstatisticswithr.com/book/anova2.html#some-data)
- [**16\.6.2** ANOVA with binary factors as a regression model](https://learningstatisticswithr.com/book/anova2.html#anova-with-binary-factors-as-a-regression-model)
- [**16\.6.3** Changing the baseline category](https://learningstatisticswithr.com/book/anova2.html#changingbaseline)
- [**16\.6.4** How to encode non binary factors as contrasts](https://learningstatisticswithr.com/book/anova2.html#how-to-encode-non-binary-factors-as-contrasts)
- [**16\.6.5** The equivalence between ANOVA and regression for non-binary factors](https://learningstatisticswithr.com/book/anova2.html#the-equivalence-between-anova-and-regression-for-non-binary-factors)
- [**16\.6.6** Degrees of freedom as parameter counting\!](https://learningstatisticswithr.com/book/anova2.html#degrees-of-freedom-as-parameter-counting)
- [**16\.6.7** A postscript](https://learningstatisticswithr.com/book/anova2.html#a-postscript)
- [**16\.7** Different ways to specify contrasts](https://learningstatisticswithr.com/book/anova2.html#contrasts)
- [**16\.7.1** Treatment contrasts](https://learningstatisticswithr.com/book/anova2.html#treatment-contrasts)
- [**16\.7.2** Helmert contrasts](https://learningstatisticswithr.com/book/anova2.html#helmert-contrasts)
- [**16\.7.3** Sum to zero contrasts](https://learningstatisticswithr.com/book/anova2.html#sum-to-zero-contrasts)
- [**16\.7.4** Viewing and setting the default contrasts in R](https://learningstatisticswithr.com/book/anova2.html#viewing-and-setting-the-default-contrasts-in-r)
- [**16\.7.5** Setting the contrasts for a single factor](https://learningstatisticswithr.com/book/anova2.html#setting-the-contrasts-for-a-single-factor)
- [**16\.7.6** Setting the contrasts for a single analysis](https://learningstatisticswithr.com/book/anova2.html#setting-the-contrasts-for-a-single-analysis)
- [**16\.8** Post hoc tests](https://learningstatisticswithr.com/book/anova2.html#posthoc2)
- [**16\.9** The method of planned comparisons](https://learningstatisticswithr.com/book/anova2.html#plannedcomparisons)
- [**16\.10** Factorial ANOVA 3: unbalanced designs](https://learningstatisticswithr.com/book/anova2.html#unbalancedanova)
- [**16\.10.1** The coffee data](https://learningstatisticswithr.com/book/anova2.html#the-coffee-data)
- [**16\.10.2** “Standard ANOVA” does not exist for unbalanced designs](https://learningstatisticswithr.com/book/anova2.html#standard-anova-does-not-exist-for-unbalanced-designs)
- [**16\.10.3** Type I sum of squares](https://learningstatisticswithr.com/book/anova2.html#type-i-sum-of-squares)
- [**16\.10.4** Type III sum of squares](https://learningstatisticswithr.com/book/anova2.html#type-iii-sum-of-squares)
- [**16\.10.5** Type II sum of squares](https://learningstatisticswithr.com/book/anova2.html#type-ii-sum-of-squares)
- [**16\.10.6** Effect sizes (and non-additive sums of squares)](https://learningstatisticswithr.com/book/anova2.html#effect-sizes-and-non-additive-sums-of-squares)
- [**16\.11** Summary](https://learningstatisticswithr.com/book/anova2.html#summary-14)
- [Part VI. Endings, alternatives and prospects](https://learningstatisticswithr.com/book/part-vi-endings-alternatives-and-prospects.html)
- [**17** Bayesian statistics](https://learningstatisticswithr.com/book/bayes.html)
- [**17\.1** Probabilistic reasoning by rational agents](https://learningstatisticswithr.com/book/bayes.html#basicbayes)
- [**17\.1.1** Priors: what you believed before](https://learningstatisticswithr.com/book/bayes.html#priors-what-you-believed-before)
- [**17\.1.2** Likelihoods: theories about the data](https://learningstatisticswithr.com/book/bayes.html#likelihoods-theories-about-the-data)
- [**17\.1.3** The joint probability of data and hypothesis](https://learningstatisticswithr.com/book/bayes.html#the-joint-probability-of-data-and-hypothesis)
- [**17\.1.4** Updating beliefs using Bayes’ rule](https://learningstatisticswithr.com/book/bayes.html#updating-beliefs-using-bayes-rule)
- [**17\.2** Bayesian hypothesis tests](https://learningstatisticswithr.com/book/bayes.html#bayesianhypothesistests)
- [**17\.2.1** The Bayes factor](https://learningstatisticswithr.com/book/bayes.html#the-bayes-factor)
- [**17\.2.2** Interpreting Bayes factors](https://learningstatisticswithr.com/book/bayes.html#interpreting-bayes-factors)
- [**17\.3** Why be a Bayesian?](https://learningstatisticswithr.com/book/bayes.html#whybayes)
- [**17\.3.1** Statistics that mean what you think they mean](https://learningstatisticswithr.com/book/bayes.html#statistics-that-mean-what-you-think-they-mean)
- [**17\.4** Evidentiary standards you can believe](https://learningstatisticswithr.com/book/bayes.html#evidentiary-standards-you-can-believe)
- [**17\.5** The \\(p\\)\-value is a lie.](https://learningstatisticswithr.com/book/bayes.html#the-p-value-is-a-lie.)
- [**17\.5.1** Is it really this bad?](https://learningstatisticswithr.com/book/bayes.html#is-it-really-this-bad)
- [**17\.6** Bayesian analysis of contingency tables](https://learningstatisticswithr.com/book/bayes.html#bayescontingency)
- [**17\.6.1** The orthodox text](https://learningstatisticswithr.com/book/bayes.html#the-orthodox-text)
- [**17\.6.2** The Bayesian test](https://learningstatisticswithr.com/book/bayes.html#the-bayesian-test)
- [**17\.6.3** Writing up the results](https://learningstatisticswithr.com/book/bayes.html#writing-up-the-results)
- [**17\.6.4** Other sampling plans](https://learningstatisticswithr.com/book/bayes.html#other-sampling-plans)
- [**17\.7** Bayesian \\(t\\)\-tests](https://learningstatisticswithr.com/book/bayes.html#ttestbf)
- [**17\.7.1** Independent samples \\(t\\)\-test](https://learningstatisticswithr.com/book/bayes.html#independent-samples-t-test)
- [**17\.7.2** Paired samples \\(t\\)\-test](https://learningstatisticswithr.com/book/bayes.html#paired-samples-t-test)
- [**17\.8** Bayesian regression](https://learningstatisticswithr.com/book/bayes.html#bayesregression)
- [**17\.8.1** A quick refresher](https://learningstatisticswithr.com/book/bayes.html#a-quick-refresher)
- [**17\.8.2** The Bayesian version](https://learningstatisticswithr.com/book/bayes.html#the-bayesian-version)
- [**17\.8.3** Finding the best model](https://learningstatisticswithr.com/book/bayes.html#finding-the-best-model)
- [**17\.8.4** Extracting Bayes factors for all included terms](https://learningstatisticswithr.com/book/bayes.html#extracting-bayes-factors-for-all-included-terms)
- [**17\.9** Bayesian ANOVA](https://learningstatisticswithr.com/book/bayes.html#bayesanova)
- [**17\.9.1** A quick refresher](https://learningstatisticswithr.com/book/bayes.html#a-quick-refresher-1)
- [**17\.9.2** The Bayesian version](https://learningstatisticswithr.com/book/bayes.html#the-bayesian-version-1)
- [**17\.9.3** Constructing Bayesian Type II tests](https://learningstatisticswithr.com/book/bayes.html#constructing-bayesian-type-ii-tests)
- [**17\.10** Summary](https://learningstatisticswithr.com/book/bayes.html#summary-15)
- [Epilogue](https://learningstatisticswithr.com/book/epilogue.html)
- [The undiscovered statistics](https://learningstatisticswithr.com/book/epilogue.html#the-undiscovered-statistics)
- [Omissions within the topics covered](https://learningstatisticswithr.com/book/epilogue.html#omissions-within-the-topics-covered)
- [Statistical models missing from the book](https://learningstatisticswithr.com/book/epilogue.html#statistical-models-missing-from-the-book)
- [Other ways of doing inference](https://learningstatisticswithr.com/book/epilogue.html#other-ways-of-doing-inference)
- [Miscellaneous topics](https://learningstatisticswithr.com/book/epilogue.html#miscellaneous-topics)
- [Learning the basics, and learning them in R](https://learningstatisticswithr.com/book/epilogue.html#learning-the-basics-and-learning-them-in-r)
- [References](https://learningstatisticswithr.com/book/references.html)
- [Published with bookdown](https://learningstatisticswithr.com/book/)
# [Learning statistics with R: A tutorial for psychology students and other beginners. (Version 0.6.1)](https://learningstatisticswithr.com/book/)
# Chapter 17 Bayesian statistics
> *In our reasonings concerning matter of fact, there are all imaginable degrees of assurance, from the highest certainty to the lowest species of moral evidence. A wise man, therefore, proportions his belief to the evidence.* – David Hume[254](https://learningstatisticswithr.com/book/bayes.html#fn254).
The ideas I’ve presented to you in this book describe inferential statistics from the frequentist perspective. I’m not alone in doing this. In fact, almost every textbook given to undergraduate psychology students presents the opinions of the frequentist statistician as *the* theory of inferential statistics, the one true way to do things. I have taught this way for practical reasons. The frequentist view of statistics dominated the academic field of statistics for most of the 20th century, and this dominance is even more extreme among applied scientists. It was and is current practice among psychologists to use frequentist methods. Because frequentist methods are ubiquitous in scientific papers, every student of statistics needs to understand those methods, otherwise they will be unable to make sense of what those papers are saying! Unfortunately – in my opinion at least – the current practice in psychology is often misguided, and the reliance on frequentist methods is partly to blame. In this chapter I explain why I think this, and provide an introduction to Bayesian statistics, an approach that I think is generally superior to the orthodox approach.
This chapter comes in two parts. In Sections [17\.1](https://learningstatisticswithr.com/book/bayes.html#basicbayes) through [17\.3](https://learningstatisticswithr.com/book/bayes.html#whybayes) I talk about what Bayesian statistics are all about, covering the basic mathematical rules for how it works as well as an explanation for why I think the Bayesian approach is so useful. Afterwards, I provide a brief overview of how you can do Bayesian versions of chi-square tests (Section [17\.6](https://learningstatisticswithr.com/book/bayes.html#bayescontingency)), \\(t\\)\-tests (Section [17\.7](https://learningstatisticswithr.com/book/bayes.html#ttestbf)), regression (Section [17\.8](https://learningstatisticswithr.com/book/bayes.html#bayesregression)) and ANOVA (Section [17\.9](https://learningstatisticswithr.com/book/bayes.html#bayesanova)).
## 17\.1 Probabilistic reasoning by rational agents
From a Bayesian perspective, statistical inference is all about *belief revision*. I start out with a set of candidate hypotheses \\(h\\) about the world. I don’t know which of these hypotheses is true, but do I have some beliefs about which hypotheses are plausible and which are not. When I observe the data \\(d\\), I have to revise those beliefs. If the data are consistent with a hypothesis, my belief in that hypothesis is strengthened. If the data inconsistent with the hypothesis, my belief in that hypothesis is weakened. That’s it! At the end of this section I’ll give a precise description of how Bayesian reasoning works, but first I want to work through a simple example in order to introduce the key ideas. Consider the following reasoning problem:
> *I’m carrying an umbrella. Do you think it will rain?*
In this problem, I have presented you with a single piece of data (\\(d =\\) I’m carrying the umbrella), and I’m asking you to tell me your beliefs about whether it’s raining. You have two possible ***hypotheses***, \\(h\\): either it rains today or it does not. How should you solve this problem?
### 17\.1.1 Priors: what you believed before
The first thing you need to do ignore what I told you about the umbrella, and write down your pre-existing beliefs about rain. This is important: if you want to be honest about how your beliefs have been revised in the light of new evidence, then you *must* say something about what you believed before those data appeared! So, what might you believe about whether it will rain today? You probably know that I live in Australia, and that much of Australia is hot and dry. And in fact you’re right: the city of Adelaide where I live has a Mediterranean climate, very similar to southern California, southern Europe or northern Africa. I’m writing this in January, and so you can assume it’s the middle of summer. In fact, you might have decided to take a quick look on Wikipedia[255](https://learningstatisticswithr.com/book/bayes.html#fn255) and discovered that Adelaide gets an average of 4.4 days of rain across the 31 days of January. Without knowing anything else, you might conclude that the probability of January rain in Adelaide is about 15%, and the probability of a dry day is 85%. If this is really what you believe about Adelaide rainfall (and now that I’ve told it to you, I’m betting that this really *is* what you believe) then what I have written here is your ***prior distribution***, written \\(P(h)\\):
| Hypothesis | Degree of Belief |
|---|---|
| Rainy day | 0\.15 |
| Dry day | 0\.85 |
### 17\.1.2 Likelihoods: theories about the data
To solve the reasoning problem, you need a theory about my behaviour. When does Dan carry an umbrella? You might guess that I’m not a complete idiot,[256](https://learningstatisticswithr.com/book/bayes.html#fn256) and I try to carry umbrellas only on rainy days. On the other hand, you also know that I have young kids, and you wouldn’t be all that surprised to know that I’m pretty forgetful about this sort of thing. Let’s suppose that on rainy days I remember my umbrella about 30% of the time (I really am awful at this). But let’s say that on dry days I’m only about 5% likely to be carrying an umbrella. So you might write out a little table like this:
| Hypothesis | Umbrella | No umbrella |
|---|---|---|
| Rainy day | 0\.30 | 0\.70 |
| Dry day | 0\.05 | 0\.95 |
It’s important to remember that each cell in this table describes your beliefs about what data \\(d\\) will be observed, *given* the truth of a particular hypothesis \\(h\\). This “conditional probability” is written \\(P(d\|h)\\), which you can read as “the probability of \\(d\\) given \\(h\\)”. In Bayesian statistics, this is referred to as ***likelihood*** of data \\(d\\) given hypothesis \\(h\\).[257](https://learningstatisticswithr.com/book/bayes.html#fn257)
### 17\.1.3 The joint probability of data and hypothesis
At this point, all the elements are in place. Having written down the priors and the likelihood, you have all the information you need to do Bayesian reasoning. The question now becomes, *how* do we use this information? As it turns out, there’s a very simple equation that we can use here, but it’s important that you understand why we use it, so I’m going to try to build it up from more basic ideas.
Let’s start out with one of the rules of probability theory. I listed it way back in Table [9\.1](https://learningstatisticswithr.com/book/probability.html#tab:probrules), but I didn’t make a big deal out of it at the time and you probably ignored it. The rule in question is the one that talks about the probability that *two* things are true. In our example, you might want to calculate the probability that today is rainy (i.e., hypothesis \\(h\\) is true) *and* I’m carrying an umbrella (i.e., data \\(d\\) is observed). The ***joint probability*** of the hypothesis and the data is written \\(P(d,h)\\), and you can calculate it by multiplying the prior \\(P(h)\\) by the likelihood \\(P(d\|h)\\). Mathematically, we say that: \\\[ P(d,h) = P(d\|h) P(h) \\\]
So, what is the probability that today is a rainy day *and* I remember to carry an umbrella? As we discussed earlier, the prior tells us that the probability of a rainy day is 15%, and the likelihood tells us that the probability of me remembering my umbrella on a rainy day is 30%. So the probability that both of these things are true is calculated by multiplying the two:
\\\[ \\begin{array} P(\\mbox{rainy}, \\mbox{umbrella}) & = & P(\\mbox{umbrella} \| \\mbox{rainy}) \\times P(\\mbox{rainy}) \\\\ & = & 0.30 \\times 0.15 \\\\ & = & 0.045 \\end{array} \\\]
In other words, before being told anything about what actually happened, you think that there is a 4.5% probability that today will be a rainy day and that I will remember an umbrella. However, there are of course *four* possible things that could happen, right? So let’s repeat the exercise for all four. If we do that, we end up with the following table:
| | Umbrella | No-umbrella |
|---|---|---|
| Rainy | 0\.045 | 0\.105 |
| Dry | 0\.0425 | 0\.8075 |
This table captures all the information about which of the four possibilities are likely. To really get the full picture, though, it helps to add the row totals and column totals. That gives us this table:
| | Umbrella | No-umbrella | Total |
|---|---|---|---|
| Rainy | 0\.0450 | 0\.1050 | 0\.15 |
| Dry | 0\.0425 | 0\.8075 | 0\.85 |
| Total | 0\.0875 | 0\.9125 | 1 |
This is a very useful table, so it’s worth taking a moment to think about what all these numbers are telling us. First, notice that the row sums aren’t telling us anything new at all. For example, the first row tells us that if we ignore all this umbrella business, the chance that today will be a rainy day is 15%. That’s not surprising, of course: that’s our prior. The important thing isn’t the number itself: rather, the important thing is that it gives us some confidence that our calculations are sensible! Now take a look at the column sums, and notice that they tell us something that we haven’t explicitly stated yet. In the same way that the row sums tell us the probability of rain, the column sums tell us the probability of me carrying an umbrella. Specifically, the first column tells us that on average (i.e., ignoring whether it’s a rainy day or not), the probability of me carrying an umbrella is 8.75%. Finally, notice that when we sum across all four logically-possible events, everything adds up to 1. In other words, what we have written down is a proper probability distribution defined over all possible combinations of data and hypothesis.
Now, because this table is so useful, I want to make sure you understand what all the elements correspond to, and how they written:
| | Umbrella | No-umbrella | |
|---|---|---|---|
| Rainy | \\(P\\)(Umbrella, Rainy) | \\(P\\)(No-umbrella, Rainy) | \\(P\\)(Rainy) |
| Dry | \\(P\\)(Umbrella, Dry) | \\(P\\)(No-umbrella, Dry) | \\(P\\)(Dry) |
| | \\(P\\)(Umbrella) | \\(P\\)(No-umbrella) | |
Finally, let’s use “proper” statistical notation. In the rainy day problem, the data corresponds to the observation that I do or do not have an umbrella. So we’ll let \\(d\_1\\) refer to the possibility that you observe me carrying an umbrella, and \\(d\_2\\) refers to you observing me not carrying one. Similarly, \\(h\_1\\) is your hypothesis that today is rainy, and \\(h\_2\\) is the hypothesis that it is not. Using this notation, the table looks like this:
| | \\(d\_1\\) | \\(d\_2\\) | |
|---|---|---|---|
| \\(h\_1\\) | \\(P(h\_1, d\_1)\\) | \\(P(h\_1, d\_2)\\) | \\(P(h\_1)\\) |
| \\(h\_2\\) | \\(P(h\_2, d\_1)\\) | \\(P(h\_2, d\_2)\\) | \\(P(h\_2)\\) |
| | \\(P(d\_1)\\) | \\(P(d\_2)\\) | |
### 17\.1.4 Updating beliefs using Bayes’ rule
The table we laid out in the last section is a very powerful tool for solving the rainy day problem, because it considers all four logical possibilities and states exactly how confident you are in each of them before being given any data. It’s now time to consider what happens to our beliefs when we are actually given the data. In the rainy day problem, you are told that I really *am* carrying an umbrella. This is something of a surprising event: according to our table, the probability of me carrying an umbrella is only 8.75%. But that makes sense, right? A guy carrying an umbrella on a summer day in a hot dry city is pretty unusual, and so you really weren’t expecting that. Nevertheless, the problem tells you that it is true. No matter how unlikely you thought it was, you must now adjust your beliefs to accommodate the fact that you now *know* that I have an umbrella.[258](https://learningstatisticswithr.com/book/bayes.html#fn258) To reflect this new knowledge, our *revised* table must have the following numbers:
| | Umbrella | No-umbrella |
|---|---|---|
| Rainy | | 0 |
| Dry | | 0 |
| Total | 1 | 0 |
In other words, the facts have eliminated any possibility of “no umbrella”, so we have to put zeros into any cell in the table that implies that I’m not carrying an umbrella. Also, you know for a fact that I am carrying an umbrella, so the column sum on the left must be 1 to correctly describe the fact that \\(P(\\mbox{umbrella})=1\\).
What two numbers should we put in the empty cells? Again, let’s not worry about the maths, and instead think about our intuitions. When we wrote out our table the first time, it turned out that those two cells had almost identical numbers, right? We worked out that the joint probability of “rain and umbrella” was 4.5%, and the joint probability of “dry and umbrella” was 4.25%. In other words, before I told you that I am in fact carrying an umbrella, you’d have said that these two events were almost identical in probability, yes? But notice that *both* of these possibilities are consistent with the fact that I actually am carrying an umbrella. From the perspective of these two possibilities, very little has changed. I hope you’d agree that it’s *still* true that these two possibilities are equally plausible. So what we expect to see in our final table is some numbers that preserve the fact that “rain and umbrella” is *slightly* more plausible than “dry and umbrella”, while still ensuring that numbers in the table add up. Something like this, perhaps?
| | Umbrella | No-umbrella |
|---|---|---|
| Rainy | 0\.514 | 0 |
| Dry | 0\.486 | 0 |
| Total | 1 | 0 |
What this table is telling you is that, after being told that I’m carrying an umbrella, you believe that there’s a 51.4% chance that today will be a rainy day, and a 48.6% chance that it won’t. That’s the answer to our problem! The ***posterior probability*** of rain \\(P(h\|d)\\) given that I am carrying an umbrella is 51.4%
How did I calculate these numbers? You can probably guess. To work out that there was a 0.514 probability of “rain”, all I did was take the 0.045 probability of “rain and umbrella” and divide it by the 0.0875 chance of “umbrella”. This produces a table that satisfies our need to have everything sum to 1, and our need not to interfere with the relative plausibility of the two events that are actually consistent with the data. To say the same thing using fancy statistical jargon, what I’ve done here is divide the joint probability of the hypothesis and the data \\(P(d,h)\\) by the ***marginal probability*** of the data \\(P(d)\\), and this is what gives us the posterior probability of the hypothesis *given* that we know the data have been observed. To write this as an equation:[259](https://learningstatisticswithr.com/book/bayes.html#fn259) \\\[ P(h \| d) = \\frac{P(d,h)}{P(d)} \\\]
However, remember what I said at the start of the last section, namely that the joint probability \\(P(d,h)\\) is calculated by multiplying the prior \\(P(h)\\) by the likelihood \\(P(d\|h)\\). In real life, the things we actually know how to write down are the priors and the likelihood, so let’s substitute those back into the equation. This gives us the following formula for the posterior probability:
\\\[ P(h \| d) = \\frac{P(d\|h) P(h)}{P(d)} \\\]
And this formula, folks, is known as ***Bayes’ rule***. It describes how a learner starts out with prior beliefs about the plausibility of different hypotheses, and tells you how those beliefs should be revised in the face of data. In the Bayesian paradigm, all statistical inference flows from this one simple rule.
## 17\.2 Bayesian hypothesis tests
In Chapter [11](https://learningstatisticswithr.com/book/hypothesistesting.html#hypothesistesting) I described the orthodox approach to hypothesis testing. It took an entire chapter to describe, because null hypothesis testing is a very elaborate contraption that people find very hard to make sense of. In contrast, the Bayesian approach to hypothesis testing is incredibly simple. Let’s pick a setting that is closely analogous to the orthodox scenario. There are two hypotheses that we want to compare, a null hypothesis \\(h\_0\\) and an alternative hypothesis \\(h\_1\\). Prior to running the experiment we have some beliefs \\(P(h)\\) about which hypotheses are true. We run an experiment and obtain data \\(d\\). Unlike frequentist statistics Bayesian statistics does allow to talk about the probability that the null hypothesis is true. Better yet, it allows us to calculate the ***posterior probability of the null hypothesis***, using Bayes’ rule:
\\\[ P(h\_0 \| d) = \\frac{P(d\|h\_0) P(h\_0)}{P(d)} \\\]
This formula tells us exactly how much belief we should have in the null hypothesis after having observed the data \\(d\\). Similarly, we can work out how much belief to place in the alternative hypothesis using essentially the same equation. All we do is change the subscript:
\\\[ P(h\_1 \| d) = \\frac{P(d\|h\_1) P(h\_1)}{P(d)} \\\]
It’s all so simple that I feel like an idiot even bothering to write these equations down, since all I’m doing is copying Bayes rule from the previous section.[260](https://learningstatisticswithr.com/book/bayes.html#fn260)
### 17\.2.1 The Bayes factor
In practice, most Bayesian data analysts tend not to talk in terms of the raw posterior probabilities \\(P(h\_0\|d)\\) and \\(P(h\_1\|d)\\). Instead, we tend to talk in terms of the ***posterior odds*** ratio. Think of it like betting. Suppose, for instance, the posterior probability of the null hypothesis is 25%, and the posterior probability of the alternative is 75%. The alternative hypothesis is three times as probable as the null, so we say that the *odds* are 3:1 in favour of the alternative. Mathematically, all we have to do to calculate the posterior odds is divide one posterior probability by the other:
\\\[ \\frac{P(h\_1 \| d)}{P(h\_0 \| d)} = \\frac{0.75}{0.25} = 3 \\\]
Or, to write the same thing in terms of the equations above: \\\[ \\frac{P(h\_1 \| d)}{P(h\_0 \| d)} = \\frac{P(d\|h\_1)}{P(d\|h\_0)} \\times \\frac{P(h\_1)}{P(h\_0)} \\\]
Actually, this equation is worth expanding on. There are three different terms here that you should know. On the left hand side, we have the posterior odds, which tells you what you believe about the relative plausibilty of the null hypothesis and the alternative hypothesis *after* seeing the data. On the right hand side, we have the ***prior odds***, which indicates what you thought *before* seeing the data. In the middle, we have the ***Bayes factor***, which describes the amount of evidence provided by the data: \\\[ \\begin{array}{ccccc}\\displaystyle \\frac{P(h\_1 \| d)}{P(h\_0 \| d)} &=& \\displaystyle\\frac{P(d\|h\_1)}{P(d\|h\_0)} &\\times& \\displaystyle\\frac{P(h\_1)}{P(h\_0)} \\\\\[6pt\] \\\\\[-2pt\] \\uparrow && \\uparrow && \\uparrow \\\\\[6pt\] \\mbox{Posterior odds} && \\mbox{Bayes factor} && \\mbox{Prior odds} \\end{array} \\\] The Bayes factor (sometimes abbreviated as ***BF***) has a special place in the Bayesian hypothesis testing, because it serves a similar role to the \\(p\\)\-value in orthodox hypothesis testing: it quantifies the strength of evidence provided by the data, and as such it is the Bayes factor that people tend to report when running a Bayesian hypothesis test. The reason for reporting Bayes factors rather than posterior odds is that different researchers will have different priors. Some people might have a strong bias to believe the null hypothesis is true, others might have a strong bias to believe it is false. Because of this, the polite thing for an applied researcher to do is report the Bayes factor. That way, anyone reading the paper can multiply the Bayes factor by their own *personal* prior odds, and they can work out for themselves what the posterior odds would be. In any case, by convention we like to pretend that we give equal consideration to both the null hypothesis and the alternative, in which case the prior odds equals 1, and the posterior odds becomes the same as the Bayes factor.
### 17\.2.2 Interpreting Bayes factors
One of the really nice things about the Bayes factor is the numbers are inherently meaningful. If you run an experiment and you compute a Bayes factor of 4, it means that the evidence provided by your data corresponds to betting odds of 4:1 in favour of the alternative. However, there have been some attempts to quantify the standards of evidence that would be considered meaningful in a scientific context. The two most widely used are from Jeffreys ([1961](https://learningstatisticswithr.com/book/bayes.html#ref-Jeffreys1961)) and Kass and Raftery ([1995](https://learningstatisticswithr.com/book/bayes.html#ref-Kass1995)). Of the two, I tend to prefer the Kass and Raftery ([1995](https://learningstatisticswithr.com/book/bayes.html#ref-Kass1995)) table because it’s a bit more conservative. So here it is:
| Bayes factor | Interpretation |
|---|---|
| 1 - 3 | Negligible evidence |
| 3 - 20 | Positive evidence |
| 20 - 150 | Strong evidence |
| \$\>\$150 | Very strong evidence |
And to be perfectly honest, I think that even the Kass and Raftery standards are being a bit charitable. If it were up to me, I’d have called the “positive evidence” category “weak evidence”. To me, anything in the range 3:1 to 20:1 is “weak” or “modest” evidence at best. But there are no hard and fast rules here: what counts as strong or weak evidence depends entirely on how conservative you are, and upon the standards that your community insists upon before it is willing to label a finding as “true”.
In any case, note that all the numbers listed above make sense if the Bayes factor is greater than 1 (i.e., the evidence favours the alternative hypothesis). However, one big practical advantage of the Bayesian approach relative to the orthodox approach is that it also allows you to quantify evidence *for* the null. When that happens, the Bayes factor will be less than 1. You can choose to report a Bayes factor less than 1, but to be honest I find it confusing. For example, suppose that the likelihood of the data under the null hypothesis \\(P(d\|h\_0)\\) is equal to 0.2, and the corresponding likelihood \\(P(d\|h\_0)\\) under the alternative hypothesis is 0.1. Using the equations given above, Bayes factor here would be:
\\\[ \\mbox{BF} = \\frac{P(d\|h\_1)}{P(d\|h\_0)} = \\frac{0.1}{0.2} = 0.5 \\\]
Read literally, this result tells is that the evidence in favour of the alternative is 0.5 to 1. I find this hard to understand. To me, it makes a lot more sense to turn the equation “upside down”, and report the amount op evidence in favour of the *null*. In other words, what we calculate is this:
\\\[ \\mbox{BF}^\\prime = \\frac{P(d\|h\_0)}{P(d\|h\_1)} = \\frac{0.2}{0.1} = 2 \\\]
And what we would report is a Bayes factor of 2:1 in favour of the null. Much easier to understand, and you can interpret this using the table above.
## 17\.3 Why be a Bayesian?
Up to this point I’ve focused exclusively on the logic underpinning Bayesian statistics. We’ve talked about the idea of “probability as a degree of belief”, and what it implies about how a rational agent should reason about the world. The question that you have to answer for yourself is this: how do *you* want to do your statistics? Do you want to be an orthodox statistician, relying on sampling distributions and \\(p\\)\-values to guide your decisions? Or do you want to be a Bayesian, relying on Bayes factors and the rules for rational belief revision? And to be perfectly honest, I can’t answer this question for you. Ultimately it depends on what you think is right. It’s your call, and your call alone. That being said, I can talk a little about why *I* prefer the Bayesian approach.
### 17\.3.1 Statistics that mean what you think they mean
> *You keep using that word. I do not think it means what you think it means*
> – Inigo Montoya, The Princess Bride[261](https://learningstatisticswithr.com/book/bayes.html#fn261)
To me, one of the biggest advantages to the Bayesian approach is that it answers the right questions. Within the Bayesian framework, it is perfectly sensible and allowable to refer to “the probability that a hypothesis is true”. You can even try to calculate this probability. Ultimately, isn’t that what you *want* your statistical tests to tell you? To an actual human being, this would seem to be the whole *point* of doing statistics: to determine what is true and what isn’t. Any time that you aren’t exactly sure about what the truth is, you should use the language of probability theory to say things like “there is an 80% chance that Theory A is true, but a 20% chance that Theory B is true instead”.
This seems so obvious to a human, yet it is explicitly forbidden within the orthodox framework. To a frequentist, such statements are a nonsense because “the theory is true” is not a repeatable event. A theory is true or it is not, and no probabilistic statements are allowed, no matter how much you might want to make them. There’s a reason why, back in Section [11\.5](https://learningstatisticswithr.com/book/hypothesistesting.html#pvalue), I repeatedly warned you *not* to interpret the \\(p\\)\-value as the probability of that the null hypothesis is true. There’s a reason why almost every textbook on statstics is forced to repeat that warning. It’s because people desperately *want* that to be the correct interpretation. Frequentist dogma notwithstanding, a lifetime of experience of teaching undergraduates and of doing data analysis on a daily basis suggests to me that most actual humans thing that “the probability that the hypothesis is true” is not only meaningful, it’s the thing we care *most* about. It’s such an appealing idea that even trained statisticians fall prey to the mistake of trying to interpret a \\(p\\)\-value this way. For example, here is a quote from an official Newspoll report in 2013, explaining how to interpret their (frequentist) data analysis:[262](https://learningstatisticswithr.com/book/bayes.html#fn262)
> Throughout the report, where relevant, statistically significant changes have been noted. All significance tests have been based on the 95 percent level of confidence. **This means that if a change is noted as being statistically significant, there is a 95 percent probability that a real change has occurred**, and is not simply due to chance variation. (emphasis added)
Nope! That’s *not* what \\(p\<.05\\) means. That’s *not* what 95% confidence means to a frequentist statistician. The bolded section is just plain wrong. Orthodox methods cannot tell you that “there is a 95% chance that a real change has occurred”, because this is not the kind of event to which frequentist probabilities may be assigned. To an ideological frequentist, this sentence should be meaningless. Even if you’re a more pragmatic frequentist, it’s still the wrong definition of a \\(p\\)\-value. It is simply not an allowed or correct thing to say if you want to rely on orthodox statistical tools.
On the other hand, let’s suppose you are a Bayesian. Although the bolded passage is the wrong definition of a \\(p\\)\-value, it’s pretty much exactly what a Bayesian means when they say that the posterior probability of the alternative hypothesis is greater than 95%. And here’s the thing. If the Bayesian posterior is actually thing you *want* to report, why are you even trying to use orthodox methods? If you want to make Bayesian claims, all you have to do is be a Bayesian and use Bayesian tools.
Speaking for myself, I found this to be a the most liberating thing about switching to the Bayesian view. Once you’ve made the jump, you no longer have to wrap your head around counterinuitive definitions of \\(p\\)\-values. You don’t have to bother remembering why you can’t say that you’re 95% confident that the true mean lies within some interval. All you have to do is be honest about what you believed before you ran the study, and then report what you learned from doing it. Sounds nice, doesn’t it? To me, this is the big promise of the Bayesian approach: you do the analysis you really want to do, and express what you really believe the data are telling you.
## 17\.4 Evidentiary standards you can believe
> *If \[\\(p\\)\] is below .02 it is strongly indicated that the \[null\] hypothesis fails to account for the whole of the facts. We shall not often be astray if we draw a conventional line at .05 and consider that \[smaller values of \\(p\\)\] indicate a real discrepancy.*
> – Sir Ronald Fisher ([1925](https://learningstatisticswithr.com/book/bayes.html#ref-Fisher1925))
Consider the quote above by Sir Ronald Fisher, one of the founders of what has become the orthodox approach to statistics. If anyone has ever been entitled to express an opinion about the intended function of \\(p\\)\-values, it’s Fisher. In this passage, taken from his classic guide *Statistical Methods for Research Workers*, he’s pretty clear about what it means to reject a null hypothesis at \\(p\<.05\\). In his opinion, if we take \\(p\<.05\\) to mean there is “a real effect”, then “we shall not often be astray”. This view is hardly unusual: in my experience, most practitioners express views very similar to Fisher’s. In essence, the \\(p\<.05\\) convention is assumed to represent a fairly stringent evidentiary standard.
Well, how true is that? One way to approach this question is to try to convert \\(p\\)\-values to Bayes factors, and see how the two compare. It’s not an easy thing to do because a \\(p\\)\-value is a fundamentally different kind of calculation to a Bayes factor, and they don’t measure the same thing. However, there have been some attempts to work out the relationship between the two, and it’s somewhat surprising. For example, Johnson ([2013](https://learningstatisticswithr.com/book/bayes.html#ref-Johnson2013)) presents a pretty compelling case that (for \\(t\\)\-tests at least) the \\(p\<.05\\) threshold corresponds roughly to a Bayes factor of somewhere between 3:1 and 5:1 in favour of the alternative. If that’s right, then Fisher’s claim is a bit of a stretch. Let’s suppose that the null hypothesis is true about half the time (i.e., the prior probability of \\(H\_0\\) is 0.5), and we use those numbers to work out the posterior probability of the null hypothesis given that it has been rejected at \\(p\<.05\\). Using the data from Johnson ([2013](https://learningstatisticswithr.com/book/bayes.html#ref-Johnson2013)), we see that if you reject the null at \\(p\<.05\\), you’ll be correct about 80% of the time. I don’t know about you, but in my opinion an evidentiary standard that ensures you’ll be wrong on 20% of your decisions isn’t good enough. The fact remains that, quite contrary to Fisher’s claim, if you reject at \\(p\<.05\\) you shall quite often go astray. It’s not a very stringent evidentiary threshold at all.
## 17\.5 The \\(p\\)\-value is a lie.
> *The cake is a lie.*
> *The cake is a lie.*
> *The cake is a lie.*
> *The cake is a lie.*
> – Portal[263](https://learningstatisticswithr.com/book/bayes.html#fn263)
Okay, at this point you might be thinking that the real problem is not with orthodox statistics, just the \\(p\<.05\\) standard. In one sense, that’s true. The recommendation that Johnson ([2013](https://learningstatisticswithr.com/book/bayes.html#ref-Johnson2013)) gives is not that “everyone must be a Bayesian now”. Instead, the suggestion is that it would be wiser to shift the conventional standard to something like a \\(p\<.01\\) level. That’s not an unreasonable view to take, but in my view the problem is a little more severe than that. In my opinion, there’s a fairly big problem built into the way most (but not all) orthodox hypothesis tests are constructed. They are grossly naive about how humans actually do research, and because of this most \\(p\\)\-values are wrong.
Sounds like an absurd claim, right? Well, consider the following scenario. You’ve come up with a really exciting research hypothesis and you design a study to test it. You’re very diligent, so you run a power analysis to work out what your sample size should be, and you run the study. You run your hypothesis test and out pops a \\(p\\)\-value of 0.072. Really bloody annoying, right?
What should you do? Here are some possibilities:
1. You conclude that there is no effect, and try to publish it as a null result
2. You guess that there might be an effect, and try to publish it as a “borderline significant” result
3. You give up and try a new study
4. You collect some more data to see if the \\(p\\) value goes up or (preferably!) drops below the “magic” criterion of \\(p\<.05\\)
Which would *you* choose? Before reading any further, I urge you to take some time to think about it. Be honest with yourself. But don’t stress about it too much, because you’re screwed no matter what you choose. Based on my own experiences as an author, reviewer and editor, as well as stories I’ve heard from others, here’s what will happen in each case:
- Let’s start with option 1. If you try to publish it as a null result, the paper will struggle to be published. Some reviewers will think that \\(p=.072\\) is not really a null result. They’ll argue it’s borderline significant. Other reviewers will agree it’s a null result, but will claim that even though some null results *are* publishable, yours isn’t. One or two reviewers might even be on your side, but you’ll be fighting an uphill battle to get it through.
- Okay, let’s think about option number 2. Suppose you try to publish it as a borderline significant result. Some reviewers will claim that it’s a null result and should not be published. Others will claim that the evidence is ambiguous, and that you should collect more data until you get a clear significant result. Again, the publication process does not favour you.
- Given the difficulties in publishing an “ambiguous” result like \\(p=.072\\), option number 3 might seem tempting: give up and do something else. But that’s a recipe for career suicide. If you give up and try a new project else every time you find yourself faced with ambiguity, your work will never be published. And if you’re in academia without a publication record you can lose your job. So that option is out.
- It looks like you’re stuck with option 4. You don’t have conclusive results, so you decide to collect some more data and re-run the analysis. Seems sensible, but unfortunately for you, if you do this all of your \\(p\\)\-values are now incorrect. *All* of them. Not just the \\(p\\)\-values that you calculated for *this* study. All of them. All the \\(p\\)\-values you calculated in the past and all the \\(p\\)\-values you will calculate in the future. Fortunately, no-one will notice. You’ll get published, and you’ll have lied.
Wait, what? How can that last part be true? I mean, it sounds like a perfectly reasonable strategy doesn’t it? You collected some data, the results weren’t conclusive, so now what you want to do is collect more data until the the results *are* conclusive. What’s wrong with that?
Honestly, there’s nothing wrong with it. It’s a reasonable, sensible and rational thing to do. In real life, this is exactly what every researcher does. Unfortunately, the theory of null hypothesis testing as I described it in Chapter [11](https://learningstatisticswithr.com/book/hypothesistesting.html#hypothesistesting) *forbids* you from doing this.[264](https://learningstatisticswithr.com/book/bayes.html#fn264) The reason is that the theory assumes that the experiment is finished and all the data are in. And because it assumes the experiment is over, it only considers *two* possible decisions. If you’re using the conventional \\(p\<.05\\) threshold, those decisions are:
| Outcome | Action |
|---|---|
| \\(p\\) less than .05 | Reject the null |
| \\(p\\) greater than .05 | Retain the null |
What *you’re* doing is adding a third possible action to the decision making problem. Specifically, what you’re doing is using the \\(p\\)\-value itself as a reason to justify continuing the experiment. And as a consequence you’ve transformed the decision-making procedure into one that looks more like this:
| Outcome | Action |
|---|---|
| \\(p\\) less than .05 | Stop the experiment and reject the null |
| \\(p\\) between .05 and .1 | Continue the experiment |
| \\(p\\) greater than .1 | Stop the experiment and retain the null |
The “basic” theory of null hypothesis testing isn’t built to handle this sort of thing, not in the form I described back in Chapter [11](https://learningstatisticswithr.com/book/hypothesistesting.html#hypothesistesting). If you’re the kind of person who would choose to “collect more data” in real life, it implies that you are *not* making decisions in accordance with the rules of null hypothesis testing. Even if you happen to arrive at the same decision as the hypothesis test, you aren’t following the decision *process* it implies, and it’s this failure to follow the process that is causing the problem.[265](https://learningstatisticswithr.com/book/bayes.html#fn265) Your \\(p\\)\-values are a lie.
Worse yet, they’re a lie in a dangerous way, because they’re all *too small*. To give you a sense of just how bad it can be, consider the following (worst case) scenario. Imagine you’re a really super-enthusiastic researcher on a tight budget who didn’t pay any attention to my warnings above. You design a study comparing two groups. You desperately want to see a significant result at the \\(p\<.05\\) level, but you really don’t want to collect any more data than you have to (because it’s expensive). In order to cut costs, you start collecting data, but every time a new observation arrives you run a \\(t\\)\-test on your data. If the \\(t\\)\-tests says \\(p\<.05\\) then you stop the experiment and report a significant result. If not, you keep collecting data. You keep doing this until you reach your pre-defined spending limit for this experiment. Let’s say that limit kicks in at \\(N=1000\\) observations. As it turns out, the truth of the matter is that there is no real effect to be found: the null hypothesis is true. So, what’s the chance that you’ll make it to the end of the experiment and (correctly) conclude that there is no effect? In an ideal world, the answer here should be 95%. After all, the whole *point* of the \\(p\<.05\\) criterion is to control the Type I error rate at 5%, so what we’d hope is that there’s only a 5% chance of falsely rejecting the null hypothesis in this situation. However, there’s no guarantee that will be true. You’re breaking the rules: you’re running tests repeatedly, “peeking” at your data to see if you’ve gotten a significant result, and all bets are off.
Figure 17.1: How badly can things go wrong if you re-run your tests every time new data arrive? If you are a frequentist, the answer is “very wrong”.
So how bad is it? The answer is shown as the solid black line in Figure [17\.1](https://learningstatisticswithr.com/book/bayes.html#fig:type1), and it’s *astoundingly* bad. If you peek at your data after every single observation, there is a 49% chance that you will make a Type I error. That’s, um, quite a bit bigger than the 5% that it’s supposed to be. By way of comparison, imagine that you had used the following strategy. Start collecting data. Every single time an observation arrives, run a *Bayesian* \\(t\\)\-test (Section [17\.7](https://learningstatisticswithr.com/book/bayes.html#ttestbf) and look at the Bayes factor. I’ll assume that Johnson ([2013](https://learningstatisticswithr.com/book/bayes.html#ref-Johnson2013)) is right, and I’ll treat a Bayes factor of 3:1 as roughly equivalent to a \\(p\\)\-value of .05.[266](https://learningstatisticswithr.com/book/bayes.html#fn266) This time around, our trigger happy researcher uses the following procedure: if the Bayes factor is 3:1 or more in favour of the null, stop the experiment and retain the null. If it is 3:1 or more in favour of the alternative, stop the experiment and reject the null. Otherwise continue testing. Now, just like last time, let’s assume that the null hypothesis is true. What happens? As it happens, I ran the simulations for this scenario too, and the results are shown as the dashed line in Figure [17\.1](https://learningstatisticswithr.com/book/bayes.html#fig:type1). It turns out that the Type I error rate is much much lower than the 49% rate that we were getting by using the orthodox \\(t\\)\-test.
In some ways, this is remarkable. The entire *point* of orthodox null hypothesis testing is to control the Type I error rate. Bayesian methods aren’t actually designed to do this at all. Yet, as it turns out, when faced with a “trigger happy” researcher who keeps running hypothesis tests as the data come in, the Bayesian approach is much more effective. Even the 3:1 standard, which most Bayesians would consider unacceptably lax, is much safer than the \\(p\<.05\\) rule.
### 17\.5.1 Is it really this bad?
The example I gave in the previous section is a pretty extreme situation. In real life, people don’t run hypothesis tests every time a new observation arrives. So it’s not fair to say that the \\(p\<.05\\) threshold “really” corresponds to a 49% Type I error rate (i.e., \\(p=.49\\)). But the fact remains that if you want your \\(p\\)\-values to be honest, then you either have to switch to a completely different way of doing hypothesis tests, or you must enforce a strict rule: *no peeking*. You are *not* allowed to use the data to decide when to terminate the experiment. You are *not* allowed to look at a “borderline” \\(p\\)\-value and decide to collect more data. You aren’t even allowed to change your data analyis strategy after looking at data. You are strictly required to follow these rules, otherwise the \\(p\\)\-values you calculate will be nonsense.
And yes, these rules are surprisingly strict. As a class exercise a couple of years back, I asked students to think about this scenario. Suppose you started running your study with the intention of collecting \\(N=80\\) people. When the study starts out you follow the rules, refusing to look at the data or run any tests. But when you reach \\(N=50\\) your willpower gives in… and you take a peek. Guess what? You’ve got a significant result! Now, sure, you know you *said* that you’d keep running the study out to a sample size of \\(N=80\\), but it seems sort of pointless now, right? The result is significant with a sample size of \\(N=50\\), so wouldn’t it be wasteful and inefficient to keep collecting data? Aren’t you tempted to stop? Just a little? Well, keep in mind that if you do, your Type I error rate at \\(p\<.05\\) just ballooned out to 8%. When you report \\(p\<.05\\) in your paper, what you’re *really* saying is \\(p\<.08\\). That’s how bad the consequences of “just one peek” can be.
Now consider this … the scientific literature is filled with \\(t\\)\-tests, ANOVAs, regressions and chi-square tests. When I wrote this book I didn’t pick these tests arbitrarily. The reason why these four tools appear in most introductory statistics texts is that these are the bread and butter tools of science. None of these tools include a correction to deal with “data peeking”: they all assume that you’re not doing it. But how realistic is that assumption? In real life, how many people do you think have “peeked” at their data before the experiment was finished and adapted their subsequent behaviour after seeing what the data looked like? Except when the sampling procedure is fixed by an external constraint, I’m guessing the answer is “most people have done it”. If that has happened, you can infer that the reported \\(p\\)\-values are wrong. Worse yet, because we don’t know what decision process they actually followed, we have no way to know what the \\(p\\)\-values *should* have been. You can’t compute a \\(p\\)\-value when you don’t know the decision making procedure that the researcher used. And so the reported \\(p\\)\-value remains a lie.
Given all of the above, what is the take home message? It’s not that Bayesian methods are foolproof. If a researcher is determined to cheat, they can always do so. Bayes’ rule cannot stop people from lying, nor can it stop them from rigging an experiment. That’s not my point here. My point is the same one I made at the very beginning of the book in Section [1\.1](https://learningstatisticswithr.com/book/why-do-we-learn-statistics.html#whywhywhy): the reason why we run statistical tests is to protect us from ourselves. And the reason why “data peeking” is such a concern is that it’s so tempting, *even for honest researchers*. A theory for statistical inference has to acknowledge this. Yes, you might try to defend \\(p\\)\-values by saying that it’s the fault of the researcher for not using them properly. But to my mind that misses the point. A theory of statistical inference that is so completely naive about humans that it doesn’t even consider the possibility that the researcher might *look at their own data* isn’t a theory worth having. In essence, my point is this:
> *Good laws have their origins in bad morals.*
> – Ambrosius Macrobius[267](https://learningstatisticswithr.com/book/bayes.html#fn267)
Good rules for statistical testing have to acknowledge human frailty. None of us are without sin. None of us are beyond temptation. A good system for statistical inference should still work even when it is used by actual humans. Orthodox null hypothesis testing does not.[268](https://learningstatisticswithr.com/book/bayes.html#fn268)
## 17\.6 Bayesian analysis of contingency tables
Time to change gears. Up to this point I’ve been talking about what Bayesian inference is and why you might consider using it. I now want to briefly describe how to do Bayesian versions of various statistical tests. The discussions in the next few sections are not as detailed as I’d like, but I hope they’re enough to help you get started. So let’s begin.
The first kind of statistical inference problem I discussed in this book appeared in Chapter [12](https://learningstatisticswithr.com/book/chisquare.html#chisquare), in which we discussed categorical data analysis problems. In that chapter I talked about several different statistical problems that you might be interested in, but the one that appears most often in real life is the analysis of *contingency tables*. In this kind of data analysis situation, we have a cross-tabulation of one variable against another one, and the goal is to find out if there is some *association* between these variables. The data set I used to illustrate this problem is found in the `chapek9.Rdata` file, and it contains a single data frame `chapek9`
```
load(file.path(projecthome, "data","chapek9.Rdata"))
head(chapek9)
```
```
## species choice
## 1 robot flower
## 2 human data
## 3 human data
## 4 human data
## 5 robot data
## 6 human flower
```
In this data set, we supposedly sampled 180 beings and measured two things. First, we checked whether they were humans or robots, as captured by the `species` variable. Second, we asked them to nominate whether they most preferred flowers, puppies, or data. When we produce the cross-tabulation, we get this as the results:
```
crosstab <- xtabs( ~ species + choice, chapek9 )
crosstab
```
```
## choice
## species puppy flower data
## robot 13 30 44
## human 15 13 65
```
Surprisingly, the humans seemed to show a much stronger preference for data than the robots did. At the time we speculated that this might have been because the questioner was a large robot carrying a gun, and the humans might have been scared.
### 17\.6.1 The orthodox text
Just to refresh your memory, here’s how we analysed these data back in [Chapter@refch](mailto:Chapter@refch):chisquare. Because we want to determine if there is some *association* between `species` and `choice`, we used the `associationTest()` function in the `lsr` package to run a chi-square test of association. The results looked like this:
```
library(lsr)
associationTest( ~species + choice, chapek9 )
```
```
##
## Chi-square test of categorical association
##
## Variables: species, choice
##
## Hypotheses:
## null: variables are independent of one another
## alternative: some contingency exists between variables
##
## Observed contingency table:
## choice
## species puppy flower data
## robot 13 30 44
## human 15 13 65
##
## Expected contingency table under the null hypothesis:
## choice
## species puppy flower data
## robot 13.5 20.8 52.7
## human 14.5 22.2 56.3
##
## Test results:
## X-squared statistic: 10.722
## degrees of freedom: 2
## p-value: 0.005
##
## Other information:
## estimated effect size (Cramer's v): 0.244
```
Because we found a small \\(p\\) value (in this case \\(p\<.01\\)), we concluded that the data are inconsistent with the null hypothesis of no association, and we rejected it.
### 17\.6.2 The Bayesian test
How do we run an equivalent test as a Bayesian? Well, like every other bloody thing in statistics, there’s a lot of different ways you *could* do it. However, for the sake of everyone’s sanity, throughout this chapter I’ve decided to rely on one R package to do the work. Specifically, I’m going to use the `BayesFactor` package written by Jeff Rouder and Rich Morey, which as of this writing is in version 0.9.10.
For the analysis of contingency tables, the `BayesFactor` package contains a function called `contingencyTableBF()`. The data that you need to give to this function is the contingency table itself (i.e., the `crosstab` variable above), so you might be expecting to use a command like this:
```
library( BayesFactor ) # ...because we have to load the package
contingencyTableBF( crosstab ) # ...because that makes sense, right?
```
However, if you try this you’ll get an error message. This is because the `contingencyTestBF()` function needs one other piece of information from you: it needs to know what *sampling plan* you used to run your experiment. You can specify the sampling plan using the `sampleType` argument. So I should probably tell you what your options are! The `contingencyTableBF()` function distinguishes between four different types of experiment:
- **Fixed sample size**. Suppose that in our `chapek9` example, our experiment was designed like this: we deliberately set out to test 180 people, but we didn’t try to control the number of humans or robots, nor did we try to control the choices they made. In this design, the total number of observations \\(N\\) is fixed, but everything else is random. This is referred to as “joint multinomial” sampling, and if that’s what you did you should specify `sampleType = "jointMulti"`. In the case of the `chapek9` data, that’s actually what I had in mind when I invented the data set.
- **Fixed row (or column) totals**. A different kind of design might work like this. We decide ahead of time that we want 180 people, but we try to be a little more systematic about it. Specifically, the *experimenter* constrains it so that we get a predetermined number of humans and robots (e.g., 90 of each). In this design, *either* the row totals or the column totals are fixed, but not both. This is referred to as “independent multinomial” sampling, and if that’s what you did you should specify `sampleType = "indepMulti"`.
- **Both row and column totals fixed**. Another logical possibility is that you designed the experiment so that *both* the row totals and the column totals are fixed. This doesn’t make any sense at all in the `chapek9` example, but there are other deisgns that can work this way. Suppose that I show you a collection of 20 toys, and then given them 10 stickers that say `boy` and another 10 that say `girl`. I then give them 10 `blue` stickers and 10 `pink` stickers. I then ask you to put the stickers on the 20 toys such that every toy has a colour and every toy has a gender. No matter how you assign the stickers, the total number of pink and blue toys will be 10, as will the number of boys and girls. In this design *both* the rows and columns of the contingency table are fixed. This is referred to as “hypergeometric” sampling, and if that’s what you’ve done you should specify `sampleType = "hypergeom"`.
- **Nothing is fixed**. Finally, it might be the case that *nothing* is fixed. Not the row columns, not the column totals, and not the total sample size either. For instance, in the `chapek9` scenario, suppose what I’d done is run the study for a fixed length of *time*. By chance, it turned out that I got 180 people to turn up to study, but it could easily have been something else. This is referred to as “Poisson” sampling, and if that’s what you’ve done you should specify `sampleType="poisson"`.
Okay, so now we have enough knowledge to actually run a test. For the `chapek9` data, I implied that we designed the study such that the total sample size \\(N\\) was fixed, so we should set `sampleType = "jointMulti"`. The command that we need is,
```
library( BayesFactor )
contingencyTableBF( crosstab, sampleType = "jointMulti" )
```
```
## Bayes factor analysis
## --------------
## [1] Non-indep. (a=1) : 15.92684 ±0%
##
## Against denominator:
## Null, independence, a = 1
## ---
## Bayes factor type: BFcontingencyTable, joint multinomial
```
As with most R commands, the output initially looks suspiciously similar to utter gibberish. Fortunately, it’s actually pretty simple once you get past the initial impression. Firstly, note that the stuff at the top and bottom are irrelevant fluff. You already know that you’re doing a Bayes factor analysis. You already know that you’re analysing a contingency table, and you already know that you specified a joint multinomial sampling plan. So let’s strip that out and take a look at what’s left over:
```
[1] Non-indep. (a=1) : 15.92684 @plusorminus0%
Against denominator:
Null, independence, a = 1
```
Let’s also ignore those two `a=1` bits, since they’re technical details that you don’t need to know about at this stage.[269](https://learningstatisticswithr.com/book/bayes.html#fn269) The rest of the output is actually pretty straightforward. At the bottom, the output defines the null hypothesis for you: in this case, the null hypothesis is that there is no relationship between `species` and `choice`. Or, to put it another way, the null hypothesis is that these two variables are *independent*. Now if you look at the line above it, you might (correctly) guess that the `Non-indep.` part refers to the *alternative* hypothesis. In this case, the alternative is that there *is* a relationship between `species` and `choice`: that is, they are not independent. So the only thing left in the output is the bit that reads
```
15.92684 @plusorminus0%
```
The 15.9 part is the Bayes factor, and it’s telling you that the odds for the alternative hypothesis against the null are about 16:1. The \\(\\pm0\\%\\) part is not very interesting: essentially, all it’s telling you is that R has calculated an exact Bayes factor, so the uncertainty about the Bayes factor is 0%.[270](https://learningstatisticswithr.com/book/bayes.html#fn270) In any case, the data are telling us that we have moderate evidence for the alternative hypothesis.
### 17\.6.3 Writing up the results
When writing up the results, my experience has been that there aren’t quite so many “rules” for how you “should” report Bayesian hypothesis tests. That might change in the future if Bayesian methods become standard and some task force starts writing up style guides, but in the meantime I would suggest using some common sense. For example, I would avoid writing this:
> A Bayesian test of association found a significant result (BF=15.92)
To my mind, this write up is unclear. Even assuming that you’ve already reported the relevant descriptive statistics, there are a number of things I am unhappy with. First, the concept of “statistical significance” is pretty closely tied with \\(p\\)\-values, so it reads slightly strangely. Second, the “BF=15.92” part will only make sense to people who already understand Bayesian methods, and not everyone does. Third, it is somewhat unclear exactly which test was run and what software was used to do so.
On the other hand, unless precision is *extremely* important, I think that this is taking things a step too far:
> We ran a Bayesian test of association using version 0.9.10-1 of the BayesFactor package using default priors and a joint multinomial sampling plan. The resulting Bayes factor of 15.92 to 1 in favour of the alternative hypothesis indicates that there is moderately strong evidence for the non-independence of species and choice.
Everything about that passage is correct, of course. Morey and Rouder ([2015](https://learningstatisticswithr.com/book/bayes.html#ref-Morey2015)) built their Bayesian tests of association using the paper by Gunel and Dickey ([1974](https://learningstatisticswithr.com/book/bayes.html#ref-Gunel1974)), the specific test we used assumes that the experiment relied on a joint multinomial sampling plan, and indeed the Bayes factor of 15.92 is moderately strong evidence. It’s just far too wordy.
In most situations you just don’t need that much information. My preference is usually to go for something a little briefer. First, if you’re reporting multiple Bayes factor analyses in your write up, then somewhere you only need to cite the software once, at the beginning of the results section. So you might have one sentence like this:
> All analyses were conducted using the BayesFactor package in R , and unless otherwise stated default parameter values were used
Notice that I don’t bother including the version number? That’s because the citation itself includes that information (go check my reference list if you don’t believe me). There’s no need to clutter up your results with redundant information that almost no-one will actually need. When you get to the actual test you can get away with this:
> A test of association produced a Bayes factor of 16:1 in favour of a relationship between species and choice.
Short and sweet. I’ve rounded 15.92 to 16, because there’s not really any important difference between 15.92:1 and 16:1. I spelled out “Bayes factor” rather than truncating it to “BF” because not everyone knows the abbreviation. I indicated exactly what the effect is (i.e., “a relationship between species and choice”) and how strong the evidence was. I *didn’t* bother indicating whether this was “moderate” evidence or “strong” evidence, because the odds themselves tell you! There’s nothing stopping you from including that information, and I’ve done so myself on occasions, but you don’t strictly need it. Similarly, I didn’t bother to indicate that I ran the “joint multinomial” sampling plan, because I’m assuming that the method section of my write up would make clear how the experiment was designed. (I might change my mind about that if the method section was ambiguous.) Neither did I bother indicating that this was a *Bayesian* test of association: if your reader can’t work that out from the fact that you’re reporting a Bayes factor and the fact that you’re citing the `BayesFactor` package for all your analyses, then there’s no chance they’ll understand anything you’ve written. Besides, if you keep writing the word “Bayes” over and over again it starts to look stupid. Bayes Bayes Bayes Bayes Bayes. See?
### 17\.6.4 Other sampling plans
Up to this point all I’ve shown you is how to use the `contingencyTableBF()` function for the joint multinomial sampling plan (i.e., when the total sample size \\(N\\) is fixed, but nothing else is). For the Poisson sampling plan (i.e., nothing fixed), the command you need is identical except for the `sampleType` argument:
```
contingencyTableBF(crosstab, sampleType = "poisson" )
```
```
## Bayes factor analysis
## --------------
## [1] Non-indep. (a=1) : 28.20757 ±0%
##
## Against denominator:
## Null, independence, a = 1
## ---
## Bayes factor type: BFcontingencyTable, poisson
```
Notice that the Bayes factor of 28:1 here is *not* the identical to the Bayes factor of 16:1 that we obtained from the last test. The sampling plan actually does matter.
What about the design in which the row columns (or column totals) are fixed? As I mentioned earlier, this corresponds to the “independent multinomial” sampling plan. Again, you need to specify the `sampleType` argument, but this time you need to specify whether you fixed the rows or the columns. For example, suppose I deliberately sampled 87 humans and 93 robots, then I would need to indicate that the `fixedMargin` of the contingency table is the `"rows"`. So the command I would use is:
```
contingencyTableBF(crosstab, sampleType = "indepMulti", fixedMargin="rows")
```
```
## Bayes factor analysis
## --------------
## [1] Non-indep. (a=1) : 8.605897 ±0%
##
## Against denominator:
## Null, independence, a = 1
## ---
## Bayes factor type: BFcontingencyTable, independent multinomial
```
Again, the Bayes factor is different, with the evidence for the alternative dropping to a mere 9:1. As you might expect, the answers would be diffrent again if it were the columns of the contingency table that the experimental design fixed.
Finally, if we turn to hypergeometric sampling in which everything is fixed, we get…
```
contingencyTableBF(crosstab, sampleType = "hypergeom")
#Error in contingencyHypergeometric(as.matrix(data2), a) :
# hypergeometric contingency tables restricted to 2 x 2 tables; see help for contingencyTableBF()
```
… an error message. Okay, some quick reading through the help files hints that support for larger contingency tables is coming, but it’s not been implemented yet. In the meantime, let’s imagine we have data from the “toy labelling” experiment I described earlier in this section. Specifically, let’s say our data look like this:
```
toys <- data.frame(stringsAsFactors=FALSE,
gender = c("girl", "boy"),
pink = c(8, 2),
blue = c(2, 8)
)
```
The Bayesian test with hypergeometric sampling gives us this:
```
contingencyTableBF(toys, sampleType = "hypergeom")
#Bayes factor analysis
#--------------
#[1] Non-indep. (a=1) : 8.294321 @plusorminus0%
#
#Against denominator:
# Null, independence, a = 1
#---
#Bayes factor type: BFcontingencyTable, hypergeometric
```
The Bayes factor of 8:1 provides modest evidence that the labels were being assigned in a way that correlates gender with colour, but it’s not conclusive.
## 17\.7 Bayesian \\(t\\)\-tests
The second type of statistical inference problem discussed in this book is the comparison between two means, discussed in some detail in the chapter on \\(t\\)\-tests (Chapter [13](https://learningstatisticswithr.com/book/ttest.html#ttest). If you can remember back that far, you’ll recall that there are several versions of the \\(t\\)\-test. The `BayesFactor` package contains a function called `ttestBF()` that is flexible enough to run several different versions of the \\(t\\)\-test. I’ll talk a little about Bayesian versions of the independent samples \\(t\\)\-tests and the paired samples \\(t\\)\-test in this section.
### 17\.7.1 Independent samples \\(t\\)\-test
The most common type of \\(t\\)\-test is the independent samples \\(t\\)\-test, and it arises when you have data that look something like this:
```
load(file.path(projecthome, "data","harpo.Rdata"))
head(harpo)
```
```
## grade tutor
## 1 65 Anastasia
## 2 72 Bernadette
## 3 66 Bernadette
## 4 74 Anastasia
## 5 73 Anastasia
## 6 71 Bernadette
```
In this data set, we have two groups of students, those who received lessons from Anastasia and those who took their classes with Bernadette. The question we want to answer is whether there’s any difference in the grades received by these two groups of student. Back in [Chapter@refch](mailto:Chapter@refch):ttest I suggested you could analyse this kind of data using the `independentSamplesTTest()` function in the `lsr` package. For example, if you want to run a Student’s \\(t\\)\-test, you’d use a command like this:
```
independentSamplesTTest(
formula = grade ~ tutor,
data = harpo,
var.equal = TRUE
)
```
```
##
## Student's independent samples t-test
##
## Outcome variable: grade
## Grouping variable: tutor
##
## Descriptive statistics:
## Anastasia Bernadette
## mean 74.533 69.056
## std dev. 8.999 5.775
##
## Hypotheses:
## null: population means equal for both groups
## alternative: different population means in each group
##
## Test results:
## t-statistic: 2.115
## degrees of freedom: 31
## p-value: 0.043
##
## Other information:
## two-sided 95% confidence interval: [0.197, 10.759]
## estimated effect size (Cohen's d): 0.74
```
Like most of the functions that I wrote for this book, the `independentSamplesTTest()` is very wordy. It prints out a bunch of descriptive statistics and a reminder of what the null and alternative hypotheses are, before finally getting to the test results. I wrote it that way deliberately, in order to help make things a little clearer for people who are new to statistics.
Again, we obtain a \\(p\\)\-value less than 0.05, so we reject the null hypothesis.
What does the Bayesian version of the \\(t\\)\-test look like? Using the `ttestBF()` function, we can obtain a Bayesian analog of Student’s independent samples \\(t\\)\-test using the following command:
```
ttestBF( formula = grade ~ tutor, data = harpo )
```
```
## Bayes factor analysis
## --------------
## [1] Alt., r=0.707 : 1.754927 ±0%
##
## Against denominator:
## Null, mu1-mu2 = 0
## ---
## Bayes factor type: BFindepSample, JZS
```
Notice that format of this command is pretty standard. As usual we have a `formula` argument in which we specify the outcome variable on the left hand side and the grouping variable on the right. The `data` argument is used to specify the data frame containing the variables. However, notice that there’s no analog of the `var.equal` argument. This is because the `BayesFactor` package does not include an analog of the Welch test, only the Student test.[271](https://learningstatisticswithr.com/book/bayes.html#fn271) In any case, when you run this command you get this as the output:
So what does all this mean? Just as we saw with the `contingencyTableBF()` function, the output is pretty dense. But, just like last time, there’s not a lot of information here that you actually need to process. Firstly, let’s examine the bottom line. The `BFindepSample` part just tells you that you ran an independent samples \\(t\\)\-test, and the `JZS` part is technical information that is a little beyond the scope of this book.[272](https://learningstatisticswithr.com/book/bayes.html#fn272) Clearly, there’s nothing to worry about in that part. In the line above, the text `Null, mu1-mu2 = 0` is just telling you that the null hypothesis is that there are no differences between means. But you already knew that. So the only part that really matters is this line here:
```
[1] Alt., r=0.707 : 1.754927 @plusorminus0%
```
Ignore the `r=0.707` part: it refers to a technical detail that we won’t worry about in this chapter.[273](https://learningstatisticswithr.com/book/bayes.html#fn273) Instead, you should focus on the part that reads `1.754927`. This is the Bayes factor: the evidence provided by these data are about 1.8:1 in favour of the alternative.
Before moving on, it’s worth highlighting the difference between the orthodox test results and the Bayesian one. According to the orthodox test, we obtained a significant result, though only barely. Nevertheless, many people would happily accept \\(p=.043\\) as reasonably strong evidence for an effect. In contrast, notice that the Bayesian test doesn’t even reach 2:1 odds in favour of an effect, and would be considered very weak evidence at best. In my experience that’s a pretty typical outcome. Bayesian methods usually require more evidence before rejecting the null.
### 17\.7.2 Paired samples \\(t\\)\-test
Back in Section [13\.5](https://learningstatisticswithr.com/book/ttest.html#pairedsamplesttest) I discussed the `chico` data frame in which students grades were measured on two tests, and we were interested in finding out whether grades went up from test 1 to test 2. Because every student did both tests, the tool we used to analyse the data was a paired samples \\(t\\)\-test. To remind you of what the data look like, here’s the first few cases:
```
load(file.path(projecthome, "data","chico.Rdata"))
head(chico)
```
```
## id grade_test1 grade_test2
## 1 student1 42.9 44.6
## 2 student2 51.8 54.0
## 3 student3 71.7 72.3
## 4 student4 51.6 53.4
## 5 student5 63.5 63.8
## 6 student6 58.0 59.3
```
We originally analysed the data using the `pairedSamplesTTest()` function in the `lsr` package, but this time we’ll use the `ttestBF()` function from the `BayesFactor` package to do the same thing. The easiest way to do it with this data set is to use the `x` argument to specify one variable and the `y` argument to specify the other. All we need to do then is specify `paired=TRUE` to tell R that this is a paired samples test. So here’s our command:
```
ttestBF(
x = chico$grade_test1,
y = chico$grade_test2,
paired = TRUE
)
```
```
## Bayes factor analysis
## --------------
## [1] Alt., r=0.707 : 5992.05 ±0%
##
## Against denominator:
## Null, mu = 0
## ---
## Bayes factor type: BFoneSample, JZS
```
At this point, I hope you can read this output without any difficulty. The data provide evidence of about 6000:1 in favour of the alternative. We could probably reject the null with some confidence\!
## 17\.8 Bayesian regression
Okay, so now we’ve seen Bayesian equivalents to orthodox chi-square tests and \\(t\\)\-tests. What’s next? If I were to follow the same progression that I used when developing the orthodox tests you’d expect to see ANOVA next, but I think it’s a little clearer if we start with regression.
### 17\.8.1 A quick refresher
In Chapter [15](https://learningstatisticswithr.com/book/regression.html#regression) I used the `parenthood` data to illustrate the basic ideas behind regression. To remind you of what that data set looks like, here’s the first six observations:
```
load(file.path(projecthome, "data","parenthood.Rdata"))
head(parenthood)
```
```
## dan.sleep baby.sleep dan.grump day
## 1 7.59 10.18 56 1
## 2 7.91 11.66 60 2
## 3 5.14 7.92 82 3
## 4 7.71 9.61 55 4
## 5 6.68 9.75 67 5
## 6 5.99 5.04 72 6
```
Back in Chapter [15](https://learningstatisticswithr.com/book/regression.html#regression) I proposed a theory in which my grumpiness (`dan.grump`) on any given day is related to the amount of sleep I got the night before (`dan.sleep`), and possibly to the amount of sleep our baby got (`baby.sleep`), though probably not to the `day` on which we took the measurement. We tested this using a regression model. In order to estimate the regression model we used the `lm()` function, like so:
```
model <- lm(
formula = dan.grump ~ dan.sleep + day + baby.sleep,
data = parenthood
)
```
The hypothesis tests for each of the terms in the regression model were extracted using the `summary()` function as shown below:
```
summary(model)
```
```
##
## Call:
## lm(formula = dan.grump ~ dan.sleep + day + baby.sleep, data = parenthood)
##
## Residuals:
## Min 1Q Median 3Q Max
## -10.906 -2.284 -0.295 2.652 11.880
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 126.278707 3.242492 38.945 <2e-16 ***
## dan.sleep -8.969319 0.560007 -16.016 <2e-16 ***
## day -0.004403 0.015262 -0.288 0.774
## baby.sleep 0.015747 0.272955 0.058 0.954
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 4.375 on 96 degrees of freedom
## Multiple R-squared: 0.8163, Adjusted R-squared: 0.8105
## F-statistic: 142.2 on 3 and 96 DF, p-value: < 2.2e-16
```
When interpreting the results, each row in this table corresponds to one of the possible predictors. The `(Intercept)` term isn’t usually interesting, though it is highly significant. The important thing for our purposes is the fact that `dan.sleep` is significant at \\(p\<.001\\) and neither of the other variables are.
### 17\.8.2 The Bayesian version
Okay, so how do we do the same thing using the `BayesFactor` package? The easiest way is to use the `regressionBF()` function instead of `lm()`. As before, we use `formula` to indicate what the full regression model looks like, and the `data` argument to specify the data frame. So the command is:
```
regressionBF(
formula = dan.grump ~ dan.sleep + day + baby.sleep,
data = parenthood
)
```
```
##
|
| | 0%
|
|========= | 14%
|
|=================== | 29%
|
|============================ | 43%
|
|===================================== | 57%
|
|============================================== | 71%
|
|======================================================== | 86%
|
|=================================================================| 100%
```
```
## Bayes factor analysis
## --------------
## [1] dan.sleep : 1.622545e+34 ±0.01%
## [2] day : 0.2724027 ±0%
## [3] baby.sleep : 10018411 ±0%
## [4] dan.sleep + day : 1.016576e+33 ±0%
## [5] dan.sleep + baby.sleep : 9.77022e+32 ±0%
## [6] day + baby.sleep : 2340755 ±0%
## [7] dan.sleep + day + baby.sleep : 7.835625e+31 ±0%
##
## Against denominator:
## Intercept only
## ---
## Bayes factor type: BFlinearModel, JZS
```
So that’s pretty straightforward: it’s exactly what we’ve been doing throughout the book. The output, however, is a little different from what you get from `lm()`. The format of this is pretty familiar. At the bottom we have some techical rubbish, and at the top we have some information about the Bayes factors. What’s new is the fact that we seem to have *lots* of Bayes factors here. What’s all this about?
The trick to understanding this output is to recognise that if we’re interested in working out which of the 3 predictor variables are related to `dan.grump`, there are actually 8 possible regression models that could be considered. One possibility is the *intercept only model*, in which none of the three variables have an effect. At the other end of the spectrum is the *full model* in which all three variables matter. So what `regressionBF()` does is treat the *intercept only* model as the null hypothesis, and print out the Bayes factors for all other models when compared against that null. For example, if we look at line 4 in the table, we see that the evidence is about \\(10^{33}\\) to 1 in favour of the claim that a model that includes both `dan.sleep` and `day` is better than the intercept only model. Or if we look at line 1, we can see that the odds are about \\(1.6 \\times 10^{34}\\) that a model containing the `dan.sleep` variable (but no others) is better than the intercept only model.
### 17\.8.3 Finding the best model
In practice, this isn’t super helpful. In most situations the intercept only model is one that you don’t really care about at all. What I find helpful is to start out by working out which model is the *best* one, and then seeing how well all the alternatives compare to it. Here’s how you do that. In this case, it’s easy enough to see that the best model is actually the one that contains `dan.sleep` only (line 1), because it has the largest Bayes factor. However, if you’ve got a lot of possible models in the output, it’s handy to know that you can use the `head()` function to pick out the best few models. First, we have to go back and save the Bayes factor information to a variable:
```
models <- regressionBF(
formula = dan.grump ~ dan.sleep + day + baby.sleep,
data = parenthood
)
```
```
##
|
| | 0%
|
|========= | 14%
|
|=================== | 29%
|
|============================ | 43%
|
|===================================== | 57%
|
|============================================== | 71%
|
|======================================================== | 86%
|
|=================================================================| 100%
```
Let’s say I want to see the best three models. To do this, I use the `head()` function specifying `n=3`, and here’s what I get as the result:
```
head( models, n = 3)
```
```
## Bayes factor analysis
## --------------
## [1] dan.sleep : 1.622545e+34 ±0.01%
## [2] dan.sleep + day : 1.016576e+33 ±0%
## [3] dan.sleep + baby.sleep : 9.77022e+32 ±0%
##
## Against denominator:
## Intercept only
## ---
## Bayes factor type: BFlinearModel, JZS
```
This is telling us that the model in line 1 (i.e., `dan.grump ~ dan.sleep`) is the best one. That’s *almost* what I’m looking for, but it’s still comparing all the models against the intercept only model. That seems silly. What I’d like to know is how big the difference is between the best model and the other good models. For that, there’s this trick:
```
head( models/max(models), n = 3)
```
```
## Bayes factor analysis
## --------------
## [1] dan.sleep : 1 ±0%
## [2] dan.sleep + day : 0.0626532 ±0.01%
## [3] dan.sleep + baby.sleep : 0.0602154 ±0.01%
##
## Against denominator:
## dan.grump ~ dan.sleep
## ---
## Bayes factor type: BFlinearModel, JZS
```
Notice the bit at the bottom showing that the “denominator” has changed. What that means is that the Bayes factors are now comparing each of those 3 models listed against the `dan.grump ~ dan.sleep` model. Obviously, the Bayes factor in the first line is exactly 1, since that’s just comparing the best model to itself. More to the point, the other two Bayes factors are both less than 1, indicating that they’re all worse than that model. The Bayes factors of 0.06 to 1 imply that the odds for the best model over the second best model are about 16:1. You can work this out by simple arithmetic (i.e., \\(0.06 / 1 \\approx 16\\)), but the other way to do it is to directly compare the models. To see what I mean, here’s the original output:
```
models
```
```
## Bayes factor analysis
## --------------
## [1] dan.sleep : 1.622545e+34 ±0.01%
## [2] day : 0.2724027 ±0%
## [3] baby.sleep : 10018411 ±0%
## [4] dan.sleep + day : 1.016576e+33 ±0%
## [5] dan.sleep + baby.sleep : 9.77022e+32 ±0%
## [6] day + baby.sleep : 2340755 ±0%
## [7] dan.sleep + day + baby.sleep : 7.835625e+31 ±0%
##
## Against denominator:
## Intercept only
## ---
## Bayes factor type: BFlinearModel, JZS
```
The best model corresponds to row 1 in this table, and the second best model corresponds to row 4. All you have to do to compare these two models is this:
```
models[1] / models[4]
```
```
## Bayes factor analysis
## --------------
## [1] dan.sleep : 15.96088 ±0.01%
##
## Against denominator:
## dan.grump ~ dan.sleep + day
## ---
## Bayes factor type: BFlinearModel, JZS
```
And there you have it. You’ve found the regression model with the highest Bayes factor (i.e., `dan.grump ~ dan.sleep`), and you know that the evidence for that model over the next best alternative (i.e., `dan.grump ~ dan.sleep + day`) is about 16:1.
### 17\.8.4 Extracting Bayes factors for all included terms
Okay, let’s say you’ve settled on a specific regression model. What Bayes factors should you report? In this example, I’m going to pretend that you decided that `dan.grump ~ dan.sleep + baby.sleep` is the model you think is best. Sometimes it’s sensible to do this, even when it’s not the one with the highest Bayes factor. Usually this happens because you have a substantive theoretical reason to prefer one model over the other. However, in this case I’m doing it because I want to use a model with more than one predictor as my example\!
Having figured out which model you prefer, it can be really useful to call the `regressionBF()` function and specifying `whichModels="top"`. You use your “preferred” model as the `formula` argument, and then the output will show you the Bayes factors that result when you try to drop predictors from this model:
```
regressionBF(
formula = dan.grump ~ dan.sleep + baby.sleep,
data = parenthood,
whichModels = "top"
)
```
```
##
|
| | 0%
|
|====================== | 33%
|
|=========================================== | 67%
|
|=================================================================| 100%
```
```
## Bayes factor top-down analysis
## --------------
## When effect is omitted from dan.sleep + baby.sleep , BF is...
## [1] Omit baby.sleep : 16.60705 ±0.01%
## [2] Omit dan.sleep : 1.025403e-26 ±0.01%
##
## Against denominator:
## dan.grump ~ dan.sleep + baby.sleep
## ---
## Bayes factor type: BFlinearModel, JZS
```
Okay, so now you can see the results a bit more clearly. The Bayes factor when you try to drop the `dan.sleep` predictor is about \\(10^{-26}\\), which is very strong evidence that you *shouldn’t* drop it. On the other hand, the Bayes factor actually goes up to 17 if you drop `baby.sleep`, so you’d usually say that’s pretty strong evidence for dropping that one.
## 17\.9 Bayesian ANOVA
As you can tell, the `BayesFactor` package is pretty flexible, and it can do Bayesian versions of pretty much everything in this book. In fact, it can do a few other neat things that I haven’t covered in the book at all. However, I have to stop somewhere, and so there’s only one other topic I want to cover: Bayesian ANOVA.
### 17\.9.1 A quick refresher
As with the other examples, I think it’s useful to start with a reminder of how I discussed ANOVA earlier in the book. First, let’s remind ourselves of what the data were. The example I used originally is the `clin.trial` data frame, which looks like this
```
load(file.path(projecthome, "data","clinicaltrial.Rdata"))
head(clin.trial)
```
```
## drug therapy mood.gain
## 1 placebo no.therapy 0.5
## 2 placebo no.therapy 0.3
## 3 placebo no.therapy 0.1
## 4 anxifree no.therapy 0.6
## 5 anxifree no.therapy 0.4
## 6 anxifree no.therapy 0.2
```
To run our orthodox analysis in earlier chapters we used the `aov()` function to do all the heavy lifting. In Chapter [16](https://learningstatisticswithr.com/book/anova2.html#anova2) I recommended using the `Anova()` function from the `car` package to produce the ANOVA table, because it uses Type II tests by default. If you’ve forgotten what “Type II tests” are, it might be a good idea to re-read Section [16\.10](https://learningstatisticswithr.com/book/anova2.html#unbalancedanova), because it will become relevant again in a moment. In any case, here’s what our analysis looked like:
```
library(car)
model <- aov( mood.gain ~ drug * therapy, data = clin.trial )
Anova(model)
```
```
## Anova Table (Type II tests)
##
## Response: mood.gain
## Sum Sq Df F value Pr(>F)
## drug 3.4533 2 31.7143 1.621e-05 ***
## therapy 0.4672 1 8.5816 0.01262 *
## drug:therapy 0.2711 2 2.4898 0.12460
## Residuals 0.6533 12
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
```
That’s pretty clearly showing us evidence for a main effect of `drug` at \\(p\<.001\\), an effect of `therapy` at \\(p\<.05\\) and no interaction.
### 17\.9.2 The Bayesian version
How do we do the same thing using Bayesian methods? The `BayesFactor` package contains a function called `anovaBF()` that does this for you. It uses a pretty standard `formula` and `data` structure, so the command should look really familiar. Just like we did with regression, it will be useful to save the output to a variable:
```
models <- anovaBF(
formula = mood.gain ~ drug * therapy,
data = clin.trial
)
```
```
##
|
| | 0%
|
|================ | 25%
|
|================================ | 50%
|
|================================================= | 75%
|
|=================================================================| 100%
```
The output is quite different to the traditional ANOVA, but it’s not too bad once you understand what you’re looking for. Let’s take a look:
```
models
```
This looks very similar to the output we obtained from the `regressionBF()` function, and with good reason. Remember what I said back in Section [16\.6](https://learningstatisticswithr.com/book/anova2.html#anovalm): under the hood, ANOVA is no different to regression, and both are just different examples of a linear model. Becasue of this, the `anovaBF()` reports the output in much the same way. For instance, if we want to identify the best model we could use the same commands that we used in the last section. One variant that I find quite useful is this:
```
models/max(models)
```
```
## Bayes factor analysis
## --------------
## [1] drug : 0.3436062 ±1.34%
## [2] therapy : 0.001022285 ±1.34%
## [3] drug + therapy : 1 ±0%
## [4] drug + therapy + drug:therapy : 0.954423 ±1.64%
##
## Against denominator:
## mood.gain ~ drug + therapy
## ---
## Bayes factor type: BFlinearModel, JZS
```
By “dividing” the `models` output by the best model (i.e., `max(models)`), what R is doing is using the best model (which in this case is `drugs + therapy`) as the denominator, which gives you a pretty good sense of how close the competitors are. For instance, the model that contains the interaction term is almost as good as the model without the interaction, since the Bayes factor is 0.98. In other words, the data do not clearly indicate whether there is or is not an interaction.
### 17\.9.3 Constructing Bayesian Type II tests
Okay, that’s all well and good, you might be thinking, but what do I report as the alternative to the \\(p\\)\-value? In the classical ANOVA table, you get a single \\(p\\)\-value for every predictor in the model, so you can talk about the significance of each effect. What’s the Bayesian analog of this?
It’s a good question, but the answer is tricky. Remember what I said in Section [16\.10](https://learningstatisticswithr.com/book/anova2.html#unbalancedanova) about ANOVA being complicated. Even in the classical version of ANOVA there are several different “things” that ANOVA might correspond to. Specifically, I discussed how you get different \\(p\\)\-values depending on whether you use Type I tests, Type II tests or Type III tests. To work out which Bayes factor is analogous to “the” \\(p\\)\-value in a classical ANOVA, you need to work out which version of ANOVA you want an analog for. For the purposes of this section, I’ll assume you want Type II tests, because those are the ones I think are most sensible in general. As I discussed back in Section [16\.10](https://learningstatisticswithr.com/book/anova2.html#unbalancedanova), Type II tests for a two-way ANOVA are reasonably straightforward, but if you have forgotten that section it wouldn’t be a bad idea to read it again before continuing.
Assuming you’ve had a refresher on Type II tests, let’s have a look at how to pull them from the Bayes factor table. Suppose we want to test the main effect of `drug`. The null hypothesis for this test corresponds to a model that includes an effect of `therapy`, but no effect of `drug`. The alternative hypothesis is the model that includes both. In other words, what we want is the Bayes factor corresponding to this comparison:
```
knitr::kable(tibble::tribble(
~V1, ~V2,
"Null model:", "`mood.gain ~ therapy`",
"Alternative model:", "`mood.gain ~ therapy + drug`"
), col.names = c("", ""))
```
| | |
|---|---|
| Null model: | `mood.gain ~ therapy` |
| Alternative model: | `mood.gain ~ therapy + drug` |
As it happens, we can read the answer to this straight off the table because it corresponds to a comparison between the model in line 2 of the table and the model in line 3: the Bayes factor in this case represents evidence *for* the null of 0.001 to 1. Or, more helpfully, the odds are about 1000 to 1 against the null.
The main effect of `therapy` can be calculated in much the same way. In this case, the null model is the one that contains only an effect of drug, and the alternative is the model that contains both. So the relevant comparison is between lines 2 and 1 in the table. The odds in favour of the null here are only 0.35 to 1. Again, I find it useful to frame things the other way around, so I’d refer to this as evidence of about 3 to 1 in favour of an effect of `therapy`.
Finally, in order to test an interaction effect, the null model here is one that contains both main effects but no interaction. The alternative model adds the interaction. That is:
```
knitr::kable(tibble::tribble(
~V1, ~V2,
"Null model:", "`mood.gain ~ drug + therapy`",
"Alternative model:", "`mood.gain ~ drug + therapy + drug:therapy`"
), col.names = c("", ""))
```
| | |
|---|---|
| Null model: | `mood.gain ~ drug + therapy` |
| Alternative model: | `mood.gain ~ drug + therapy + drug:therapy` |
If we look those two models up in the table, we see that this comparison is between the models on lines 3 and 4 of the table. The odds of 0.98 to 1 imply that these two models are fairly evenly matched.
You might be thinking that this is all pretty laborious, and I’ll concede that’s true. At some stage I might consider adding a function to the `lsr` package that would automate this process and construct something like a “Bayesian Type II ANOVA table” from the output of the `anovaBF()` function. However, I haven’t had time to do this yet, nor have I made up my mind about whether it’s really a good idea to do this. In the meantime, I thought I should show you the trick for how I do this in practice. The command that I use when I want to grab the right Bayes factors for a Type II ANOVA is this one:
```
max(models)/models
```
```
## denominator
## numerator drug therapy drug + therapy
## drug + therapy 2.910308 978.2007 1
## denominator
## numerator drug + therapy + drug:therapy
## drug + therapy 1.047753
```
The output isn’t quite so pretty as the last one, but the nice thing is that you can read off everything you need. The best model is `drug + therapy`, so all the other models are being compared to that. What’s the Bayes factor *for* the main effect of `drug`? The relevant null hypothesis is the one that contains only `therapy`, and the Bayes factor in question is 954:1. The main effect of `therapy` is weaker, and the evidence here is only 2.8:1. Finally, the evidence *against* an interaction is very weak, at 1.01:1.
Reading the results off this table is sort of counterintuitive, because you have to read off the answers from the “wrong” part of the table. For instance, the evidence for an effect of `drug` can be read from the column labelled `therapy`, which is pretty damned weird. To be fair to the authors of the package, I don’t think they ever intended for the `anovaBF()` function to be used this way. My understanding[274](https://learningstatisticswithr.com/book/bayes.html#fn274) is that their view is simply that you should find the best model and report that model: there’s no inherent reason why a Bayesian ANOVA should try to follow the exact same design as an orthodox ANOVA.[275](https://learningstatisticswithr.com/book/bayes.html#fn275)
In any case, if you know what you’re looking for, you can look at this table and then report the results of the Bayesian analysis in a way that is pretty closely analogous to how you’d report a regular Type II ANOVA. As I mentioned earlier, there’s still no convention on how to do that, but I usually go for something like this:
> A Bayesian Type II ANOVA found evidence for main effects of drug (Bayes factor: 954:1) and therapy (Bayes factor: 3:1), but no clear evidence for or against an interaction (Bayes factor: 1:1).
## 17\.10 Summary
The first half of this chapter was focused primarily on the theoretical underpinnings of Bayesian statistics. I introduced the mathematics for how Bayesian inference works (Section [17\.1](https://learningstatisticswithr.com/book/bayes.html#basicbayes)), and gave a very basic overview of how Bayesian hypothesis testing is typically done (Section [17\.2](https://learningstatisticswithr.com/book/bayes.html#bayesianhypothesistests)). Finally, I devoted some space to talking about why I think Bayesian methods are worth using (Section [17\.3](https://learningstatisticswithr.com/book/bayes.html#whybayes).
The second half of the chapter was a lot more practical, and focused on tools provided by the `BayesFactor` package. Specifically, I talked about using the `contingencyTableBF()` function to do Bayesian analogs of chi-square tests (Section [17\.6](https://learningstatisticswithr.com/book/bayes.html#bayescontingency), the `ttestBF()` function to do Bayesian \\(t\\)\-tests, (Section [17\.7](https://learningstatisticswithr.com/book/bayes.html#ttestbf)), the `regressionBF()` function to do Bayesian regressions, and finally the `anovaBF()` function for Bayesian ANOVA.
If you’re interested in learning more about the Bayesian approach, there are many good books you could look into. John Kruschke’s book *Doing Bayesian Data Analysis* is a pretty good place to start (Kruschke [2011](https://learningstatisticswithr.com/book/bayes.html#ref-Kruschke2011)), and is a nice mix of theory and practice. His approach is a little different to the “Bayes factor” approach that I’ve discussed here, so you won’t be covering the same ground. If you’re a cognitive psychologist, you might want to check out Michael Lee and E.J. Wagenmakers’ book *Bayesian Cognitive Modeling* (Lee and Wagenmakers [2014](https://learningstatisticswithr.com/book/bayes.html#ref-Lee2014)). I picked these two because I think they’re especially useful for people in my discipline, but there’s a lot of good books out there, so look around\!
### References
Jeffreys, Harold. 1961. *The Theory of Probability*. 3rd ed. Oxford.
Kass, Robert E., and Adrian E. Raftery. 1995. “Bayes Factors.” *Journal of the American Statistical Association* 90: 773–95.
Fisher, R. 1925. *Statistical Methods for Research Workers*. Edinburgh, UK: Oliver; Boyd.
Johnson, Valen E. 2013. “Revised Standards for Statistical Evidence.” *Proceedings of the National Academy of Sciences*, no. 48: 19313–7.
Morey, Richard D., and Jeffrey N. Rouder. 2015. *BayesFactor: Computation of Bayes Factors for Common Designs*. <http://CRAN.R-project.org/package=BayesFactor>.
Gunel, Erdogan, and James Dickey. 1974. “Bayes Factors for Independence in Contingency Tables.” *Biometrika*, 545–57.
Kruschke, J. K. 2011. *Doing Bayesian Data Analysis: A Tutorial with R and BUGS*. Burlington, MA: Academic Press.
Lee, Michael D, and Eric-Jan Wagenmakers. 2014. *Bayesian Cognitive Modeling: A Practical Course*. Cambridge University Press.
***
1. <http://en.wikiquote.org/wiki/David_Hume>[↩](https://learningstatisticswithr.com/book/bayes.html#fnref254)
2. <http://en.wikipedia.org/wiki/Climate_of_Adelaide>[↩](https://learningstatisticswithr.com/book/bayes.html#fnref255)
3. It’s a leap of faith, I know, but let’s run with it okay?[↩](https://learningstatisticswithr.com/book/bayes.html#fnref256)
4. Um. I hate to bring this up, but some statisticians would object to me using the word “likelihood” here. The problem is that the word “likelihood” has a very specific meaning in frequentist statistics, and it’s not quite the same as what it means in Bayesian statistics. As far as I can tell, Bayesians didn’t originally have any agreed upon name for the likelihood, and so it became common practice for people to use the frequentist terminology. This wouldn’t have been a problem, except for the fact that the way that Bayesians use the word turns out to be quite different to the way frequentists do. This isn’t the place for yet another lengthy history lesson, but to put it crudely: when a Bayesian says “*a* likelihood function” they’re usually referring one of the *rows* of the table. When a frequentist says the same thing, they’re referring to the same table, but to them “*a* likelihood function” almost always refers to one of the *columns*. This distinction matters in some contexts, but it’s not important for our purposes.[↩](https://learningstatisticswithr.com/book/bayes.html#fnref257)
5. If we were being a bit more sophisticated, we could extend the example to accommodate the possibility that I’m lying about the umbrella. But let’s keep things simple, shall we?[↩](https://learningstatisticswithr.com/book/bayes.html#fnref258)
6. You might notice that this equation is actually a restatement of the same basic rule I listed at the start of the last section. If you multiply both sides of the equation by \\(P(d)\\), then you get \\(P(d) P(h\| d) = P(d,h)\\), which is the rule for how joint probabilities are calculated. So I’m not actually introducing any “new” rules here, I’m just using the same rule in a different way.[↩](https://learningstatisticswithr.com/book/bayes.html#fnref259)
7. Obviously, this is a highly simplified story. All the complexity of real life Bayesian hypothesis testing comes down to how you calculate the likelihood \\(P(d\|h)\\) when the hypothesis \\(h\\) is a complex and vague thing. I’m not going to talk about those complexities in this book, but I do want to highlight that although this simple story is true as far as it goes, real life is messier than I’m able to cover in an introductory stats textbook.[↩](https://learningstatisticswithr.com/book/bayes.html#fnref260)
8. <http://www.imdb.com/title/tt0093779/quotes>. I should note in passing that I’m not the first person to use this quote to complain about frequentist methods. Rich Morey and colleagues had the idea first. I’m shamelessly stealing it because it’s such an awesome pull quote to use in this context and I refuse to miss any opportunity to quote *The Princess Bride*.[↩](https://learningstatisticswithr.com/book/bayes.html#fnref261)
9. <http://about.abc.net.au/reports-publications/appreciation-survey-summary-report-2013/>[↩](https://learningstatisticswithr.com/book/bayes.html#fnref262)
10. <http://knowyourmeme.com/memes/the-cake-is-a-lie>[↩](https://learningstatisticswithr.com/book/bayes.html#fnref263)
11. In the interests of being completely honest, I should acknowledge that not all orthodox statistical tests that rely on this silly assumption. There are a number of *sequential analysis* tools that are sometimes used in clinical trials and the like. These methods are built on the assumption that data are analysed as they arrive, and these tests aren’t horribly broken in the way I’m complaining about here. However, sequential analysis methods are constructed in a very different fashion to the “standard” version of null hypothesis testing. They don’t make it into any introductory textbooks, and they’re not very widely used in the psychological literature. The concern I’m raising here is valid for every single orthodox test I’ve presented so far, and for almost every test I’ve seen reported in the papers I read.[↩](https://learningstatisticswithr.com/book/bayes.html#fnref264)
12. A related problem: <http://xkcd.com/1478/>[↩](https://learningstatisticswithr.com/book/bayes.html#fnref265)
13. Some readers might wonder why I picked 3:1 rather than 5:1, given that Johnson ([2013](https://learningstatisticswithr.com/book/bayes.html#ref-Johnson2013)) suggests that \\(p=.05\\) lies somewhere in that range. I did so in order to be charitable to the \\(p\\)\-value. If I’d chosen a 5:1 Bayes factor instead, the results would look even better for the Bayesian approach.[↩](https://learningstatisticswithr.com/book/bayes.html#fnref266)
14. <http://www.quotationspage.com/quotes/Ambrosius_Macrobius/>[↩](https://learningstatisticswithr.com/book/bayes.html#fnref267)
15. Okay, I just *know* that some knowledgeable frequentists will read this and start complaining about this section. Look, I’m not dumb. I absolutely know that if you adopt a sequential analysis perspective you can avoid these errors within the orthodox framework. I also know that you can explictly design studies with interim analyses in mind. So yes, in one sense I’m attacking a “straw man” version of orthodox methods. However, the straw man that I’m attacking is the one that *is used by almost every single practitioner*. If it ever reaches the point where sequential methods become the norm among experimental psychologists and I’m no longer forced to read 20 extremely dubious ANOVAs a day, I promise I’ll rewrite this section and dial down the vitriol. But until that day arrives, I stand by my claim that *default* Bayes factor methods are much more robust in the face of data analysis practices as they exist in the real world. *Default* orthodox methods suck, and we all know it.[↩](https://learningstatisticswithr.com/book/bayes.html#fnref268)
16. If you’re desperate to know, you can find all the gory details in Gunel and Dickey (1974). However, that’s a pretty technical paper. The help documentation to the `contingencyTableBF()` gives this explanation: “the argument `priorConcentration` indexes the expected deviation from the null hypothesis under the alternative, and corresponds to Gunel and Dickey’s (1974) \\(a\\) parameter.” As I write this I’m about halfway through the Gunel and Dickey paper, and I agree that setting \\(a=1\\) is a pretty sensible default choice, since it corresponds to an assumption that you have very little *a priori* knowledge about the contingency table.[↩](https://learningstatisticswithr.com/book/bayes.html#fnref269)
17. In some of the later examples, you’ll see that this number is not always 0%. This is because the `BayesFactor` package often has to run some simulations to compute approximate Bayes factors. So the answers you get won’t always be identical when you run the command a second time. That’s why the output of these functions tells you what the margin for error is.[↩](https://learningstatisticswithr.com/book/bayes.html#fnref270)
18. Apparently this omission is deliberate. I have this vague recollection that I spoke to Jeff Rouder about this once, and his opinion was that when homogeneity of variance is violated the results of a \\(t\\)\-test are uninterpretable. I can see the argument for this, but I’ve never really held a strong opinion myself. (Jeff, if you never said that, I’m sorry)[↩](https://learningstatisticswithr.com/book/bayes.html#fnref271)
19. Just in case you’re interested: the “JZS” part of the output relates to how the Bayesian test expresses the prior uncertainty about the variance \\(\\sigma^2\\), and it’s short for the names of three people: “Jeffreys Zellner Siow”. See Rouder et al. ([2009](https://learningstatisticswithr.com/book/bayes.html#ref-Rouder2009)) for details.[↩](https://learningstatisticswithr.com/book/bayes.html#fnref272)
20. Again, in case you care … the null hypothesis here specifies an effect size of 0, since the two means are identical. The alternative hypothesis states that there *is* an effect, but it doesn’t specify exactly how big the effect will be. The \\(r\\) value here relates to how big the effect is expected to be according to the alternative. You can type `?ttestBF` to get more details.[↩](https://learningstatisticswithr.com/book/bayes.html#fnref273)
21. Again, guys, sorry if I’ve misread you.[↩](https://learningstatisticswithr.com/book/bayes.html#fnref274)
22. I don’t even disagree with them: it’s *not* at all obvious why a Bayesian ANOVA should reproduce (say) the same set of model comparisons that the Type II testing strategy uses. It’s precisely because of the fact that I haven’t really come to any strong conclusions that I haven’t added anything to the `lsr` package to make Bayesian Type II tests easier to produce.[↩](https://learningstatisticswithr.com/book/bayes.html#fnref275) |
| Readable Markdown | null |
| Shard | 143 (laksa) |
| Root Hash | 3668643058299424143 |
| Unparsed URL | com,learningstatisticswithr!/book/bayes.html s443 |