PP5000 problem set 1: univariate-stats

Choose one of the following 3 options

Solve exercise 2611121 and exercise 10812121614
Solve exercise 16610101020
Solve exercise 1200110030

Items marked OPTIONAL are extensions, you can skip them and still receive full credit.

The deadline is Tuesday, September 01, 19:30.

You are allowed to discuss exercises and to collaborate, but this is an individual assignment.

Upload a pdf document with your answers. Add brief comments, particularly if you reviewed your answers and found them suspicious or implausible. Don’t overthink formatting and cosmetics: copy-pasting console output or taking screenshots is fine for this assignment.

Acknowledge any assistance received, whether from humans or from AIs.

Alongside your pdf answers, but separate from them, you must provide replication code sufficient to produce every result, table, plot… in your pdf answers.

Exercise 2611121: comparing `2d6+1` and `1d12+1`

Start your code with:

import random
import pandas as pd
random.seed(2611121)
n_realizations=4000

Generate n_realizations realizations of:

2d6+1, i.e. the result of rolling two d6, adding the faces, and adding 1

and n_realizations realizations of:

1d12+1, i.e. the result of rolling 1 d12, and adding 1

Store the results in lists named respectively: realizations_2d6_plus_1 and realizations_1d12_plus_1

Put both lists in the same pandas dataframe.

Produce a table with 2 columns named 2d6_plus_1 and 1d12_plus_1 and rows displaying:

minimum
maximum
mean
variance

as well as 5 rows for the frequency of:

realization >= 8
realization >= 9
realization >= 10
realization >= 11
realization >= 12

10812121614 Hypothesis testing — Reverse the hypotheses

The goal of this exercise is to follow the lecture notes’ example of hypothesis testing, except we flip the hypotheses. The null hypothesis is now H0: the d6 is rigged, and the alternative hypothesis is now H1: the d6 is fair. The test statistic is the sample mean.

Start your code with:

import random
random.seed(10812121614)

n_simulations = 2000
alpha = 0.10
sample_size = 50

interval_widths_fair = [12, 12, 12, 12, 12, 12]
interval_widths_rigged = [10, 8, 12, 12, 16, 14]

def jitter(x, scale):
    return x + random.uniform(-scale, scale)

Apply a small jitter (1e-6) whenever you compute a test statistic. The jitter ensures that no two samples have exactly the same value of the test statistic.

Distribution of the test statistic and decision rule.

Draw n_simulations samples assuming H0 and compute the test statistic for each sample.
Display the distribution of the test statistic under H0 as a histogram.
Draw n_simulations samples assuming H1 and compute the test statistic for each sample.
Display the distribution of the test statistic under H1 as a histogram on the same graph.
Choose the decision rule (the rejection region) such that when H0 is true, type I errors occur with probability alpha.

30 trials of testing a d6 of unknown type.

Perform 30 trials: in each trial, the type of d6 is chosen randomly (chosen to be fair or rigged with equal probability), then a sample of size sample_size is generated, the test statistic is computed, then a decision to reject or not the null hypothesis is taken. For each trial, print:

true type of the d6
realization of the test statistic in this sample
p-value of the realization of the test statistic assuming H0 is true
Decision (either reject H0 or fail to reject H0)
Whether the decision was correct, a type I error, or a type II error

From trials to simulations

Using the samples drawn previously, compute the frequency of type I errors.
Using the samples drawn previously, compute the frequency of type II errors.
Give the empirical power of the test.

16610101020 Hypothesis testing — Mean vs Variance (OPTIONAL: vs log likelihood ratio)

The goal of this exercise is to consider a different type of rigged d6 and different test statistics. The null hypothesis is H0: the d6 is fair, and the alternative hypothesis is H1: the d6 is rigged.

Start your code with:

import random
random.seed(16610101020)

n_simulations = 3000
alpha = 0.10
sample_size = 50

interval_widths_fair = [12, 12, 12, 12, 12, 12]
interval_widths_rigged = [16, 6, 10, 10, 10, 20]

def jitter(x, scale):
    return x + random.uniform(-scale, scale)

Apply a small jitter (1e-6) whenever you compute a test statistic. The jitter ensures that no two samples have exactly the same value of the test statistic.

Follow the hypothesis testing procedure outlined below, using the sample mean as the test statistic (as in the lecture notes). Then, repeat the procedure, with a different test statistic: the sample variance. OPTIONAL: repeat the procedure with yet another test statistic: the log-likelihood ratio.

Remark: you can compute multiple test statistics on the same sample, rather than draw new samples when you change the test statistic. At the end, construct a dataframe with one row per test statistic, and with columns for:

decision rule (cutoff)
frequency of type I errors (empirical alpha)
frequency of type II errors
empirical power of the test statistic

Sort the dataframe by decreasing empirical power and print the dataframe.

Hypothesis testing procedure

Distribution of the test statistic and decision rule.

Draw n_simulations samples assuming H0 and compute the test statistic for each sample.
Display the distribution of the test statistic under H0 as a histogram.
Draw n_simulations samples assuming H1 and compute the test statistic for each sample.
Display the distribution of the test statistic under H1 as a histogram on the same graph.
Choose the decision rule (the rejection region) such that when H0 is true, type I errors occur with probability alpha.

From trials to simulations

Using the fair d6 samples drawn previously, compute the frequency of type I errors.
Using the rigged d6 samples drawn previously, compute the frequency of type II errors.
Give the empirical power of the test.

Exercise 1200110030 — Pooled virus testing

Source: Cyganowski, S., Kloeden, P. E., and Ombach, J. (2002), From Elementary Probability to Stochastic Differential Equations with Maple, Springer.

Start your code with:

import random
random.seed(1200110030)

n_simulations = 1400
p_infected = 0.01
n_population=1200

A population of n_population individuals is screened for a rare virus. The probability that a given individual is infected is p_infection, independent across individuals. The testing procedure is excellent (no false positives or false negatives), but expensive, so in an effort to reduce costs, the testing center practices pooled testing:

divide the population into groups of 30 individuals (pools).
take half the blood sample (5mL out of 10mL) of every individual within a group and mix them together.
test the mixed blood
if the mixed blood test is negative, the group is “clean”.
if the mixed blood test is positive, test separately and individually the remaining 5mL of blood from every individual in the group.

One simulation of one testing campaign.

Simulate one testing campaign: generate a sample with n_population individuals, each with a status (is_infected=True or is_infected=False) generated randomly using probability p_infected, and assume that individuals 0 through 29 form one group, then individuals 30 through 59 another group, etc.

Print:

the number of clean groups (tests performed of mixed blood which were negative)
the number of groups with at least one infected individual (tests performed of mixed blood which were positive)
the total number of group tests performed
the total number of individual tests performed
the total number of tests performed

Multiple simulations

Simulate a testing campaign n_simulations times: in each simulation, generate a new sample, then conduct the testing campaign. For each simulation, store the total number of tests performed.

Print descriptive statistics on the total number of tests:

average
minimum
25th percentile
median
75th percentile
maximum

Choosing the group size optimally

We are interested in whether dividing the population into groups of 30 is optimal, or whether a different number would reduce the expected cost of testing the entire population. We want to simulate testing campaigns assuming different group sizes: 2, 6, 8, 10, 12, 14, 16, 20, 25, 32, 40, 48.

For each possible group size, simulate the testing campaign n_simulations times.

Construct a dataframe where each row is a group size, and different columns display:

the minimum number of tests across simulations
the average number of tests across simulations
the maximum number of tests across simulations

Print the group size with the lowest average number of tests.

OPTIONAL: run the code again with a lower probability of infection p_infected = 0.005. How does the optimal group size change?

--- title: "PP5000 problem set 1: univariate-stats" # Title of the document as shown in the output format: html: code-fold: false # Always show code (no hiding or collapsing) code-tools: true # Enable copy / view tools for code blocks code-appendix: false # Do not automatically gather code at the end, we do that directly via the parser tool math: mathjax # forces MathJax, which is more robust with Greek letters. execute: echo: true # Display all code in the output eval: true # Execute all code blocks error: true # Show errors instead of stopping or hiding them --- **Choose one of the following 3 options** - Solve exercise 2611121 and exercise 10812121614 - Solve exercise 16610101020 - Solve exercise 1200110030 Items marked **OPTIONAL** are extensions, you can skip them and still receive full credit. The deadline is **Tuesday, September 01, 19:30**. You are allowed to discuss exercises and to collaborate, but this is an individual assignment. Upload a pdf document with your answers. Add brief comments, particularly if you reviewed your answers and found them suspicious or implausible. Don't overthink formatting and cosmetics: copy-pasting console output or taking screenshots is fine for this assignment. Acknowledge any assistance received, whether from humans or from AIs. Alongside your pdf answers, but separate from them, you must provide replication code sufficient to produce every result, table, plot... in your pdf answers. ___ ___ ___ ## Exercise 2611121: comparing `2d6+1` and `1d12+1` Start your code with: ```{python} import random import pandas as pd random.seed(2611121) n_realizations=4000 ``` Generate `n_realizations` realizations of: `2d6+1`, i.e. the result of rolling two d6, adding the faces, and adding 1 and `n_realizations` realizations of: `1d12+1`, i.e. the result of rolling 1 d12, and adding 1 Store the results in lists named respectively: `realizations_2d6_plus_1` and `realizations_1d12_plus_1` Put both lists in the same `pandas` dataframe. Produce a table with 2 columns named `2d6_plus_1` and `1d12_plus_1` and rows displaying: - minimum - maximum - mean - variance as well as 5 rows for the frequency of: - realization >= 8 - realization >= 9 - realization >= 10 - realization >= 11 - realization >= 12 ___ ___ ___ ## 10812121614 Hypothesis testing — Reverse the hypotheses The goal of this exercise is to follow the lecture notes' example of hypothesis testing, except we flip the hypotheses. The null hypothesis is now H0: the d6 is rigged, and the alternative hypothesis is now H1: the d6 is fair. The test statistic is the sample mean. Start your code with: ```{python} import random random.seed(10812121614) n_simulations = 2000 alpha = 0.10 sample_size = 50 interval_widths_fair = [12, 12, 12, 12, 12, 12] interval_widths_rigged = [10, 8, 12, 12, 16, 14] def jitter(x, scale): return x + random.uniform(-scale, scale) ``` Apply a small **jitter** (`1e-6`) whenever you compute a test statistic. The jitter ensures that no two samples have exactly the same value of the test statistic. ### Distribution of the test statistic and decision rule. - Draw `n_simulations` samples assuming H0 and compute the test statistic for each sample. - Display the distribution of the test statistic under H0 as a histogram. - Draw `n_simulations` samples assuming H1 and compute the test statistic for each sample. - Display the distribution of the test statistic under H1 as a histogram **on the same graph**. - Choose the decision rule (the rejection region) such that when H0 is true, type I errors occur with probability `alpha`. ___ ### 30 trials of testing a d6 of unknown type. Perform 30 **trials**: in each trial, the type of d6 is chosen randomly (chosen to be fair or rigged with equal probability), then a sample of size `sample_size` is generated, the test statistic is computed, then a decision to reject or not the null hypothesis is taken. For each trial, print: - true type of the d6 - realization of the test statistic in this sample - p-value of the realization of the test statistic assuming H0 is true - Decision (either `reject H0` or `fail to reject H0`) - Whether the decision was correct, a type I error, or a type II error ___ ### From trials to simulations - Using the samples drawn previously, compute the frequency of type I errors. - Using the samples drawn previously, compute the frequency of type II errors. - Give the empirical **power** of the test. ___ ___ ___ ## 16610101020 Hypothesis testing — Mean vs Variance (OPTIONAL: vs log likelihood ratio) The goal of this exercise is to consider a different type of rigged d6 and different test statistics. The null hypothesis is H0: the d6 is fair, and the alternative hypothesis is H1: the d6 is rigged. Start your code with: ```{python} import random random.seed(16610101020) n_simulations = 3000 alpha = 0.10 sample_size = 50 interval_widths_fair = [12, 12, 12, 12, 12, 12] interval_widths_rigged = [16, 6, 10, 10, 10, 20] def jitter(x, scale): return x + random.uniform(-scale, scale) ``` Apply a small **jitter** (`1e-6`) whenever you compute a test statistic. The jitter ensures that no two samples have exactly the same value of the test statistic. Follow the hypothesis testing procedure outlined below, using the **sample mean** as the test statistic (as in the lecture notes). Then, repeat the procedure, with a different test statistic: the **sample variance**. **OPTIONAL**: repeat the procedure with yet another test statistic: the **log-likelihood ratio**. **Remark**: you can compute multiple test statistics on the same sample, rather than draw new samples when you change the test statistic. At the end, construct a dataframe with one row per test statistic, and with columns for: - decision rule (cutoff) - frequency of type I errors (empirical alpha) - frequency of type II errors - empirical power of the test statistic Sort the dataframe by decreasing empirical power and print the dataframe. ::: {.callout-tip title="Hypothesis testing procedure"} ### Distribution of the test statistic and decision rule. - Draw `n_simulations` samples assuming H0 and compute the test statistic for each sample. - Display the distribution of the test statistic under H0 as a histogram. - Draw `n_simulations` samples assuming H1 and compute the test statistic for each sample. - Display the distribution of the test statistic under H1 as a histogram **on the same graph**. - Choose the decision rule (the rejection region) such that when H0 is true, type I errors occur with probability `alpha`. ___ ### From trials to simulations - Using the fair d6 samples drawn previously, compute the frequency of type I errors. - Using the rigged d6 samples drawn previously, compute the frequency of type II errors. - Give the empirical **power** of the test. ::: ___ ___ ___ ### Exercise 1200110030 — Pooled virus testing *Source:* Cyganowski, S., Kloeden, P. E., and Ombach, J. (2002), *From Elementary Probability to Stochastic Differential Equations with Maple*, Springer. Start your code with: ```{python} import random random.seed(1200110030) n_simulations = 1400 p_infected = 0.01 n_population=1200 ``` A population of `n_population` individuals is screened for a rare virus. The probability that a given individual is infected is `p_infection`, independent across individuals. The testing procedure is excellent (no false positives or false negatives), but expensive, so in an effort to reduce costs, the testing center practices pooled testing: 1. divide the population into groups of 30 individuals (pools). 2. take half the blood sample (5mL out of 10mL) of every individual within a group and mix them together. 3. test the mixed blood 4. if the mixed blood test is negative, the group is "clean". 5. if the mixed blood test is positive, test separately and individually the remaining 5mL of blood from every individual in the group. ### One simulation of one testing campaign. Simulate one testing campaign: generate a sample with `n_population` individuals, each with a status (`is_infected=True` or `is_infected=False`) generated randomly using probability `p_infected`, and assume that individuals 0 through 29 form one group, then individuals 30 through 59 another group, etc. Print: - the number of clean groups (tests performed of mixed blood which were negative) - the number of groups with at least one infected individual (tests performed of mixed blood which were positive) - the total number of group tests performed - the total number of individual tests performed - the total number of tests performed ___ ### Multiple simulations Simulate a testing campaign `n_simulations` times: in each simulation, generate a new sample, then conduct the testing campaign. For each simulation, store the **total number of tests** performed. Print descriptive statistics on the **total number of tests**: - average - minimum - 25th percentile - median - 75th percentile - maximum ___ ### Choosing the group size optimally We are interested in whether dividing the population into groups of 30 is optimal, or whether a different number would reduce the expected cost of testing the entire population. We want to simulate testing campaigns assuming different group sizes: 2, 6, 8, 10, 12, 14, 16, 20, 25, 32, 40, 48. For each possible group size, simulate the testing campaign `n_simulations` times. Construct a dataframe where each row is a group size, and different columns display: - the minimum number of tests across simulations - the average number of tests across simulations - the maximum number of tests across simulations Print the group size with the **lowest average number of tests**. **OPTIONAL**: run the code again with a lower probability of infection `p_infected = 0.005`. How does the optimal group size change?

Exercise 2611121: comparing 2d6+1 and 1d12+1

10812121614 Hypothesis testing — Reverse the hypotheses

Distribution of the test statistic and decision rule.

30 trials of testing a d6 of unknown type.

From trials to simulations

16610101020 Hypothesis testing — Mean vs Variance (OPTIONAL: vs log likelihood ratio)

Distribution of the test statistic and decision rule.

From trials to simulations

Exercise 1200110030 — Pooled virus testing

One simulation of one testing campaign.

Multiple simulations

Choosing the group size optimally

Exercise 2611121: comparing `2d6+1` and `1d12+1`