Evaluating Evidence and Making Decisions using Bayesian Statistics

# Evaluating Evidence and Making Decisions using Bayesian Statistics

## ISCoP Conference 2021

### Mattan S. Ben-Shachar

.right[
.big[<svg aria-hidden="true" role="img" viewBox="0 0 512 512" style="height:1em;width:1em;vertical-align:-0.125em;margin-left:auto;margin-right:auto;font-size:inherit;fill:currentColor;overflow:visible;fill-opacity:white;position:relative;"><path d="M326.612 185.391c59.747 59.809 58.927 155.698.36 214.59-.11.12-.24.25-.36.37l-67.2 67.2c-59.27 59.27-155.699 59.262-214.96 0-59.27-59.26-59.27-155.7 0-214.96l37.106-37.106c9.84-9.84 26.786-3.3 27.294 10.606.648 17.722 3.826 35.527 9.69 52.721 1.986 5.822.567 12.262-3.783 16.612l-13.087 13.087c-28.026 28.026-28.905 73.66-1.155 101.96 28.024 28.579 74.086 28.749 102.325.51l67.2-67.19c28.191-28.191 28.073-73.757 0-101.83-3.701-3.694-7.429-6.564-10.341-8.569a16.037 16.037 0 0 1-6.947-12.606c-.396-10.567 3.348-21.456 11.698-29.806l21.054-21.055c5.521-5.521 14.182-6.199 20.584-1.731a152.482 152.482 0 0 1 20.522 17.197zM467.547 44.449c-59.261-59.262-155.69-59.27-214.96 0l-67.2 67.2c-.12.12-.25.25-.36.37-58.566 58.892-59.387 154.781.36 214.59a152.454 152.454 0 0 0 20.521 17.196c6.402 4.468 15.064 3.789 20.584-1.731l21.054-21.055c8.35-8.35 12.094-19.239 11.698-29.806a16.037 16.037 0 0 0-6.947-12.606c-2.912-2.005-6.64-4.875-10.341-8.569-28.073-28.073-28.191-73.639 0-101.83l67.2-67.19c28.239-28.239 74.3-28.069 102.325.51 27.75 28.3 26.872 73.934-1.155 101.96l-13.087 13.087c-4.35 4.35-5.769 10.79-3.783 16.612 5.864 17.194 9.042 34.999 9.69 52.721.509 13.906 17.454 20.446 27.294 10.606l37.106-37.106c59.271-59.259 59.271-155.699.001-214.959z"/></svg> [.white[tinyurl.com/ISCoP-2021-bayes]](https://tinyurl.com/ISCoP-2021-bayes)]  
.big[<svg aria-hidden="true" role="img" viewBox="0 0 512 512" style="height:1em;width:1em;vertical-align:-0.125em;margin-left:auto;margin-right:auto;font-size:inherit;fill:currentColor;overflow:visible;fill-opacity:white;position:relative;"><path d="M512 0C460.22 3.56 96.44 38.2 71.01 287.61c-3.09 26.66-4.84 53.44-5.99 80.24l178.87-178.69c6.25-6.25 16.4-6.25 22.65 0s6.25 16.38 0 22.63L7.04 471.03c-9.38 9.37-9.38 24.57 0 33.94 9.38 9.37 24.59 9.37 33.98 0l57.13-57.07c42.09-.14 84.15-2.53 125.96-7.36 53.48-5.44 97.02-26.47 132.58-56.54H255.74l146.79-48.88c11.25-14.89 21.37-30.71 30.45-47.12h-81.14l106.54-53.21C500.29 132.86 510.19 26.26 512 0z"/></svg> <svg aria-hidden="true" role="img" viewBox="0 0 640 512" style="height:1em;width:1.25em;vertical-align:-0.125em;margin-left:auto;margin-right:auto;font-size:inherit;fill:currentColor;overflow:visible;fill-opacity:white;position:relative;"><path d="M593.8 59.1H46.2C20.7 59.1 0 79.8 0 105.2v301.5c0 25.5 20.7 46.2 46.2 46.2h547.7c25.5 0 46.2-20.7 46.1-46.1V105.2c0-25.4-20.7-46.1-46.2-46.1zM338.5 360.6H277v-120l-61.5 76.9-61.5-76.9v120H92.3V151.4h61.5l61.5 76.9 61.5-76.9h61.5v209.2zm135.3 3.1L381.5 256H443V151.4h61.5V256H566z"/></svg> [.white[Raw rmarkdown file]](https://github.com/mattansb/bayesian-evidence/blob/main/bayesian-evidence.Rmd)]  
.big[Presented at February 23, 2021]  
.small[(updated June 09, 2021)]
]

---

# About Me

## Mattan S. Ben-Shachar

### PhD Student + Stats Lover + R Developer

---

background-color: var(--myred)
class: inverse

.big[
<svg aria-hidden="true" role="img" viewBox="0 0 512 512" style="height:1em;width:1em;vertical-align:-0.125em;margin-left:auto;margin-right:auto;font-size:inherit;fill:currentColor;overflow:visible;fill-opacity:white;position:relative;"><path d="M512 199.652c0 23.625-20.65 43.826-44.8 43.826h-99.851c16.34 17.048 18.346 49.766-6.299 70.944 14.288 22.829 2.147 53.017-16.45 62.315C353.574 425.878 322.654 448 272 448c-2.746 0-13.276-.203-16-.195-61.971.168-76.894-31.065-123.731-38.315C120.596 407.683 112 397.599 112 385.786V214.261l.002-.001c.011-18.366 10.607-35.889 28.464-43.845 28.886-12.994 95.413-49.038 107.534-77.323 7.797-18.194 21.384-29.084 40-29.092 34.222-.014 57.752 35.098 44.119 66.908-3.583 8.359-8.312 16.67-14.153 24.918H467.2c23.45 0 44.8 20.543 44.8 43.826zM96 200v192c0 13.255-10.745 24-24 24H24c-13.255 0-24-10.745-24-24V200c0-13.255 10.745-24 24-24h48c13.255 0 24 10.745 24 24zM68 368c0-11.046-8.954-20-20-20s-20 8.954-20 20 8.954 20 20 20 20-8.954 20-20z"/></svg> You use statistical models (in `R`)
- ANOVAs, regression
- Maybe some mixed models
]

.big[
<svg aria-hidden="true" role="img" viewBox="0 0 512 512" style="height:1em;width:1em;vertical-align:-0.125em;margin-left:auto;margin-right:auto;font-size:inherit;fill:currentColor;overflow:visible;fill-opacity:white;position:relative;"><path d="M512 199.652c0 23.625-20.65 43.826-44.8 43.826h-99.851c16.34 17.048 18.346 49.766-6.299 70.944 14.288 22.829 2.147 53.017-16.45 62.315C353.574 425.878 322.654 448 272 448c-2.746 0-13.276-.203-16-.195-61.971.168-76.894-31.065-123.731-38.315C120.596 407.683 112 397.599 112 385.786V214.261l.002-.001c.011-18.366 10.607-35.889 28.464-43.845 28.886-12.994 95.413-49.038 107.534-77.323 7.797-18.194 21.384-29.084 40-29.092 34.222-.014 57.752 35.098 44.119 66.908-3.583 8.359-8.312 16.67-14.153 24.918H467.2c23.45 0 44.8 20.543 44.8 43.826zM96 200v192c0 13.255-10.745 24-24 24H24c-13.255 0-24-10.745-24-24V200c0-13.255 10.745-24 24-24h48c13.255 0 24 10.745 24 24zM68 368c0-11.046-8.954-20-20-20s-20 8.954-20 20 8.954 20 20 20 20-8.954 20-20z"/></svg> You've heard about (and maybe even used) Bayes factors
]

---

<br>

<br>

---

# Outline

- What is a Bayesian model?

- How to Bayes, even?

- Why to Bayes? (aka "Why is this better than how I currently model?")

- Demo: Building a Bayesian model
  - Posterior Estimates
  
--

- **Evaluating Evidence and Making Decisions using Bayesian Statistics**

*Let us begin...*

---

# It's all About the <br> <s>Bass</s> Bayesian Modeling

---

# What *is* a Bayesian model?

A Bayesian model is a statistical model where you use **probability** to represent **all uncertainty** within the model, both the uncertainty regarding the output but also the uncertainty regarding the input (aka parameters) to the model<sup>1</sup>...

.footnote[
[1] Bååth (2015). *From [stackexchange](https://stats.stackexchange.com/a/129712/293056)*
]

???

- "uncertainty regarding the output" = how un/certain we are about our predictions.
- "uncertainty regarding the input" = how un/certain we are about our parameters.

... where probability expresses *a degree of belief* in an event.

---

# How to Bayes? <svg aria-hidden="true" role="img" viewBox="0 0 576 512" style="height:1em;width:1.12em;vertical-align:-0.125em;margin-left:auto;margin-right:auto;font-size:inherit;fill:currentColor;overflow:visible;fill-opacity:white;position:relative;"><path d="M208 0c-29.9 0-54.7 20.5-61.8 48.2-.8 0-1.4-.2-2.2-.2-35.3 0-64 28.7-64 64 0 4.8.6 9.5 1.7 14C52.5 138 32 166.6 32 200c0 12.6 3.2 24.3 8.3 34.9C16.3 248.7 0 274.3 0 304c0 33.3 20.4 61.9 49.4 73.9-.9 4.6-1.4 9.3-1.4 14.1 0 39.8 32.2 72 72 72 4.1 0 8.1-.5 12-1.2 9.6 28.5 36.2 49.2 68 49.2 39.8 0 72-32.2 72-72V64c0-35.3-28.7-64-64-64zm368 304c0-29.7-16.3-55.3-40.3-69.1 5.2-10.6 8.3-22.3 8.3-34.9 0-33.4-20.5-62-49.7-74 1-4.5 1.7-9.2 1.7-14 0-35.3-28.7-64-64-64-.8 0-1.5.2-2.2.2C422.7 20.5 397.9 0 368 0c-35.3 0-64 28.6-64 64v376c0 39.8 32.2 72 72 72 31.8 0 58.4-20.7 68-49.2 3.9.7 7.9 1.2 12 1.2 39.8 0 72-32.2 72-72 0-4.8-.5-9.5-1.4-14.1 29-12 49.4-40.6 49.4-73.9z"/></svg>

---

### .blue[A Prior]

A probability distribution representing your prior *belief* about the probability of possible values each parameter can take.

> *"Sounds too subjective to be used in Science!"*  
>     .center[\- You (2021)?]

In real life applications, you would be hard-pressed to just use whatever prior you like - you would need to somehow **justify your prior** (which requires domain specific knowledge).

Watch also [Bürkner (2018). *Why not to be afraid of priors (too much)*](https://www.youtube.com/watch?v=Uz9r8eV2erQ)

]

Similar to how you must also justify and use a reasonable likelihood function.

---

### .orange[A Likelihood Function]

What process best describes the (conditional) data generation process?

For example:
- A .orange[binomial] likelihood function for **binary** data
- A .orange[Poisson] likelihood function for **count** data
- A .orange[cumulative multinomial] likelihood function for **ordinal** data
- An .orange[inverse Gaussian / ex-Gaussian / [other]] likelihood function for **reaction times**
- ...
- A .orange[Gaussian] likelihood function for **conditionally normal** data

The likelihood function tells us the *probability of observing our data given the value(s) of some parameter(s)*.

This function is used to ***update the priors***, resulting in ***The Posterior***...

---

### Prior + Likelihhod = .green[Posterior]

This is that whole pesky *Bayes' Rule* thing everyone keeps going on about:

.content-box-green[
`$$\overbrace{P(\theta|Data)}^{\text{Posterior}} = \frac{\overbrace{P(Data|\theta)}^{\text{Likelihood}} \times \overbrace{P(\theta)}^{\text{Prior}}}{P(Data)}$$`

In words:

The **posterior probability** of some parameter `$\theta$` having a value of `$x$`, is a equal to probability of the observed data occurring if that were the value of `$\theta$` (**the likelihood**), normalized by our **prior belief** that `$\theta$` can have a value of `$x$`.
]

---

![](index_files/figure-html/normal_prior1a-1.png)

---

![](index_files/figure-html/normal_prior1b-1.png)

---

![](index_files/figure-html/normal_prior1c-1.png)

![](index_files/figure-html/normal_prior2-1.png)

![](index_files/figure-html/normal_prior3-1.png)

![](index_files/figure-html/normal_prior4-1.png)

---

???

You will usually see some bell-curve-like prior, but priors can, in theory, take any weird shape you can think of.

![](index_files/figure-html/weird_prior1-1.png)

???

This is also called a horse-shoe prior - used for regularization.

![](index_files/figure-html/weird_prior2-1.png)

![](index_files/figure-html/weird_prior3-1.png)

![](index_files/figure-html/weird_prior4-1.png)

???

The posterior is not only affected by how strong or weak or priors are, but also by how stong/weak or data is...

---

![](index_files/figure-html/lik_plot1-1.png)

???

Stronger data = larger sample sizes and stronger effects - lead to more specific likelihood functions.

![](index_files/figure-html/lik_plot2-1.png)![](index_files/figure-html/lik_plot2-2.png)

---

background-color: var(--myred)
class: inverse
layout: true

# Why to Bayes?

### *AKA* "Why is this better than what I currently do?"

<hr>

---

- **Speak in the language of probabilities** (*probabilitese?*).

> *There is a 0.2 (posterior) probability of the treatment alleviating more than 3 ADHD symptoms.*

> *There is a 0.85 (posterior) probability of realibility of the test being at least `$\alpha > 0.8$`.*

- **The power of Priors**

- Utilize prior knowledge - *add* the information gained from the current data to the existing corpus of knowledge.

- Not every study is *tabula rasa*.

- Use priors to prevent over-fitting (regularization via horseshoe, spike-and-slab).

???

We don't have to invent the wheel... And we can use priors to be cautious (regularization)...

---

**Fit complex models / to complex data**:

- Limiting the search space of our model's parameters to what is *a-priori* reasonable, reduces issues that plague other estimation methods.

- failed convergence, local maxima, complete separation...
  
--

- With a likelihood function and a prior, you can add endless complexity to your model (even allow `$n<p$`).
  
  - Easily model heteroscedasticity,
  - Model individual differences in ICC in HLM,
  - Easily obtain CIs for random effects,
  - ...

- **Some types of models cannot practically be analyzed using frequentists methods** 🤷 [(Rouder & Lu, 2005)](https://twitter.com/Nate__Haines/status/1360227275711668228)

---

# Demo

## Let's get our hands dirty...

.right[
See the full analysis script [here <svg aria-hidden="true" role="img" viewBox="0 0 512 512" style="height:1em;width:1em;vertical-align:-0.125em;margin-left:auto;margin-right:auto;font-size:inherit;fill:currentColor;overflow:visible;fill-opacity:#bf94e4;position:relative;"><path d="M326.612 185.391c59.747 59.809 58.927 155.698.36 214.59-.11.12-.24.25-.36.37l-67.2 67.2c-59.27 59.27-155.699 59.262-214.96 0-59.27-59.26-59.27-155.7 0-214.96l37.106-37.106c9.84-9.84 26.786-3.3 27.294 10.606.648 17.722 3.826 35.527 9.69 52.721 1.986 5.822.567 12.262-3.783 16.612l-13.087 13.087c-28.026 28.026-28.905 73.66-1.155 101.96 28.024 28.579 74.086 28.749 102.325.51l67.2-67.19c28.191-28.191 28.073-73.757 0-101.83-3.701-3.694-7.429-6.564-10.341-8.569a16.037 16.037 0 0 1-6.947-12.606c-.396-10.567 3.348-21.456 11.698-29.806l21.054-21.055c5.521-5.521 14.182-6.199 20.584-1.731a152.482 152.482 0 0 1 20.522 17.197zM467.547 44.449c-59.261-59.262-155.69-59.27-214.96 0l-67.2 67.2c-.12.12-.25.25-.36.37-58.566 58.892-59.387 154.781.36 214.59a152.454 152.454 0 0 0 20.521 17.196c6.402 4.468 15.064 3.789 20.584-1.731l21.054-21.055c8.35-8.35 12.094-19.239 11.698-29.806a16.037 16.037 0 0 0-6.947-12.606c-2.912-2.005-6.64-4.875-10.341-8.569-28.073-28.073-28.191-73.639 0-101.83l67.2-67.19c28.239-28.239 74.3-28.069 102.325.51 27.75 28.3 26.872 73.934-1.155 101.96l-13.087 13.087c-4.35 4.35-5.769 10.79-3.783 16.612 5.864 17.194 9.042 34.999 9.69 52.721.509 13.906 17.454 20.446 27.294 10.606l37.106-37.106c59.271-59.259 59.271-155.699.001-214.959z"/></svg>](https://mattansb.github.io/bayesian-evidence/files/full%20analysis%20script.nb.html)
]

???

Due to time constraints (fitting Bayesian models does take some time), I will walk you through the process of model fitting, exploration, and inference.

---

We will be looking at a regression model,

but the tools from this demo can be applied to Bayesian [SEM](https://faculty.missouri.edu/~merklee/blavaan/), IRT, SDT, [etc](https://cran.r-project.org/view=Bayesian)...

---

## The Data

Thirty 4 year old children completed **Flanker's task**. .small[(real data.)]

<img src="img/flanker_fishC.png" width="80%"/>
]
]

<img src="img/flanker_fishI.png" width="80%"/>
]
]

<img src="img/flanker_fishN.png" width="40%"/>
]

We will be examining their **Interference** (Incongruent - Neutral) and **Facilitation** (Neutral - Congruent) effects, **controlling for age** (in months).

---

We will be working in **`R`** with the following packages:

- `brms` for Bayesian Regression Models with *Stan*.

- *Stan* is a probabilistic programming language

- `emmeans` for extracting estimates / contrasts / slopes from the model.

- `bayestestR` for descriptive and inferential statistics.

- Plots are made with `ggplot2` + `patchwork` + `tidybayes` + `ggdist` + `see`.

.footnote[
[See other package versions and packages used <svg aria-hidden="true" role="img" viewBox="0 0 512 512" style="height:1em;width:1em;vertical-align:-0.125em;margin-left:auto;margin-right:auto;font-size:inherit;fill:currentColor;overflow:visible;fill-opacity:#bf94e4;position:relative;"><path d="M326.612 185.391c59.747 59.809 58.927 155.698.36 214.59-.11.12-.24.25-.36.37l-67.2 67.2c-59.27 59.27-155.699 59.262-214.96 0-59.27-59.26-59.27-155.7 0-214.96l37.106-37.106c9.84-9.84 26.786-3.3 27.294 10.606.648 17.722 3.826 35.527 9.69 52.721 1.986 5.822.567 12.262-3.783 16.612l-13.087 13.087c-28.026 28.026-28.905 73.66-1.155 101.96 28.024 28.579 74.086 28.749 102.325.51l67.2-67.19c28.191-28.191 28.073-73.757 0-101.83-3.701-3.694-7.429-6.564-10.341-8.569a16.037 16.037 0 0 1-6.947-12.606c-.396-10.567 3.348-21.456 11.698-29.806l21.054-21.055c5.521-5.521 14.182-6.199 20.584-1.731a152.482 152.482 0 0 1 20.522 17.197zM467.547 44.449c-59.261-59.262-155.69-59.27-214.96 0l-67.2 67.2c-.12.12-.25.25-.36.37-58.566 58.892-59.387 154.781.36 214.59a152.454 152.454 0 0 0 20.521 17.196c6.402 4.468 15.064 3.789 20.584-1.731l21.054-21.055c8.35-8.35 12.094-19.239 11.698-29.806a16.037 16.037 0 0 0-6.947-12.606c-2.912-2.005-6.64-4.875-10.341-8.569-28.073-28.073-28.191-73.639 0-101.83l67.2-67.19c28.239-28.239 74.3-28.069 102.325.51 27.75 28.3 26.872 73.934-1.155 101.96l-13.087 13.087c-4.35 4.35-5.769 10.79-3.783 16.612 5.864 17.194 9.042 34.999 9.69 52.721.509 13.906 17.454 20.446 27.294 10.606l37.106-37.106c59.271-59.259 59.271-155.699.001-214.959z"/></svg>](https://mattansb.github.io/bayesian-evidence/files/full%20analysis%20script.nb.html#setup)
]

---

## Building a Bayesian Model

---

```r
m_flanker <- brm(
* RT ~ Congruency + age_mo + (Congruency | id),
* data = child_flanker,
  prior = 
    # Two parameters for Congruency 
    set_prior("student_t(3, 0, 100)", class = "b",
              coef = c("Congruency1", "Congruency2")) +
    # Slope of age_mo
    set_prior("student_t(3, 0, 1000)", class = "b",
              coef = "age_mo"),
  family = gaussian())
```

We will be fitting an hierarchical linear model - predicting (single trial) RTs from `Congruency` (I, N, C) which is nested within each child (`id`) - controlling for the children's age (in months, `age_mo`).

This is essentially a repeated measures ANCOVA.

.content-box-red[
Note: For `Congruency` I've used ***orthonormal* dummy-coding**. This is important, but ⏳! Read more about that [here](https://easystats.github.io/bayestestR/articles/bayes_factors.html#contr_bayes).
]

---

```r
m_flanker <- brm(
  RT ~ Congruency + age_mo + (Congruency | id), 
  data = child_flanker, 
* prior =
*   # Two parameters for Congruency
*   set_prior("student_t(3, 0, 100)", class = "b",
*             coef = c("Congruency1", "Congruency2")) +
*   # Slope of age_mo
*   set_prior("student_t(3, 0, 1000)", class = "b",
*             coef = "age_mo"),
  family = gaussian())
```

For our fixed effects, we will be somewhat conservative and use a scaled *t*(3)-prior centered on 0. This prior has the benefit of the scaling factor giving the range where 60% of the prior's mass is.

- In adults, the Flanker effect is about 20-50ms. Here we have 4yo - reasonable (?) that any differences between means (effect) would be **~100ms**, which we will use as our scaling factor ([Jonkman et al, 1999](https://doi.org/10.1111/1469-8986.3640419)).

- Prior on effect of age - no idea. We will use a weakly informative prior scaled to 1000ms/month (covering a very large range of possible effects).

---

.content-box-red[
Notes:
  - By default, `brms` sets **flat** (*diffused, extremely uninformative*) priors for fixed effects.
  - We can also set a prior of `sigma` (error variance), and many others. See more options with `brms::get_prior()`.
]

---

```r
m_flanker <- brm(
  RT ~ Congruency + age_mo + (Congruency | id), 
  data = child_flanker, 
  prior = 
    # Two parameters for Congruency 
    set_prior("student_t(3, 0, 100)", class = "b",
              coef = c("Congruency1", "Congruency2")) + 
    # Slope of age_mo 
    set_prior("student_t(3, 0, 1000)", class = "b",
              coef = "age_mo"), 
* family = gaussian())
```

We will be using a Gaussian likelihood function of `$RT \sim N(\mu_i, \sigma^2)$`, where `$\mu_i =a + \sum b_j X_{ij}$`.

AKA, a boring linear regression.

---

### Prior & Posterior Checks

???

These checks just make sure that our model is very generally reasonable, and that we've done a good job of sampling from the posterior...

.footnote[
But you can find them, and more, in [the full analysis script <svg aria-hidden="true" role="img" viewBox="0 0 512 512" style="height:1em;width:1em;vertical-align:-0.125em;margin-left:auto;margin-right:auto;font-size:inherit;fill:currentColor;overflow:visible;fill-opacity:#bf94e4;position:relative;"><path d="M326.612 185.391c59.747 59.809 58.927 155.698.36 214.59-.11.12-.24.25-.36.37l-67.2 67.2c-59.27 59.27-155.699 59.262-214.96 0-59.27-59.26-59.27-155.7 0-214.96l37.106-37.106c9.84-9.84 26.786-3.3 27.294 10.606.648 17.722 3.826 35.527 9.69 52.721 1.986 5.822.567 12.262-3.783 16.612l-13.087 13.087c-28.026 28.026-28.905 73.66-1.155 101.96 28.024 28.579 74.086 28.749 102.325.51l67.2-67.19c28.191-28.191 28.073-73.757 0-101.83-3.701-3.694-7.429-6.564-10.341-8.569a16.037 16.037 0 0 1-6.947-12.606c-.396-10.567 3.348-21.456 11.698-29.806l21.054-21.055c5.521-5.521 14.182-6.199 20.584-1.731a152.482 152.482 0 0 1 20.522 17.197zM467.547 44.449c-59.261-59.262-155.69-59.27-214.96 0l-67.2 67.2c-.12.12-.25.25-.36.37-58.566 58.892-59.387 154.781.36 214.59a152.454 152.454 0 0 0 20.521 17.196c6.402 4.468 15.064 3.789 20.584-1.731l21.054-21.055c8.35-8.35 12.094-19.239 11.698-29.806a16.037 16.037 0 0 0-6.947-12.606c-2.912-2.005-6.64-4.875-10.341-8.569-28.073-28.073-28.191-73.639 0-101.83l67.2-67.19c28.239-28.239 74.3-28.069 102.325.51 27.75 28.3 26.872 73.934-1.155 101.96l-13.087 13.087c-4.35 4.35-5.769 10.79-3.783 16.612 5.864 17.194 9.042 34.999 9.69 52.721.509 13.906 17.454 20.446 27.294 10.606l37.106-37.106c59.271-59.259 59.271-155.699.001-214.959z"/></svg>](https://mattansb.github.io/bayesian-evidence/files/full%20analysis%20script.nb.html)...
]

---

## Explore the Model

---

```r
means_Congruency <-
  emmeans(m_flanker, ~ Congruency)
```

<hr>

]

---

```r
means_Congruency <-
  emmeans(m_flanker, ~ Congruency)
```

<hr>

Frequentist estimation methods (such as **OLS** or **maximum likelihood (ML)**) produce a point estimate for each parameter.

But in Bayes we have not a single value, but .green[a whole distribution of values]!

We we can either .green[present the whole distribution, *as is*]...

]

![](index_files/figure-html/plot_cong_means-1.png)

]

---

Or we can summarize the posterior distribution:

- Median (most common)
- Mean
- Maximum A Posteriori (MAP)

]

- The Highest Density Interval (HDI; most common)
- The Equal-Tailed Interval (ETI)

]

<hr>

```r
describe_posterior(means_Congruency, 
                   centrality = "median",
                   ci = 0.89, ci_method = "hdi",
                   test = NULL)
```

```
## Summary of Posterior Distribution
## 
## Parameter   |  Median |             89% CI
## ------------------------------------------
## Incongruent | 1602.49 | [1420.11, 1753.12]
## Neutral     | 1433.93 | [1286.05, 1587.07]
## Congruent   | 1490.33 | [1341.28, 1639.82]
```

---

background-color: var(--myred)
class: inverse, center, middle, title-slide

## Evaluating Evidence and Making Decisions <br> using Bayesian Statistics <br> <svg aria-hidden="true" role="img" viewBox="0 0 512 512" style="height:1em;width:1em;vertical-align:-0.125em;margin-left:auto;margin-right:auto;font-size:inherit;fill:currentColor;overflow:visible;fill-opacity:white;position:relative;"><path d="M505 442.7L405.3 343c-4.5-4.5-10.6-7-17-7H372c27.6-35.3 44-79.7 44-128C416 93.1 322.9 0 208 0S0 93.1 0 208s93.1 208 208 208c48.3 0 92.7-16.4 128-44v16.3c0 6.4 2.5 12.5 7 17l99.7 99.7c9.4 9.4 24.6 9.4 33.9 0l28.3-28.3c9.4-9.4 9.4-24.6.1-34zM208 336c-70.7 0-128-57.2-128-128 0-70.7 57.2-128 128-128 70.7 0 128 57.2 128 128 0 70.7-57.2 128-128 128z"/></svg>

---

We are limiting our discussion to evaluating evidence for **single estimates / parameters** (expected values, slopes, contrasts...).

But it is also possible to evaluating evidence for multiple parameters, with order restrictions and model comparisons. (Maybe next year...)

<hr>

We will be looking at two contrasts: the Interference and Facilitation effects:

```r
diffs_Congruency <- contrast(means_Congruency, 
                             list(Interference = c(1, -1, 0),
                                  Facilitation = c(0,  1, -1)))

describe_posterior(diffs_Congruency, test = NULL, ci = .89)
```

```
## Summary of Posterior Distribution
## 
## Parameter    | Median |            89% CI
## -----------------------------------------
## Interference | 166.48 | [  48.51, 285.05]
## Facilitation | -56.01 | [-154.24,  35.10]
```

---

### The Probability of Direction

- The maximal probability of the estimate being strictly directional (larger or smaller than 0).

- Generally ranges from 50% (no preference) to 100%.

```r
p_direction(diffs_Congruency)
```

```
## Probability of Direction
## 
## Parameter    |     pd
## ---------------------
## Interference | 98.58%
## Facilitation | 83.58%
```

]

![](index_files/figure-html/pd_plot-1.png)

]

For the Interference effect it seems like these is a high probability of direction, but not that great for the Facilitation effect ( `$p_d$` < 0.95 ).

- <b>.green[Pros]</b>: Easy to understand; Resembles the *p*-value - `$r \simeq -1$`. <sup>*</sup>

- <b>.red[Cons]</b>: like the *p*-values, a *low* `$p_d$` cannot be used to support the null.

---

### *p*-MAP

- The *density ratio* between the null and the MAP value.

- Values range from 1 (the null *is* the MAP) to ~0 (the MAP is much much more probable than the null).

```r
p_map(diffs_Congruency)
```

```
## MAP-based p-value
## 
## Parameter    | p (MAP)
## ----------------------
## Interference |   0.086
## Facilitation |   0.651
```

]

![](index_files/figure-html/pmap_plot-1.png)

]

For the Interference effect it seems like the MAP is more th 10 times more probable than the null. But for the Facilitation effect it is not even twice as probable.

- <b>.green[Pros]</b>: Closely related to LRT tests - familiar; Also closely associated with the *p*-value.

- <b>.red[Cons]</b>: Again, a *high* *p*-MAP cannot be used to support the null.

---

### *p*-ROPE

- The probability that our estimate is *basically* null.

- We first define a **Region of Practical Equivalence (ROPE)** - a range of effects that are, for any practical purposes, the same as no effect at all.

For the Congruency effects, we will define any effect that is smaller in magnitude than 30ms, to be consider to be just as good as no effect at all - so ROPE [-30, +30].

---

### *p*-ROPE

- How much of the posterior falls in the ROPE.

- Or: How much of the most probable values (e.g., those in the HDI) fall in the ROPE.

```r
rope(diffs_Congruency, 
     range = c(-30, 30), ci = 0.89)
```

```
## # Proportion of samples inside the ROPE [-30.00, 30.00]:
## 
## Parameter    | inside ROPE
## --------------------------
## Interference |      0.00 %
## Facilitation |     29.99 %
```

]

![](index_files/figure-html/therope_plot-1.png)

]

It is very improbable that the Interference effect is very small. *But* there is about a 30% that among 4 year olds, there is no Facilitation effect - ([though not very conclusive](https://easystats.github.io/bayestestR/articles/guidelines.html#significance)) we are supporting the null!

???

If you've ever heard that "Bayes is good for small samples" this is what is meant by that: that unlike frequentist methods where small samples and non-significant results leave you high and dry, Bayes allows you to same *something*, weak as it may be.

---

The *p*<sub>*d*</sub>, *p*-MAP and ROPE are **posterior based methods** - they inform us about the accumulated information in the priors + our data.

Often we are interested in **what has been *learned* in the current study, from the current data**.

E.g., by now it's clear that there exists an Interference effect in Flanker's task. But **which values of the effect are supported or contradicted by the *current* data?** Maybe our data supports the null value(s) - what can be learnt from that?

To answer these types of questions we can *compare* .blue[The Prior] to .green[The Posterior] to see .orange[what our data taught us] - what values became more / less plausible.

]

![](index_files/figure-html/support_plot-1.png)

]

We can use this information to look at different *sets of parameter values* - or **hypotheses** - E.g. `$H_{small}: \theta \in [-3, 3]$`, `$H_{positive}: \theta \in [0, \infty]$` and ask:

> Which *hypothesis* is supported *more* by the data?

---

### The Bayes Factor

This index of evidence is a ***Bayes Factor***:

- It quantifies how the prior was *updated* to the posterior.

- It compares two "hypotheses".

**Any** measure that quantifies this 👆 is a Bayes factor.

There are many different type of questions that can be answered with Bayes factors - we will be looking at two.

]

???

You maybe have used Bayes factors to compare between models... Those too have these properties.

<hr>

For technical reasons we need a model that represents *only our priors* - which we will then *compare* to the results from our updated (posterior) model.

We can do that with the `unupdate()` function:

```r
# Get the priors only ("un-update" the model).
m_flanker_prior <- unupdate(m_flanker)
```

---

### The Null-Interval Bayes Factor

---

The null-interval Bayes factor is an extension of the ROPE test;

> How has the *relative probability*<sup>[1]</sup> of the the effect being practically null changed? Does the data support or contradict the effect being null?

<hr>

The two hypotheses we will be comparing, using the same ROPE:

- `$H_0: \text{effect} \in [-30, +30]$`  
- `$H_A: \text{effect} \notin [-30, +30]$`  
  - Or: `$H_A: \text{effect} < -30$` or `$+30 < \text{effect}$`

---

```r
bayesfactor_parameters(
  diffs_Congruency, 
  prior = m_flanker_prior,
  null = c(-30, 30) # same ROPE as before
)
```

```
## Bayes Factor (Null-Interval)
## 
## Parameter    |    BF
## --------------------
## Interference |  6.02
## Facilitation | 0.518
## 
## * Evidence Against The Null: [-30.000, 30.000]
```

]

![](index_files/figure-html/bf_ROPE_plot-1.png)

]

- For the Interference effect, the ROPE has *become* relatively less probable - with the data giving 6 times more support for non-ROPE values.

- For the Facilitation effect, the ROPE has *become* relatively **more** probable - with the data giving (1/0.5 =) 2 times more support compared to the non-ROPE values.

---

### The Point-Null Bayes Factor

---

The point-null can be thought of as the null-interval Bayes factor with an infinitesimally small ROPE - that includes only one null value, exactly.

> How has the probability<sup>[1,2]</sup> of the the null value changed? Does the data support or contradict the effect being null?

This Bayes factor is also called the *Savage-Dickey density ratio*, .small[and it is analogous to a Bayes factor comparing two nested models.]

[1] Actually the density of the null.  
[2] This is also relative - if the null became more probable, necessarily the non-null values became less, and vice versa.

]

<hr>

The two hypotheses we will be comparing:

- `$H_0: \text{effect} = 0$`  
- `$H_A: \text{effect} \neq 0$`  
  - Or: `$H_A: \text{effect} < 0$` or `$0 < \text{effect}$`

---

```r
bayesfactor_parameters(
  diffs_Congruency, 
  prior = m_flanker_prior,
  null = 0
)
```

```
## Bayes Factor (Savage-Dickey density ratio)
## 
## Parameter    |    BF
## --------------------
## Interference |  5.93
## Facilitation | 0.591
## 
## * Evidence Against The Null: 0
```

]

![](index_files/figure-html/bf_point_plot-1.png)

]

- For the Interference effect, the *mass* of the posterior is shifted *away* from the null (compared to the prior) - the data giving ~6 times more support for non-null values.

- For the Facilitation effect the mass has moved *towards* 0, the data giving (1/0.6 =) 1.7 times more support compared to the non-null values.

---

```r
bayesfactor_parameters(
  diffs_Congruency, 
  prior = m_flanker_prior,
  null = 0
)
```

```
## Bayes Factor (Savage-Dickey density ratio)
## 
## Parameter    |    BF
## --------------------
## Interference |  5.93
## Facilitation | 0.591
## 
## * Evidence Against The Null: 0
```

]

```r
bayesfactor_parameters(
  diffs_Congruency, 
  prior = m_flanker_prior,
  null = c(-30, 30) # same ROPE as before
)
```

```
## Bayes Factor (Null-Interval)
## 
## Parameter    |    BF
## --------------------
## Interference |  6.02
## Facilitation | 0.518
## 
## * Evidence Against The Null: [-30.000, 30.000]
```

]

Here the point-null and the null-interval BFs gave similar results, but that need not be the case - depending on the effect size, the definition of the ROPE, the sample size, etc.

]

---

### Other Bayes Factors

- **Directional** null-interval / point-null Bayes factors

- e.g., [-30, +30] *vs* [+30, Inf]

- Bayes factor for **dividing hypotheses**

- e.g., [-Inf, 0] *vs* [0, Inf]

- **Model restricted** Bayes factors

- [Incongruent > Neutral > Congruent] vs [Incongruent ≠ Neutral ≠ Congruent]
  
- And more...

Read more about these Bayes factors [here](https://easystats.github.io/bayestestR/articles/bayes_factors.html)!

---

## Age

---

For *covariates*, we can present the posterior distribution of slopes, but be can also present a *trace plot* of slopes from the posterior.

]

![](index_files/figure-html/age_lines_plot-1.png)

]

---

Here too we can summarize the posterior distribution:

```r
slope_age <- emtrends(m_flanker, ~1, "age_mo")

describe_posterior(slope_age, test = NULL, ci = .89)
```

```
## Summary of Posterior Distribution
## 
## Parameter | Median |           89% CI
## -------------------------------------
## overall   |  22.48 | [-47.48, 101.12]
```

]

![](index_files/figure-html/slopes_standard_plot-1.png)

]

---

#### *p*-Direction & *p*-MAP

```r
p_direction(slope_age)
```

```
## Probability of Direction
## 
## Parameter |     pd
## ------------------
## overall   | 67.70%
```

]

```r
p_map(slope_age)
```

```
## MAP-based p-value
## 
## Parameter | p (MAP)
## -------------------
## overall   |   0.830
```

]

Not very decisive… (remember, these cannot be used to support the null!)

---

#### *p*-ROPE

For the ROPE - I think any effect smaller an overall change of less than 500ms a year = ~40ms a month, is practically 0 (you may disagree…):

```r
rope(slope_age, range = c(-40, 40), ci = 0.89)
```

```
## # Proportion of samples inside the ROPE [-40.00, 40.00]:
## 
## Parameter | inside ROPE
## -----------------------
## overall   |     63.89 %
```

There is about a 60% probability that the effect of age on reaction times is *practically* nothing!

Not strongly conclusive, but at the very least it is suggestive!

---

#### Bayes Factor

---

```r
bayesfactor_parameters(
  slope_age,
  prior = m_flanker_prior,
  null = 0
)
```

```
## Bayes Factor (Savage-Dickey density ratio)
## 
## Parameter |    BF
## -----------------
## overall   | 0.052
## 
## * Evidence Against The Null: 0
```

Wow! It seems that the data strongly support (by a factor of 1/0.05 = 20) the effect of age being null over it being non-null!

*But wait* - the Bayes factor measures the change from the prior to the posterior… But what was our prior here?

---

```r
bayesfactor_parameters(
  slope_age,
  prior = m_flanker_prior,
  null = 0
)
```

```
## Bayes Factor (Savage-Dickey density ratio)
## 
## Parameter |    BF
## -----------------
## overall   | 0.052
## 
## * Evidence Against The Null: 0
```

]

![](index_files/figure-html/slopes_BF_plot-1.png)

]

We used a super vague prior - which give some non-trivial probability to extreme effects!

So is it really surprising that the posterior is now, relatively closer to the *null*? ***No.***

???

In our prior, the null was very very importable - it is not therefore surprising that it became *more* probable.

.content-box-red[
With wide and uninformative enough priors, the Bayes factor will **always favor the null / ROPE**!  
DO NOT COMPUTE BAYES FACTORS WITH UNINFORMATIVE PRIORS! <sup>*</sup>
]

???

This is only an issue if one of your hypotheses is a point, or almost a point.

---

background-color: var(--myred)
class: inverse
layout: false

# Recommendations

### *What to actually report?*

We (Makowski et al., 2019) [recommend](https://easystats.github.io/bayestestR/articles/guidelines.html) reporting for inferential statistics:

- **The *p*-direction**: Easy to understand, easy to "translate" to *p*-values.

- ***p*-ROPE**: Provides information about the practical relevance of the effect, and allows to accept the null.

*If* informed priors are used,

- **Bayes factor** .small[(instead or in addition to the *p*-ROPE)]: Provides information about hypotheses supported or contradicted by the data.

---

# Summary

### *What you now know! <svg aria-hidden="true" role="img" viewBox="0 0 640 512" style="height:1em;width:1.25em;vertical-align:-0.125em;margin-left:auto;margin-right:auto;font-size:inherit;fill:currentColor;overflow:visible;fill-opacity:white;position:relative;"><path d="M622.34 153.2L343.4 67.5c-15.2-4.67-31.6-4.67-46.79 0L17.66 153.2c-23.54 7.23-23.54 38.36 0 45.59l48.63 14.94c-10.67 13.19-17.23 29.28-17.88 46.9C38.78 266.15 32 276.11 32 288c0 10.78 5.68 19.85 13.86 25.65L20.33 428.53C18.11 438.52 25.71 448 35.94 448h56.11c10.24 0 17.84-9.48 15.62-19.47L82.14 313.65C90.32 307.85 96 298.78 96 288c0-11.57-6.47-21.25-15.66-26.87.76-15.02 8.44-28.3 20.69-36.72L296.6 284.5c9.06 2.78 26.44 6.25 46.79 0l278.95-85.7c23.55-7.24 23.55-38.36 0-45.6zM352.79 315.09c-28.53 8.76-52.84 3.92-65.59 0l-145.02-44.55L128 384c0 35.35 85.96 64 192 64s192-28.65 192-64l-14.18-113.47-145.03 44.56z"/></svg>*

- What a Bayesian model *is*.

- What Bayes can give you, that no one else can.

- A taste of Bayesian model fitting with `brms`.

- The richness of inferences that can be made with Bayesian statistics.

---

# Suggested Reading

#### For Bayesian Beginners

- Makowski, D., Ben-Shachar, M. S., Chen, S. H., & Lüdecke, D. (2019). [Indices of effect existence and significance in the Bayesian framework. *Frontiers in psychology, 10*, 2767.](https://doi.org/10.3389/fpsyg.2019.02767)

- [`bayestestR` guides and articles.](https://easystats.github.io/bayestestR)

- Van de Schoot, R. et al (2021). [Bayesian statistics and modelling. *Nature Reviews Methods Primers, 1*(1), 1-26.](https://doi.org/10.1038/s43586-020-00001-2)

- [Bayesian Inference for Psychology. *Psychonomic Bulletin and Review*.](https://scholar.google.co.il/scholar?q=Bayesian+inference+for+psychology+Psychonomic+Bulletin+and+Review)

#### Books

- Kruschke, J. (2014). Doing bayesian data analysis: A tutorial with r, jags, and stan. Academic Press.

- McElreath, R. (2018). Statistical rethinking: A bayesian course with examples in r and stan. Chapman; Hall/CRC.

- [Richard's YouTube channel](https://www.youtube.com/channel/UCNJK6_DZvcMqNSzQdEkzvzA)

---

background-color: var(--myred)
class: inverse

.pull-left[
<img style="border-radius: 50%;" src="https://mattansb.github.io/CV/headshots/BrainOrange.jpg" width="150px"/><img src="img/BGU-logo-round-clear.png" width="20%" /><img src="img/lab_logo.png" width="20%" />

# Thank you!

### Follow me!

]

.pull-right[
<br><br><br><br>
.center[
<img src="img/easystats.png" width="30%" /><img src="img/bayestestR.png" width="30%" />
]

.small[
The [**`bayestestR`**](https://easystats.github.io/bayestestR) package is part of the `easystats` project. Core team members:

- Me 👋

- Dominique Makowski ([@Dom_Makowski](https://twitter.com/Dom_Makowski))

- Daniel Lüdecke ([@strengejacke](https://twitter.com/strengejacke))

- Indrajeet Patil ([@patilindrajeets](https://twitter.com/patilindrajeets))
]
]

.footnote[
*Slides created with the R package [**`xaringan`**](https://github.com/yihui/xaringan).*
]

---

[.black[@kareem_carr]](https://twitter.com/kareem_carr/status/1356986263975395333)