What this class is about

Revisiting hypothesis testing

(Frequentist) hypothesis testing framework

  • Suppose \(X\) is the number of heads out of \(n\) independent tosses.
  • Let \(p\) be the probability of getting a for this coin.

Hypotheses

\(H_0: p = 0.5\) vs. \(H_a: p > 0.5\). Note \(p_0=0.5\).
Alternative \(H_a\) is saying we believe that the coin is biased to heads.

NOTE: Alternative needs to be decided before seeing data.

Assumptions Each toss is independent with equal chance of getting a head.

Test statistic

It has been mathematically shown that \(X \sim B(n, p)\). If \(H_0\) is true, we expect the number of heads to be \(np_0\).
We observe \(n, x\) and compute the estimate, \(\widehat{p} = \frac{x}{n}\). If \(np_0 \geq 5, n(1-p_0)\geq 5\) it has been shown mathematically that this test statistic is \(\frac{(\widehat{p} - p_0)}{\sqrt{p_0(1-p_0)/n}}\) has a normal distribution.

This can be used to compute the probability of seeing more than \(x\) heads in \(n\) coin flips, the \(p\)-value. If this is small (e.g. \(<\) 0.05), it is unlikely to have happened, if \(H_0\) was true.

Computationally testing coin bias (1/3)

  • Experiment 1: I flipped the coin 10 times and this is the result:
Code 
set.seed(924)
samp10 <- sample(rep(c(head, tail), c(7, 3)))
cat(paste0(samp10, collapse = ""))

  • The result is 7 head and 3 tails. So 70% are heads.
  • Do you believe the coin is biased based on this data?

Computationally testing coin bias (2/3)

  • Experiment 2: Suppose now I flip the coin 100 times and this is the outcome:
Code 
samp100 <- sample(rep(c(head, tail), c(70, 30)))
cat(paste0(samp100, collapse = ""))

  • We observe 70 heads and 30 tails. So again 70% are heads.
  • Based on this data, do you think the coin is biased?

Computationally testing coin bias (3/3)

Calculate it

Experiment 1 (n=10)

  • We observed \(x=7\), or \(\widehat{p} = 0.7\).

  • Assuming \(H_0\) is true, we expect \(np=10\times 0.5=5\).

  • Calculate the \(P(X \geq 7)\)



Code 
sum(dbinom(7:10, 10, 0.5))
[1] 0.17

Experiment 2 (n=100)

  • We observed \(x=70\), or \(\widehat{p} = 0.7\).

  • Assuming \(H_0\) is true, we expect \(np=100\times 0.5=50\).

  • Calculate the \(P(X \geq 70)\)



Code 
sum(dbinom(70:100, 100, 0.5))
[1] 3.9e-05

Your turn: simulate coin flipping 1

Your turn

  1. Simulate flipping a fair coin 10 times using this code:
coinflips <- sample(c("H", "T"), 10, replace=TRUE)
sum(coinflips == "H")
  1. Enter your number of heads into the ED.ie response.

Your turn: simulate coin flipping 2

Your turn

  1. Simulate flipping a fair coin 100 times using this code:
coinflips <- sample(c("H", "T"), 100, replace=TRUE)
sum(coinflips == "H")
  1. Enter your number into the ED.ie response

Judicial system

Why is the null hypothesis always specific?

You need to be able to calculate the probability of something happening, if the null was true.

Evidence by test statistic computer using the observations.
Judgement by \(p\)-value (probability that the event, or more extreme events, happened assuming the null is true).

Does the test statistic have to be numerical?

Testing hypotheses using data plots

Motivation

Is this an example of a good residual plot, or does it suggest that there is a problem with the model fit?

Code 
cars_fit <- lm(dist ~ speed, data = cars)
cars_lm <- augment(cars_fit)
set.seed(1051)
ggplot(cars_lm, aes(x=speed, y=.resid)) +
  geom_point()

Visual inference concepts

  • Hypothesis testing in a visual inference framework is where:
    • the test statistic is a plot and
    • judgement is by human visual perception.

Why is the plot a test statistic? We’ll see why soon.

  • You, we, me actually do visual inference many times but generally in an informal fashion.
  • The problem with doing this is we are making an inference on whether the plot has any patterns based on a single data plot.
  • The single data plot needs to be examined in the context of what might this look like if different samples were shown.

Reading data plots requires calibration.

Visual inference steps

  1. Write your plot description using ggplot. DON’T MAKE THE DATA PLOT!
  2. Based on the description articulate the null and alternate hypotheses, e.g. \(H_A:\) the two variables have association \(\longrightarrow\) \(H_0:\) there is no association between the two variables.
  3. Decide on a method to generate null data, sample that are consistent with the null hypothesis.
  4. Make \(m-1\) null samples.
  5. Make a lineup of the plots of the null samples with the plot of the data plots randomly inserted. All plots use the same ggplot description, just different samples.
  6. Ask \(n\) human viewers to select a plot in the lineup that looks different to others. Don’t give any context, and disguise axes.

The nullabor::lineup() function can do simple lineups.

Using the same type of calculation as done for the coi flipping, compute the probability that \(x\) out of \(n\) people detected the data plot from a lineup, then

  • the visual inference p-value is given as \[P(X \geq x)\] where \(X \sim B(n, 1/m)\), and
  • the power of a lineup is estimated as \(x/n\).

Easier, still use the nullabor::pvisual() function to compute it.

Example: residual plots (1/4)

Which plot has a pattern that is most different from other plots?

Residuals from dist~speed using datasets::cars.

Code 
lm(dist ~ speed, data = cars)
  • This is a lineup of the residual plot
  • Which plot (if any) looks different from the others?
  • Why do you think it looks different?
> decrypt("clZx bKhK oL 3OHohoOL 0B")
[1] "True data in position  11"


Code 
nullabor::pvisual(2, 16, 20)
     x simulated binom
[1,] 2       0.2  0.19

Residual plots need context


You are asked to decide IF THERE IS NO PATTERN. This is hard!

  • Numerical tests (e.g. Breusch-Pagan, Shapiro-Wilk) are either insensitive or overly-sensitive
  • Visual inspection is effective but subjective and unscalable
  • Different analysts reading the same plot reach different conclusions

Residual plots are better when viewed in the context of good residual plots, where we know the assumptions of the model are satisfied.

Example: residual plots (2/4)

Code 
library(broom)
diamonds <- diamonds %>%
  mutate(lprice = log10(price),
         lcarat = log10(carat))
d_fit <- lm(lprice ~ lcarat, data=diamonds)
d_res <- augment(d_fit, diamonds)

set.seed(923)
l <- lineup(null_lm(lprice ~ lcarat,
                      method="rotate"), d_res)
ggplot(l, aes(lcarat, .resid)) + 
  geom_hline(yintercept=0, colour="grey70") +
  geom_point(alpha = 0.01) +
  geom_smooth(data=l, method = "lm", colour="orange", se=F) +
  facet_wrap(~.sample, ncol=5) +
  theme_bw() +
  theme(axis.text=element_blank(),
        axis.title=element_blank())

Residuals from log-transformed price and carat ggplot2::diamonds. Did it linearise the relationship nicely?

Code 
d_fit <- lm(lprice ~ lcarat, data=diamonds)
  • Which plot (if any) looks different from the others?
  • Why do you think it looks different?
> decrypt("clZx bKhK oL 3OHohoOL 0Q")
[1] "True data in position  15"
Code 
nullabor::pvisual(8, 12, 20)
     x simulated   binom
[1,] 8         0 1.6e-08

Example: Residual plots (3/4)



Is there a problem with the model?



Let’s see what the computer vision model suggests: Shiny web app

Residual plot (4/4)

Use the app to check how the computer vision model performs with the lineup of residual plots from the cars model.

Use the app to check how the computer vision model performs on the second example data set provided in the app.

It’s not only for residual plots

Sports analytics: basketball



Which plot is most different?

Time series: cross-currency rates



Which plot is most different?

Association: cars



Which plot is most different?

Spatial analysis: cancer incidence



Which plot is most different?





From Steff Kobakian’s Master’s thesis

Reading any plot is easier in the context of null plots

Why is a data plot a statistic?

Why is a data plot a statistic? (1/2)

  • The concept of tidy data matches elementary statistics
  • Tabular form puts variables in columns and observations in rows

\[X = \left[ \begin{array}{rrrr} X_1 & X_2 & ... & X_p \end{array} \right] \\ = \left[ \begin{array}{rrrr} X_{11} & X_{12} & ... & X_{1p} \\ X_{21} & X_{22} & ... & X_{2p} \\ \vdots & \vdots & \ddots& \vdots \\ X_{n1} & X_{n2} & ... & X_{np} \end{array} \right]\]

  • Variables can have distributions, e.g. \(X_1 \sim N(0,1), ~~X_2 \sim \text{Exp}(1) ...\)

Why is a data plot a statistic? (2/2)

  • A statistic is a function on the values of items in a sample, e.g. for \(n\) iid random variates \(\bar{X}_1=\sum_{i=1}^n X_{i1}\), \(s_1^2=\frac{1}{n-1}\sum_{i=1}^n(X_{i1}-\bar{X}_1)^2\)
  • We study the behaviour of the statistic over all possible samples of size \(n\).
  • The grammar of graphics is the mapping of (random) variables to graphical elements, making plots of data into statistics

Example 1:

ggplot(threept_sub, 
       aes(x=angle, y=r)) + 
  geom_point(alpha=0.3)


angle is mapped to the x axis

r is mapped to the y axis

Example 2:

ggplot(penguins, 
      aes(x=bl, 
          y=fl, 
          colour=species)) +
  geom_point()


bl is mapped to the x axis

fl is mapped to the y axis

species is mapped to colour

Example 3:

ggplot(aud, aes(x=date, y=rate)) + 
  geom_line()


date is mapped to the x axis

rate is mapped to the y axis

displayed as a line geom

Determining the null hypothesis

What is the null hypothesis? (1/2)

To determine the null hypothesis, you need to think about what pattern would NOT be interesting.

A

ggplot(data) + 
  geom_point(aes(x=x1, y=x2))


B

ggplot(data) + 
  geom_point(aes(x=x1, 
    y=x2, colour=cl))

C

ggplot(data) + 
  geom_histogram(aes(x=x1))


D

ggplot(data) + 
  geom_boxplot(aes(x=cl, y=x1))



🤔 Which of these plot definitions would most match to a null hypothesis stating there is no difference in the distribution between the groups?

What is the null hypothesis? (2/2)



A

\(H_o:\) no association between x1 and x2



B

\(H_o:\) no difference in association of between x1 and x2 between levels of cl

C

\(H_o:\) the distribution of x1 is XXX



D

\(H_o:\) no difference in the distribution of x1 between levels of cl

How do you generate null samples

Primary null-generating mechanisms


Null samples can be generated using two basic approaches:

  • Permutation: randomizing the order of one of the variables breaks association, but keeps marginal distributions the same.
  • Simulation: from a given distribution, or model. Assumption is that the data comes from that model.

applied to subsets, or conditioning on other variables. Simulation may require computing summary statistics from the data to use as parameter estimates.

Association: cars



Null plots generated by permuting x variable.

Time series: cross-currency rates



Nulls generated by simulating from an ARIMA model.

Beyond \(p\)-value to power

What is power?

  • A statistic is said to be more powerful than another statistic if it has a higher probability of correctly rejecting the null hypothesis when the alternative hypothesis is true.
  • The effectiveness of two plots designs for the same data can be compared by computing power from a lineup.
  • The power of a lineup is calculated as \(x/n\) where \(x\) is the number of people who detected the data plot out of \(n\) people.

Power of plot design

We’re going to play a game to determine, which of these plots is more effective for assessing difference between groups.



Each group will see only one of the next four plots for 3 second. Choose which of the plots, 1, …, 10, is most different, and record it until we are ready to collect the data.

Make a note whether it is plot 1,, 2,, 3,, 4,, 5,, 6,, 7,, 8,, 9,, 10, which is most different.

Compiling the results

Share your group, and which plot you chose.

Computing the power

Note: Different people evaluated each lineup.

Plot type \(x\) \(n\) Power
geom_point \(x_1=??\) \(n_1=??\) \(x_1 / n_1=??\)
geom_boxplot \(x_2=??\) \(n_2=??\) \(x_2 / n_2=??\)
geom_violin \(x_3=??\) \(n_3=??\) \(x_3 / n_3=??\)
ggbeeswarm::geom_quasirandom \(x_4=??\) \(n_4=??\) \(x_4 / n_4=??\)


  • The plot type with a higher power is preferable
  • You can use this framework to find the optimal plot design
  • Which is the winner?

The process for choosing a plot design

  1. Decide on your best plot designs
  2. Make lineups, using the same data, same position with each
  3. Show independent groups the lineups to evaluate
  4. Calculate the power for each, or use speed as an alternative
  5. Whichever design produces the higher power, is the one where the reader can best see the structure.

Main takeaways

  • Hypothesis tests can be built computationally, by simulating or permuting data under \(H_0\), instead of relying on textbook reference distributions.
  • A data plot can itself be a test statistic: the lineup protocol embeds the data plot among null plots generated under \(H_0\), and a viewer “rejects” \(H_0\) if they pick out the data plot.
  • Choosing how to generate null plots matters: use permutation when the null implies no association, and simulation (e.g. from a fitted model) when the null specifies a particular data-generating process.
  • Any plot — residual plots, sports analytics, time series, associations, spatial maps — is easier to read and trust when compared against null plots as a reference.
  • The effectiveness of competing plot designs can be measured empirically and quantitatively as power (\(x/n\), the proportion of viewers who spot the real data), giving an evidence-based way to choose between visualization designs.

Resources