Problem Set 1: Frequentist hypothesis testing

Paired, exact, and contingency-table tests

Null-hypothesis significance testing

A null-hypothesis significance test has four parts. First, state a null hypothesis about a population parameter or a relation between variables. Second, define a statistic that measures departure from that hypothesis. Third, derive or select the statistic’s sampling distribution under the null hypothesis. Fourth, calculate the probability, under that distribution, of a result at least as discrepant as the observed result. This probability is the p value. It is not the probability that the null hypothesis is true.

This assignment applies that sequence to five tests: the paired t test, the Wilcoxon signed-rank test, the exact sign test, Fisher’s exact test, and the chi-squared test of independence. For each test, you will derive the statistic or its null probability and then implement the calculation in R. The acoustic data become relevant in Exercises 1–4. The morphological data become relevant in Exercises 5–6.

All six exercises are required. The points sum to 100.

Assignment files

Download the complete Problem Set 1 package. Extract the package and open ps1.ipynb from the extracted ps1-assignment directory. Select the R kernel. The package contains the R notebook, the assignment source, the README, the preparation scripts, and the redistributable teaching tables.

The files are also available separately:

Use <- for assignment, use snake_case object names, and record a seed before every random calculation. Every numerical answer must be produced by visible code.

Each derivation is divided into short steps. Complete only the displayed line or lines in each step. The surrounding equations supply the notation and the rule to apply. You do not need to reproduce them.

In an R block, replace every ___ with the requested expression or value. Code without a blank is supplied as scaffolding and should be retained.

Code

# Run the notebook from the extracted ps1-assignment directory.
stopifnot(
  file.exists("README.md"),
  file.exists("scripts/prepare_data.R")
)

if (!file.exists("data/hillenbrand_vowels.csv")) {
  source("scripts/prepare_data.R")
  prepare_hillenbrand_data("data/hillenbrand_vowels.csv")
}

stopifnot(
  file.exists("data/hillenbrand_vowels.csv"),
  file.exists("data/unimorph_genitive_lexemes.csv"),
  file.exists("data/unimorph_genitive_lexeme_population.csv"),
  file.exists("data/unimorph_plural_lexemes.csv")
)

library(dplyr)
library(ggplot2)
library(tidyr)

Exercise 1: Vowel formants and paired differences (10 points)

We’ll use hypothesis tests to analyze vowel formant frequencies. Formants are acoustic correlates of vowel quality in speech production. The first formant (F1) correlates with vowel height: lower F1 tends to correspond to a higher vowel. The second formant (F2) correlates with vowel backness: higher F2 tends to correspond to a fronter vowel.

The data for Exercises 1–4 come from Hillenbrand, Getty, Clark, and Wheeler’s study of American English vowels (1995, JASA). The preparation script extracts the phonTools::h95 teaching table. Each of the 45 adult men produced the same 12 vowels in an /hVd/ context. In particular, each man produced /i/ as in heed and /ɪ/ as in hid. Because the two measurements come from the same speaker, they form a paired observation.

Generate data/hillenbrand_vowels.csv by following the README. The columns used here are speaker_id, speaker_type, vowel, f1_hz, and f2_hz.

1.1 Prepare the paired table

Load the table. Retain rows for adult men and vowels i and I. Verify that each retained speaker has one row for each vowel. Join the two vowel tables by speaker_id and speaker_type. Create one row per speaker with f1_i_hz, f1_I_hz, and

d_s=\operatorname{F1}_{s,/i/}-\operatorname{F1}_{s,/ɪ/}.

Here d_s is the observed difference for speaker s.

Code

vowels <- read.csv("data/hillenbrand_vowels.csv")

# Inspect the available columns and vowel labels.
str(vowels)
sort(unique(vowels$vowel))

adult_men <- vowels |>
  filter(speaker_type == "man")

target_vowels <- adult_men |>
  filter(vowel %in% c("i", "I"))

# This table must contain one row for every speaker-vowel combination.
pair_counts <- target_vowels |>
  count(speaker_id, vowel)
stopifnot(all(pair_counts$n == 1))

# pivot_wider() joins the two measurements by the identifier columns.
paired_f1 <- target_vowels |>
  select(speaker_id, speaker_type, vowel, f1_hz) |>
  pivot_wider(
    names_from = vowel,
    values_from = f1_hz,
    names_glue = "f1_{vowel}_hz"
  ) |>
  mutate(difference_hz = f1_i_hz - f1_I_hz)

stopifnot(nrow(paired_f1) == 45, !anyNA(paired_f1$difference_hz))

1.2 Visualize vowel space

The F1–F2 vowel space provides an acoustic representation of vowel quality. Plot F2 on the horizontal axis and F1 on the vertical axis, with both axes reversed to mirror the conventional articulatory vowel space. Use the complete adult-male subset, not only /i/ and /ɪ/.

Code

ggplot(adult_men, aes(x = f2_hz, y = f1_hz, color = vowel)) +
  geom_point(alpha = .35) +
  stat_ellipse(level = .68) +
  scale_x_reverse() +
  scale_y_reverse() +
  labs(x = "F2 (Hz)", y = "F1 (Hz)", color = "Vowel")

The ellipse summarizes the center and spread of each vowel cloud under a bivariate normal approximation. It is a descriptive layer, not a confidence region for a vowel mean.

1.3 Summarize the differences

Calculate the number of speakers, the mean, the standard deviation, the median, and the number of negative, zero, and positive differences. Plot the differences with zero marked.

Code

difference_summary <- paired_f1 |>
  summarize(
    speaker_count = n(),
    mean_difference_hz = ___,
    sd_difference_hz = ___,
    median_difference_hz = ___,
    negative_count = ___,
    zero_count = ___,
    positive_count = ___
  )
difference_summary

ggplot(paired_f1, aes(x = difference_hz)) +
  geom_dotplot(binwidth = 5) +
  geom_vline(xintercept = 0, linetype = "dashed") +
  labs(x = "F1 difference: /i/ minus /ɪ/ (Hz)", y = "Count")

Report the requested values in the subtraction order /i/ minus /ɪ/.

Your result

Replace this sentence with the sample size and requested summaries.

Exercise 2: Paired t test (20 points)

Let D_s:\Omega\to\mathbb R be the random /i/ minus /ɪ/ difference for speaker s, and let \mu_D=\mathbb E[D_s]. Under the paired normal model,

D_s\mid\mu_D,\sigma_D\sim\mathcal N(\mu_D,\sigma_D^2).

Test H_0:\mu_D=0 against H_1:\mu_D\ne0.

Required background

Review the Student t distribution, null hypotheses and test statistics, paired observations, and inference for a paired mean difference before beginning the derivation below. The two vowel measurements are paired because each speaker produced both vowels. Thus the paired t test applies the one-sample t procedure to the within-speaker differences.

2.1 Derive the estimator and its standard error

Let S be the number of speakers and \overline D=S^{-1}\sum_{s=1}^{S}D_s. Assume that D_1,\ldots,D_S are independent and identically distributed with mean \mu_D and variance \sigma_D^2.

Step 1: expectation

Linearity of expectation moves the expectation through the finite sum. Replace each \mathbb E[D_s] with \mu_D, then simplify. Complete one line.

\mathbb E[\overline D] =\frac{1}{S}\sum_{s=1}^{S}\mathbb E[D_s] =\text{enter the simplification}.

Step 2: variance

Independence allows the variance of the sum to become the sum of the variances. Replace each \operatorname{Var}(D_s) with \sigma_D^2, then simplify. Complete one line.

\operatorname{Var}(\overline D) =\frac{1}{S^2}\sum_{s=1}^{S}\operatorname{Var}(D_s) =\text{enter the simplification}.

Step 3: estimated standard error

Take the square root of the variance from Step 2. Then replace the unknown \sigma_D with the sample standard deviation S_D. Complete two lines.

\begin{aligned} \operatorname{SE}(\overline D)&=\text{enter the population standard error},\\ \widehat{\operatorname{SE}}(\overline D)&=\text{enter the estimated standard error}. \end{aligned}

2.2 Derive the test statistic

Standardize the difference between \overline D and the null value. Because \sigma_D is estimated by S_D, the null reference distribution is a Student t distribution with S-1 degrees of freedom.

Step 1: random statistic

Insert the null value 0 in the numerator and the estimated standard error from Step 3 in the denominator. Complete one line.

Your LaTeX

\mathcal T=\text{enter the standardized random statistic}.

Step 2: observed statistic

Replace the random quantities \overline D and S_D with their observed values \overline d and s_d. Complete one line.

t_{\mathrm{obs}}=\text{enter the statistic as a function of }\overline d,s_d,\text{ and }S.

2.3 Derive the two-sided p value

Define the event containing null outcomes at least as far from zero as the observed statistic. Use full set-builder notation. Then use the cumulative distribution function F_{t_{S-1}} to write its probability.

Step 1: two-sided event

Use an absolute-value comparison so that the event contains both tails. Complete one line.

Your LaTeX

B_{\mathcal T}=\{\omega\in\Omega\mid\text{enter the two-sided condition}\}.

Step 2: event probability

Symmetry makes the two tail probabilities equal. Write twice the upper-tail probability in terms of F_{t_{S-1}}. Complete one line.

\mathbb P_0(B_{\mathcal T})=\text{enter the two-sided tail probability}.

2.4 Derive the confidence interval

A 95% confidence interval is obtained by inverting the central 95% t reference region. Subtract the critical value times the estimated standard error for the lower endpoint and add it for the upper endpoint. Complete two lines.

Your LaTeX

\begin{aligned} L&=\text{enter the lower endpoint},\\ U&=\text{enter the upper endpoint}. \end{aligned}

2.5 Implement the test

Implement the paired t-test calculation from scratch to verify R’s t.test() output. Without calling t.test(), calculate \overline d, s_d, the estimated standard error, t_{\mathrm{obs}}, the two-sided p value, and the 95% confidence interval.

Code

difference_vector <- paired_f1$difference_hz
speaker_count <- length(difference_vector)
mean_difference <- ___
sd_difference <- ___
estimated_se <- ___
degrees_of_freedom <- ___
t_observed <- ___

# pt(-abs(t_observed), ...) evaluates one tail. Double it.
p_value <- 2 * pt(-abs(t_observed), df = degrees_of_freedom)

t_critical <- qt(.975, df = degrees_of_freedom)
confidence_interval <- mean_difference + c(-1, 1) * ___

c(
  estimate = mean_difference,
  standard_error = estimated_se,
  t = t_observed,
  df = degrees_of_freedom,
  p_value = p_value,
  lower_95 = confidence_interval[1],
  upper_95 = confidence_interval[2]
)

Verify the result with t.test(..., paired = TRUE).

Code

software_result <- t.test(
  paired_f1$f1_i_hz,
  paired_f1$f1_I_hz,
  paired = TRUE
)
software_result

Report \overline d, the interval, t_{\mathrm{obs}}, the degrees of freedom, and the p value.

Your result

Replace this sentence with the requested test result.

2.6 Choose your own vowel pair

Select another pair of vowels and test whether their mean F1 or F2 values differ. Use the adult-male measurements and preserve the paired speaker design. State the two vowels, the acoustic measurement, the subtraction order, and the hypotheses. Construct the speaker-level differences, calculate the paired t test from the differences, and verify it with t.test(..., paired = TRUE).

Code

vowel_1 <- "___"
vowel_2 <- "___"
measurement <- "___"  # Enter either "f1_hz" or "f2_hz".

chosen_pair <- adult_men |>
  filter(vowel %in% c(vowel_1, vowel_2)) |>
  select(speaker_id, speaker_type, vowel, all_of(measurement)) |>
  pivot_wider(names_from = vowel, values_from = all_of(measurement))

# Subtract in the declared order.
chosen_pair$difference <-
  chosen_pair[[vowel_1]] - chosen_pair[[vowel_2]]

# Use chosen_pair$difference for the plot and manual calculation.
# Verify with t.test(chosen_pair$difference, mu = 0).

Report the estimate, 95% confidence interval, t statistic, degrees of freedom, and p value. Relate the sign of the estimate to the stated subtraction order and the selected formant.

Your result

Replace this sentence with the requested vowel comparison.

Exercise 3: Wilcoxon signed-rank test (10 points)

The paired t statistic uses the numerical magnitudes of the differences. The Wilcoxon signed-rank statistic replaces each nonzero magnitude |d_s| with its rank. It thus uses the order of the magnitudes rather than their spacing.

For S nonzero differences, rank |d_1|,\ldots,|d_S| from 1 to S, assigning the mean rank to tied magnitudes. Let r_s be the rank assigned to |d_s|. Define

w_+=\sum_{s=1}^{S}r_s\mathbb I(d_s>0),

where \mathbb I(d_s>0) equals one when d_s>0 and zero otherwise. Under the signed-rank null model, the differences are independent and symmetric about zero. Conditional on the observed absolute ranks, positive and negative signs are equally probable.

3.1 Construct the statistic

Use rank(abs(difference_hz), ties.method = "average") to calculate r_s. Calculate w_+ from the ranks and signs. Display speaker_id, difference_hz, the absolute difference, its rank, and its sign indicator.

Code

signed_ranks <- paired_f1 |>
  filter(difference_hz != 0) |>
  mutate(
    absolute_difference = abs(difference_hz),
    absolute_rank = rank(absolute_difference, ties.method = "average"),
    positive_indicator = as.integer(difference_hz > 0),
    positive_rank = ___
  )

w_plus <- sum(___)

signed_ranks |>
  select(
    speaker_id,
    difference_hz,
    absolute_difference,
    absolute_rank,
    positive_indicator
  )
w_plus

3.2 Derive the null mean

Let B_s:\Omega\to\{0,1\} be the random positive-sign indicator under the null model. Then \mathbb E[B_s]=.5 and \mathcal W_+=\sum_{s=1}^{S}r_sB_s.

Step 1: move the expectation through the sum

Use linearity of expectation and replace \mathbb E_0[B_s] with .5. Complete one line.

Your LaTeX

\mathbb E_0[\mathcal W_+] =\sum_{s=1}^{S}r_s\mathbb E_0[B_s] =\text{enter the expression after substitution}.

Step 2: sum the ranks

Use \sum_{s=1}^{S}r_s=S(S+1)/2 and simplify. Complete one line.

\mathbb E_0[\mathcal W_+]=\text{enter the null mean}.

3.3 Implement the test

Run wilcox.test() on the paired F1 measurements with paired = TRUE, alternative = "two.sided", exact = FALSE, and correct = FALSE. The observed integer-valued F1 measurements produce tied absolute differences, so this call uses the large-sample signed-rank reference approximation rather than exact enumeration.

Code

signed_rank_result <- wilcox.test(
  paired_f1$f1_i_hz,
  paired_f1$f1_I_hz,
  paired = TRUE,
  alternative = "two.sided",
  exact = FALSE,
  correct = FALSE
)

signed_rank_result
c(manual_w_plus = w_plus, software_v = unname(signed_rank_result$statistic))

Report w_+, the p value, and the null hypothesis tested by the signed-rank procedure.

Your result

Replace this sentence with the requested signed-rank result.

Exercise 4: Exact sign test (10 points)

The Bernoulli distribution represents whether each difference has the predicted sign. Define

Z_s(\omega)= \begin{cases} 1 & \text{if }D_s(\omega)<0,\\ 0 & \text{if }D_s(\omega)\ge0. \end{cases}

Let \pi_D=\mathbb P(Z_s=1). Test H_0:\pi_D=.5 against H_1:\pi_D>.5. There are no zero differences in these data.

4.1 Derive the null distribution and p value

Let \mathcal K=\sum_{s=1}^{S}Z_s and k_{\mathrm{obs}}=\sum_{s=1}^{S}z_s.

Step 1: null distribution

Under H_0, the Z_s are independent \operatorname{Bernoulli}(.5) variables. Use the binomial distribution for their sum. Complete one line.

Your LaTeX

\mathcal K\mid H_0\sim\text{enter the null distribution}.

Step 2: upper-tail event

The alternative predicts more negative differences than the null. Define the event containing counts at least as large as k_{\mathrm{obs}}. Complete one line.

B_{\mathcal K}=\{\omega\in\Omega\mid\text{enter the upper-tail condition}\}.

Step 3: exact probability

Insert the binomial probability mass function and sum from the observed count through S. Complete one line.

\mathbb P_0(B_{\mathcal K})=\text{enter the binomial sum from }k_{\mathrm{obs}}\text{ to }S.

4.2 Implement the test

Calculate k_{\mathrm{obs}}, k_{\mathrm{obs}}/S, and the upper-tail probability with pbinom(..., lower.tail = FALSE). Verify the calculation with binom.test(..., alternative = "greater") and report the one-sided 95% exact confidence bound.

Code

negative_indicator <- as.integer(paired_f1$difference_hz < 0)
k_observed <- sum(negative_indicator)
speaker_count <- length(negative_indicator)
observed_proportion <- k_observed / speaker_count

# pbinom(q, lower.tail = FALSE) returns P(K > q).
# Setting q = k_observed - 1 therefore returns P(K >= k_observed).
p_value <- pbinom(
  k_observed - 1,
  size = speaker_count,
  prob = .5,
  lower.tail = FALSE
)

sign_result <- binom.test(
  k_observed,
  speaker_count,
  p = .5,
  alternative = "greater"
)

c(
  negative_count = k_observed,
  observed_proportion = observed_proportion,
  lower_95 = sign_result$conf.int[1],
  p_value = p_value
)

Report the sign count, the observed proportion, the confidence bound, and the p value.

Your result

Replace this sentence with the requested sign-test result.

Exercise 5: Fisher’s exact test (25 points)

We’ll use Fisher’s exact test to analyze German genitive allomorphy, a case of morphological variation in which the genitive singular ending on masculine and neuter nouns alternates between -s and -es. For instance, the genitive of Tag “day” may be realized as Tages or Tags. The traditional pedagogical rule states that stems ending in sibilants require -es. Syllable count has also been proposed to matter, with monosyllabic nouns tending to favor -es and polysyllabic nouns tending to favor -s.

The data come from the German UniMorph paradigm inventory described by McCarthy et al. (2020). The preparation scripts reduce the genitive records to one row per lemma. They record syllable class, whether a lemma ends in one of the strings s, z, ß, sch, or x, and whether a genitive form ending in -es is attested. These are variables produced by the scripts. UniMorph does not supply them as linguistic analyses.

The raw paradigm may contain both -s and -es forms for one lemma. Thus es_attested records attestation in the extracted paradigm rather than an exclusive suffix choice.

The teaching table contains a reproducible simple random sample of 120 lexemes. Its relevant columns are lemma, syllable_class, ends_in_sibilant, and es_attested.

Required background

Review the hypergeometric distribution and Fisher’s exact test before beginning the manual calculation below. These notes develop the 2×2 contingency table, fixed margins, and the probability-ordering definition of a two-sided p value.

5.1 Load and inspect the lexeme table

Load data/unimorph_genitive_lexemes.csv. Verify that lemma is unique. Examine the levels and marginal distributions of syllable_class, ends_in_sibilant, and es_attested.

Code

genitive_lexemes <- read.csv("data/unimorph_genitive_lexemes.csv")

str(genitive_lexemes)
stopifnot(!anyDuplicated(genitive_lexemes$lemma))

count(genitive_lexemes, syllable_class)
count(genitive_lexemes, ends_in_sibilant)
count(genitive_lexemes, es_attested)

5.2 Test syllable count and genitive -es attestation

In German genitive allomorphy, is there an association between syllable class and the attestation of an -es form? Construct the two-by-two table for syllable_class by es_attested, with monosyllabic before polysyllabic and TRUE before FALSE. Plot the table with a mosaic plot or grouped bar chart. Run fisher.test() and calculate the sample odds ratio. With this order, the odds ratio places the odds of -es attestation among monosyllabic lexemes in the numerator and the corresponding odds among polysyllabic lexemes in the denominator.

Code

syllable_table <- with(
  genitive_lexemes,
  table(
    factor(
      syllable_class,
      levels = c("monosyllabic", "polysyllabic")
    ),
    factor(es_attested, levels = c(TRUE, FALSE))
  )
)
syllable_table
mosaicplot(syllable_table, color = TRUE)

syllable_fisher <- fisher.test(syllable_table)

a_syllable <- syllable_table[1, 1]
b_syllable <- syllable_table[1, 2]
c_syllable <- syllable_table[2, 1]
d_syllable <- syllable_table[2, 2]
syllable_sample_or <- ___

syllable_fisher
syllable_sample_or

Report the observed table, sample odds ratio, conditional estimate and interval returned by fisher.test(), and p value.

Your result

Replace this sentence with the requested syllable-class result.

5.3 Construct the sibilant table

Is there an association between the string-based sibilant classification and the attestation of an -es form? Construct the table for ends_in_sibilant by es_attested, with the character value "True" before "False" in the rows and the logical value TRUE before FALSE in the columns. With this order, the sample odds ratio places the odds of -es attestation among sibilant-final lexemes in the numerator and the corresponding odds among other lexemes in the denominator.

Code

sibilant_table <- with(
  genitive_lexemes,
  table(
    factor(ends_in_sibilant, levels = c("True", "False")),
    factor(es_attested, levels = c(TRUE, FALSE))
  )
)

sibilant_table
rowSums(sibilant_table)
colSums(sibilant_table)

Report the four cells and both sets of margins.

Your result

Replace this sentence with the observed sibilant table and its margins.

5.4 Derive the conditional null probability

Fisher’s exact test tests independence in a two-by-two table while holding both observed margins fixed. Write the observed cells in row-major order as a,b,c,d. Let A be the random upper-left cell count among all tables with those margins. Under the null hypothesis of independence, A has a hypergeometric distribution.

Derive the probability of the observed table by counting the ways to select its first-column total from the two rows and dividing by the number of selections from all observations.

Step 1: selections within each row

The first row contains a+b observations and contributes a observations to the first column. The second row contains c+d observations and contributes c. Write the product of the two counts. Complete one line.

Your LaTeX

N_{\mathrm{observed}}=\text{enter the product of two binomial coefficients}.

Step 2: all selections with the fixed margins

There are a+b+c+d observations in total, and the first column contains a+c of them. Write the number of possible first-column selections. Complete one line.

N_{\mathrm{all}}=\text{enter one binomial coefficient}.

Step 3: conditional probability

Divide the count in Step 1 by the count in Step 2. Complete one line.

\mathbb P_0(A=a)=\text{enter }N_{\mathrm{observed}}/N_{\mathrm{all}}\text{ in terms of }a,b,c,d.

5.5 Implement the sibilant test

Calculate the sample odds ratio ad/(bc) and the observed-table probability with dhyper(). Run fisher.test() to obtain the two-sided p value and conditional confidence interval. Fisher’s two-sided p value sums the probabilities of all margin-compatible tables whose null probabilities are no greater than the observed table’s probability.

Code

a <- sibilant_table[1, 1]
b <- sibilant_table[1, 2]
c_cell <- sibilant_table[2, 1]
d_cell <- sibilant_table[2, 2]

sibilant_sample_or <- ___

# dhyper(x, m, n, k) uses the upper-left cell, the two row totals,
# and the first-column total in that order.
observed_table_probability <- dhyper(
  x = a,
  m = a + b,
  n = c_cell + d_cell,
  k = a + c_cell
)

sibilant_fisher <- fisher.test(sibilant_table)

tibble(
  predictor = c("syllable class", "sibilant ending"),
  sample_odds_ratio = c(syllable_sample_or, sibilant_sample_or)
)

observed_table_probability
sibilant_fisher

Report the odds-ratio orientation, the sample odds ratio, the conditional estimate returned by fisher.test(), its confidence interval, and the p value. Put the sibilant and syllable-class sample odds ratios in the same table.

Your result

Replace this sentence with the requested sibilant result and odds-ratio comparison.

5.6 Estimate sampling distributions across lexicon sizes

A sampling distribution is the distribution of a statistic across repeated samples from a fixed population. For this exercise, treat data/unimorph_genitive_lexeme_population.csv as the finite population.

For each n\in\{12,20,40,80,120\}, draw 500 simple random samples without replacement. Use one recorded seed. For every sample, record the odds ratio, Fisher p value, confidence-interval endpoints, whether any cell is zero, and whether either margin is zero. Use explicit factor levels when constructing every table.

A zero cell is a table entry equal to zero. A zero margin occurs when an entire row or column is absent. fisher.test() should not be called when a margin is zero. Return NA_real_ for its p value and interval in that case. A sample odds ratio may be zero or infinite when a cell is zero. Retain that value and exclude it only from summaries that require a finite logarithm.

Code

genitive_population <- read.csv(
  "data/unimorph_genitive_lexeme_population.csv"
)

analyze_lexeme_sample <- function(sample_size, replicate_id) {
  sampled_lexemes <- genitive_population |>
    slice_sample(n = sample_size, replace = FALSE)

  sampled_table <- with(
    sampled_lexemes,
    table(
      factor(ends_in_sibilant, levels = c("True", "False")),
      factor(es_attested, levels = c(TRUE, FALSE))
    )
  )

  zero_cell <- any(sampled_table == 0)
  zero_margin <- any(rowSums(sampled_table) == 0) ||
    any(colSums(sampled_table) == 0)

  a <- sampled_table[1, 1]
  b <- sampled_table[1, 2]
  c_cell <- sampled_table[2, 1]
  d_cell <- sampled_table[2, 2]

  # Preserve zero and infinite sample odds ratios.
  sample_or <- if (b * c_cell == 0) {
    if (a * d_cell > 0) Inf else NA_real_
  } else {
    ___
  }

  if (zero_margin) {
    fisher_p <- NA_real_
    ci_lower <- NA_real_
    ci_upper <- NA_real_
  } else {
    sampled_fisher <- fisher.test(sampled_table)
    fisher_p <- sampled_fisher$p.value
    ci_lower <- sampled_fisher$conf.int[1]
    ci_upper <- sampled_fisher$conf.int[2]
  }

  tibble(
    sample_size = sample_size,
    replicate_id = replicate_id,
    zero_cell = zero_cell,
    zero_margin = zero_margin,
    odds_ratio = sample_or,
    fisher_p = fisher_p,
    ci_lower = ci_lower,
    ci_upper = ci_upper
  )
}

set.seed(___)
sampling_plan <- expand_grid(
  sample_size = c(12, 20, 40, 80, 120),
  replicate_id = seq_len(500)
)

sampling_results <- bind_rows(Map(
  analyze_lexeme_sample,
  sampling_plan$sample_size,
  sampling_plan$replicate_id
))

For each n, calculate the proportion of samples with a zero cell, the proportion with a zero margin, the proportion with a nonfinite odds ratio, the median finite log odds ratio, and the median finite confidence-interval width on the log scale. Plot the sampling distributions of the finite log odds ratios by n.

Code

sampling_summary <- sampling_results |>
  mutate(
    finite_log_or = is.finite(odds_ratio) & odds_ratio > 0,
    finite_log_interval =
      is.finite(ci_lower) & ci_lower > 0 &
      is.finite(ci_upper) & ci_upper > 0
  ) |>
  group_by(sample_size) |>
  summarize(
    zero_cell_rate = mean(zero_cell),
    zero_margin_rate = mean(zero_margin),
    nonfinite_or_rate = mean(!finite_log_or),
    median_log_or = median(log(odds_ratio[finite_log_or])),
    median_log_interval_width = median(
      log(ci_upper[finite_log_interval]) -
        log(ci_lower[finite_log_interval])
    ),
    .groups = "drop"
  )

sampling_summary

sampling_results |>
  filter(is.finite(odds_ratio), odds_ratio > 0) |>
  ggplot(aes(x = factor(sample_size), y = log(odds_ratio))) +
  geom_boxplot() +
  labs(x = "Lexemes sampled", y = "Finite log sample odds ratio")

Report the summary table.

Your result

Replace this sentence with the sample-size summary table.

Exercise 6: Chi-squared test of independence (25 points)

We’ll use the chi-squared test to analyze German plural allomorphy. German plural marking combines several suffix classes with possible stem changes. The classes represented in this teaching table include -e as in Tag → Tage, -en as in Frau → Frauen, -er as in Kind → Kinder, -s as in Auto → Autos, and forms for which the preparation rule finds none of these overt suffixes.

The second UniMorph table has one row per noun lemma. The preparation rule retains one plural form and classifies its ending as e, en, er, s, or ø. Here ø means that the rule found none of the listed overt suffixes. It is a label produced by the script, not a claim that the plural has no exponence or that it contains umlaut.

The linguistic question is whether grammatical gender and the derived plural-suffix class are independent in the lexeme inventory represented by the table. This creates a contingency table larger than 2×2.

Required background

Review the chi-squared test of independence and the chi-squared distribution before beginning the derivations below. These notes develop the independence factorization, expected counts, Pearson’s statistic, degrees of freedom, and the large-sample reference distribution.

Load data/unimorph_plural_lexemes.csv. The relevant columns are lemma, gender, and plural_suffix. Construct the observed three-by-five table with gender ordered FEM, MASC, NEUT and suffix ordered e, en, er, s, ø.

Code

plural_lexemes <- read.csv("data/unimorph_plural_lexemes.csv")
stopifnot(!anyDuplicated(plural_lexemes$lemma))

plural_table <- with(
  plural_lexemes,
  table(
    factor(gender, levels = c("FEM", "MASC", "NEUT")),
    factor(plural_suffix, levels = c("e", "en", "er", "s", "ø"))
  )
)
plural_table

6.1 Derive the expected counts

Under independence, the expected count in cell (i,j) is determined by its row total, column total, and the grand total. Follow the derivation in the chi-squared-test notes.

Step 1: independence probability

Under the null hypothesis, the joint cell probability is the product of its row and column probabilities. Replace both probabilities with the corresponding observed marginal proportions. Complete one line.

Your LaTeX

\widehat p_{ij} =\widehat p_{i+}\widehat p_{+j} =\text{enter the product of observed marginal proportions}.

Step 2: expected count

Multiply the probability in Step 1 by the grand total O_{++} and cancel one factor of O_{++}. Complete one line.

E_{ij}=\text{enter the expected-count formula}.

Step 3: feminine-by-en cell

Read the feminine row total, the -en column total, and the grand total from plural_table. Substitute those three numbers. Complete one line.

E_{\mathrm{FEM},\mathrm{en}} =\text{enter the numerical substitution and result}.

6.2 Derive the statistic and reference distribution

Pearson’s statistic sums the standardized squared difference between every observed and expected count.

Step 1: statistic

Sum the cell contribution (O_{ij}-E_{ij})^2/E_{ij} over all rows and columns. Complete one line.

Your LaTeX

\mathcal X^2=\text{enter Pearson's statistic}.

Step 2: degrees of freedom

An R\times C table has RC cells. Fixed row and column margins leave (R-1)(C-1) free cells. Substitute R=3 and C=5. Complete two lines.

\begin{aligned} \nu&=(R-1)(C-1),\\ \nu&=\text{enter the numerical degrees of freedom}. \end{aligned}

Step 3: upper-tail event

Larger values of the statistic represent larger departures from independence. Complete one line.

B_{\mathcal X^2}=\{\omega\in\Omega\mid\text{enter the upper-tail condition}\}.

Step 4: reference probability

Write the upper-tail probability using the cumulative distribution function F_{\chi^2_\nu}. Complete one line.

\mathbb P_0(B_{\mathcal X^2})=\text{enter the chi-squared tail probability}.

6.3 Implement the test

Compute the chi-squared statistic manually to verify R’s output. Without calling chisq.test(), calculate the complete expected table, every cell contribution, x^2_{\mathrm{obs}}, the degrees of freedom, and the p value. Verify that every expected count is at least five.

Code

row_totals <- rowSums(plural_table)
column_totals <- colSums(plural_table)
grand_total <- sum(plural_table)

expected_counts <- outer(row_totals, column_totals) / ___
cell_contributions <- ___
x_squared_observed <- sum(cell_contributions)
degrees_of_freedom <- ___
p_value <- pchisq(
  x_squared_observed,
  df = degrees_of_freedom,
  lower.tail = FALSE
)

expected_counts
cell_contributions
c(
  x_squared = x_squared_observed,
  df = degrees_of_freedom,
  p_value = p_value,
  minimum_expected = min(expected_counts)
)

Verify the statistic, degrees of freedom, p value, and expected table with chisq.test(..., correct = FALSE). Calculate the Pearson residuals (O_{ij}-E_{ij})/\sqrt{E_{ij}} and identify the two cells with the largest absolute residuals.

Code

software_result <- chisq.test(plural_table, correct = FALSE)
software_result

stopifnot(all(expected_counts >= 5))

pearson_residuals <- ___
residual_table <- as.data.frame(as.table(pearson_residuals)) |>
  mutate(absolute_residual = abs(Freq)) |>
  arrange(desc(absolute_residual))

slice_head(residual_table, n = 2)

Report x^2_{\mathrm{obs}}, the degrees of freedom, the p value, the minimum expected count, and the two largest absolute residuals with their signs.

Your result

Replace this sentence with the requested chi-squared-test result.

Submission checklist

The notebook runs from a fresh R session after the data-preparation steps.
Every derivation is entered in its LaTeX block.
Every computation is entered in the R block immediately following its instruction.
Every requested result is entered in its result block.
The paired table is joined by identifiers rather than aligned by row position.
Every hypothesis names its parameter or independence relation.
Every p value is defined relative to its null reference distribution.
Every random calculation has a recorded seed.

--- title: "Problem Set 1: Frequentist hypothesis testing" subtitle: "Paired, exact, and contingency-table tests" format: html: toc: true execute: enabled: false --- ## Null-hypothesis significance testing A [null-hypothesis significance test](../../statistical-inference/null-hypotheses-and-test-statistics.qmd) has four parts. First, state a null hypothesis about a population parameter or a relation between variables. Second, define a statistic that measures departure from that hypothesis. Third, derive or select the statistic's sampling distribution under the null hypothesis. Fourth, calculate the probability, under that distribution, of a result at least as discrepant as the observed result. This probability is the p value. It is not the probability that the null hypothesis is true. This assignment applies that sequence to five tests: the paired t test, the Wilcoxon signed-rank test, the exact sign test, Fisher's exact test, and the chi-squared test of independence. For each test, you will derive the statistic or its null probability and then implement the calculation in R. The acoustic data become relevant in Exercises 1–4. The morphological data become relevant in Exercises 5–6. All six exercises are required. The points sum to 100. ::: {.callout-note title="Assignment files"} [Download the complete Problem Set 1 package](../../downloads/ps1-assignment.zip){download="ps1-assignment.zip"}. Extract the package and open `ps1.ipynb` from the extracted `ps1-assignment` directory. Select the R kernel. The package contains the R notebook, the assignment source, the README, the preparation scripts, and the redistributable teaching tables. The files are also available separately: - [R notebook](ps1.ipynb){download="ps1.ipynb"} - [README and setup instructions](README.md) - [`prepare_data.R`](scripts/prepare_data.R){download="prepare_data.R"} - [`fetch_unimorph_data.py`](scripts/fetch_unimorph_data.py){download="fetch_unimorph_data.py"} - [`prepare_unimorph_tables.R`](scripts/prepare_unimorph_tables.R){download="prepare_unimorph_tables.R"} ::: Use `<-` for assignment, use `snake_case` object names, and record a seed before every random calculation. Every numerical answer must be produced by visible code. Each derivation is divided into short steps. Complete only the displayed line or lines in each step. The surrounding equations supply the notation and the rule to apply. You do not need to reproduce them. In an R block, replace every `___` with the requested expression or value. Code without a blank is supplied as scaffolding and should be retained. ```{r} # Run the notebook from the extracted ps1-assignment directory. stopifnot( file.exists("README.md"), file.exists("scripts/prepare_data.R") ) if (!file.exists("data/hillenbrand_vowels.csv")) { source("scripts/prepare_data.R") prepare_hillenbrand_data("data/hillenbrand_vowels.csv") } stopifnot( file.exists("data/hillenbrand_vowels.csv"), file.exists("data/unimorph_genitive_lexemes.csv"), file.exists("data/unimorph_genitive_lexeme_population.csv"), file.exists("data/unimorph_plural_lexemes.csv") ) library(dplyr) library(ggplot2) library(tidyr) ``` ## Exercise 1: Vowel formants and paired differences (10 points) We'll use hypothesis tests to analyze **vowel formant frequencies**. Formants are acoustic correlates of vowel quality in speech production. The first formant (F1) correlates with vowel height: lower F1 tends to correspond to a higher vowel. The second formant (F2) correlates with vowel backness: higher F2 tends to correspond to a fronter vowel. The data for Exercises 1–4 come from Hillenbrand, Getty, Clark, and Wheeler's study of American English vowels ([1995, JASA](https://doi.org/10.1121/1.411872)). The preparation script extracts the `phonTools::h95` teaching table. Each of the 45 adult men produced the same 12 vowels in an /hVd/ context. In particular, each man produced /i/ as in *heed* and /ɪ/ as in *hid*. Because the two measurements come from the same speaker, they form a [paired observation](../../statistical-inference/paired-observations.qmd). Generate `data/hillenbrand_vowels.csv` by following the [README](README.md). The columns used here are `speaker_id`, `speaker_type`, `vowel`, `f1_hz`, and `f2_hz`. ### 1.1 Prepare the paired table Load the table. Retain rows for adult men and vowels `i` and `I`. Verify that each retained speaker has one row for each vowel. Join the two vowel tables by `speaker_id` and `speaker_type`. Create one row per speaker with `f1_i_hz`, `f1_I_hz`, and $$ d_s=\operatorname{F1}_{s,/i/}-\operatorname{F1}_{s,/ɪ/}. $$ Here $d_s$ is the observed difference for speaker $s$. ```{r} vowels <- read.csv("data/hillenbrand_vowels.csv") # Inspect the available columns and vowel labels. str(vowels) sort(unique(vowels$vowel)) adult_men <- vowels |> filter(speaker_type == "man") target_vowels <- adult_men |> filter(vowel %in% c("i", "I")) # This table must contain one row for every speaker-vowel combination. pair_counts <- target_vowels |> count(speaker_id, vowel) stopifnot(all(pair_counts$n == 1)) # pivot_wider() joins the two measurements by the identifier columns. paired_f1 <- target_vowels |> select(speaker_id, speaker_type, vowel, f1_hz) |> pivot_wider( names_from = vowel, values_from = f1_hz, names_glue = "f1_{vowel}_hz" ) |> mutate(difference_hz = f1_i_hz - f1_I_hz) stopifnot(nrow(paired_f1) == 45, !anyNA(paired_f1$difference_hz)) ``` ### 1.2 Visualize vowel space The **F1–F2 vowel space** provides an acoustic representation of vowel quality. Plot F2 on the horizontal axis and F1 on the vertical axis, with both axes reversed to mirror the conventional articulatory vowel space. Use the complete adult-male subset, not only /i/ and /ɪ/. ```{r} ggplot(adult_men, aes(x = f2_hz, y = f1_hz, color = vowel)) + geom_point(alpha = .35) + stat_ellipse(level = .68) + scale_x_reverse() + scale_y_reverse() + labs(x = "F2 (Hz)", y = "F1 (Hz)", color = "Vowel") ``` The ellipse summarizes the center and spread of each vowel cloud under a bivariate normal approximation. It is a descriptive layer, not a confidence region for a vowel mean. ### 1.3 Summarize the differences Calculate the number of speakers, the mean, the standard deviation, the median, and the number of negative, zero, and positive differences. Plot the differences with zero marked. ```{r} difference_summary <- paired_f1 |> summarize( speaker_count = n(), mean_difference_hz = ___, sd_difference_hz = ___, median_difference_hz = ___, negative_count = ___, zero_count = ___, positive_count = ___ ) difference_summary ggplot(paired_f1, aes(x = difference_hz)) + geom_dotplot(binwidth = 5) + geom_vline(xintercept = 0, linetype = "dashed") + labs(x = "F1 difference: /i/ minus /ɪ/ (Hz)", y = "Count") ``` Report the requested values in the subtraction order /i/ minus /ɪ/. ::: {.callout-note title="Your result"} Replace this sentence with the sample size and requested summaries. ::: ## Exercise 2: Paired t test (20 points) Let $D_s:\Omega\to\mathbb R$ be the random /i/ minus /ɪ/ difference for speaker $s$, and let $\mu_D=\mathbb E[D_s]$. Under the [paired normal model](../../statistical-inference/paired-mean-inference.qmd), $$ D_s\mid\mu_D,\sigma_D\sim\mathcal N(\mu_D,\sigma_D^2). $$ Test $H_0:\mu_D=0$ against $H_1:\mu_D\ne0$. ::: {.callout-note title="Required background"} Review [the Student t distribution](../../random-variables-and-distributions/student-t-distribution.qmd), [null hypotheses and test statistics](../../statistical-inference/null-hypotheses-and-test-statistics.qmd), [paired observations](../../statistical-inference/paired-observations.qmd), and [inference for a paired mean difference](../../statistical-inference/paired-mean-inference.qmd) before beginning the derivation below. The two vowel measurements are paired because each speaker produced both vowels. Thus the paired t test applies the one-sample t procedure to the within-speaker differences. ::: ### 2.1 Derive the estimator and its standard error Let $S$ be the number of speakers and $\overline D=S^{-1}\sum_{s=1}^{S}D_s$. Assume that $D_1,\ldots,D_S$ are independent and identically distributed with mean $\mu_D$ and variance $\sigma_D^2$. #### Step 1: expectation [Linearity of expectation](../../random-variables-and-distributions/linearity-of-expectation.qmd) moves the expectation through the finite sum. Replace each $\mathbb E[D_s]$ with $\mu_D$, then simplify. Complete one line. $$ \mathbb E[\overline D] =\frac{1}{S}\sum_{s=1}^{S}\mathbb E[D_s] =\text{enter the simplification}. $$ #### Step 2: variance Independence allows the variance of the sum to become the sum of the variances. Replace each $\operatorname{Var}(D_s)$ with $\sigma_D^2$, then simplify. Complete one line. $$ \operatorname{Var}(\overline D) =\frac{1}{S^2}\sum_{s=1}^{S}\operatorname{Var}(D_s) =\text{enter the simplification}. $$ #### Step 3: estimated standard error Take the square root of the variance from Step 2. Then replace the unknown $\sigma_D$ with the sample standard deviation $S_D$. Complete two lines. $$ \begin{aligned} \operatorname{SE}(\overline D)&=\text{enter the population standard error},\\ \widehat{\operatorname{SE}}(\overline D)&=\text{enter the estimated standard error}. \end{aligned} $$ ### 2.2 Derive the test statistic Standardize the difference between $\overline D$ and the null value. Because $\sigma_D$ is estimated by $S_D$, the null reference distribution is a [Student t distribution](../../random-variables-and-distributions/student-t-distribution.qmd) with $S-1$ degrees of freedom. #### Step 1: random statistic Insert the null value $0$ in the numerator and the estimated standard error from Step 3 in the denominator. Complete one line. **Your LaTeX** $$ \mathcal T=\text{enter the standardized random statistic}. $$ #### Step 2: observed statistic Replace the random quantities $\overline D$ and $S_D$ with their observed values $\overline d$ and $s_d$. Complete one line. $$ t_{\mathrm{obs}}=\text{enter the statistic as a function of }\overline d,s_d,\text{ and }S. $$ ### 2.3 Derive the two-sided p value Define the event containing null outcomes at least as far from zero as the observed statistic. Use full set-builder notation. Then use the cumulative distribution function $F_{t_{S-1}}$ to write its probability. #### Step 1: two-sided event Use an absolute-value comparison so that the event contains both tails. Complete one line. **Your LaTeX** $$ B_{\mathcal T}=\{\omega\in\Omega\mid\text{enter the two-sided condition}\}. $$ #### Step 2: event probability Symmetry makes the two tail probabilities equal. Write twice the upper-tail probability in terms of $F_{t_{S-1}}$. Complete one line. $$ \mathbb P_0(B_{\mathcal T})=\text{enter the two-sided tail probability}. $$ ### 2.4 Derive the confidence interval A 95% [confidence interval](../../statistical-inference/confidence-intervals.qmd) is obtained by inverting the central 95% t reference region. Subtract the critical value times the estimated standard error for the lower endpoint and add it for the upper endpoint. Complete two lines. **Your LaTeX** $$ \begin{aligned} L&=\text{enter the lower endpoint},\\ U&=\text{enter the upper endpoint}. \end{aligned} $$ ### 2.5 Implement the test Implement the paired t-test calculation from scratch to verify R's `t.test()` output. Without calling `t.test()`, calculate $\overline d$, $s_d$, the estimated standard error, $t_{\mathrm{obs}}$, the two-sided p value, and the 95% confidence interval. ```{r} difference_vector <- paired_f1$difference_hz speaker_count <- length(difference_vector) mean_difference <- ___ sd_difference <- ___ estimated_se <- ___ degrees_of_freedom <- ___ t_observed <- ___ # pt(-abs(t_observed), ...) evaluates one tail. Double it. p_value <- 2 * pt(-abs(t_observed), df = degrees_of_freedom) t_critical <- qt(.975, df = degrees_of_freedom) confidence_interval <- mean_difference + c(-1, 1) * ___ c( estimate = mean_difference, standard_error = estimated_se, t = t_observed, df = degrees_of_freedom, p_value = p_value, lower_95 = confidence_interval[1], upper_95 = confidence_interval[2] ) ``` Verify the result with `t.test(..., paired = TRUE)`. ```{r} software_result <- t.test( paired_f1$f1_i_hz, paired_f1$f1_I_hz, paired = TRUE ) software_result ``` Report $\overline d$, the interval, $t_{\mathrm{obs}}$, the degrees of freedom, and the p value. ::: {.callout-note title="Your result"} Replace this sentence with the requested test result. ::: ### 2.6 Choose your own vowel pair Select another pair of vowels and test whether their mean F1 or F2 values differ. Use the adult-male measurements and preserve the paired speaker design. State the two vowels, the acoustic measurement, the subtraction order, and the hypotheses. Construct the speaker-level differences, calculate the paired t test from the differences, and verify it with `t.test(..., paired = TRUE)`. ```{r} vowel_1 <- "___" vowel_2 <- "___" measurement <- "___" # Enter either "f1_hz" or "f2_hz". chosen_pair <- adult_men |> filter(vowel %in% c(vowel_1, vowel_2)) |> select(speaker_id, speaker_type, vowel, all_of(measurement)) |> pivot_wider(names_from = vowel, values_from = all_of(measurement)) # Subtract in the declared order. chosen_pair$difference <- chosen_pair[[vowel_1]] - chosen_pair[[vowel_2]] # Use chosen_pair$difference for the plot and manual calculation. # Verify with t.test(chosen_pair$difference, mu = 0). ``` Report the estimate, 95% confidence interval, t statistic, degrees of freedom, and p value. Relate the sign of the estimate to the stated subtraction order and the selected formant. ::: {.callout-note title="Your result"} Replace this sentence with the requested vowel comparison. ::: ## Exercise 3: Wilcoxon signed-rank test (10 points) The paired t statistic uses the numerical magnitudes of the differences. The Wilcoxon signed-rank statistic replaces each nonzero magnitude $|d_s|$ with its rank. It thus uses the order of the magnitudes rather than their spacing. For $S$ nonzero differences, rank $|d_1|,\ldots,|d_S|$ from 1 to $S$, assigning the mean rank to tied magnitudes. Let $r_s$ be the rank assigned to $|d_s|$. Define $$ w_+=\sum_{s=1}^{S}r_s\mathbb I(d_s>0), $$ where $\mathbb I(d_s>0)$ equals one when $d_s>0$ and zero otherwise. Under the signed-rank null model, the differences are independent and symmetric about zero. Conditional on the observed absolute ranks, positive and negative signs are equally probable. ### 3.1 Construct the statistic Use `rank(abs(difference_hz), ties.method = "average")` to calculate $r_s$. Calculate $w_+$ from the ranks and signs. Display `speaker_id`, `difference_hz`, the absolute difference, its rank, and its sign indicator. ```{r} signed_ranks <- paired_f1 |> filter(difference_hz != 0) |> mutate( absolute_difference = abs(difference_hz), absolute_rank = rank(absolute_difference, ties.method = "average"), positive_indicator = as.integer(difference_hz > 0), positive_rank = ___ ) w_plus <- sum(___) signed_ranks |> select( speaker_id, difference_hz, absolute_difference, absolute_rank, positive_indicator ) w_plus ``` ### 3.2 Derive the null mean Let $B_s:\Omega\to\{0,1\}$ be the random positive-sign indicator under the null model. Then $\mathbb E[B_s]=.5$ and $\mathcal W_+=\sum_{s=1}^{S}r_sB_s$. #### Step 1: move the expectation through the sum Use linearity of expectation and replace $\mathbb E_0[B_s]$ with $.5$. Complete one line. **Your LaTeX** $$ \mathbb E_0[\mathcal W_+] =\sum_{s=1}^{S}r_s\mathbb E_0[B_s] =\text{enter the expression after substitution}. $$ #### Step 2: sum the ranks Use $\sum_{s=1}^{S}r_s=S(S+1)/2$ and simplify. Complete one line. $$ \mathbb E_0[\mathcal W_+]=\text{enter the null mean}. $$ ### 3.3 Implement the test Run `wilcox.test()` on the paired F1 measurements with `paired = TRUE`, `alternative = "two.sided"`, `exact = FALSE`, and `correct = FALSE`. The observed integer-valued F1 measurements produce tied absolute differences, so this call uses the large-sample signed-rank reference approximation rather than exact enumeration. ```{r} signed_rank_result <- wilcox.test( paired_f1$f1_i_hz, paired_f1$f1_I_hz, paired = TRUE, alternative = "two.sided", exact = FALSE, correct = FALSE ) signed_rank_result c(manual_w_plus = w_plus, software_v = unname(signed_rank_result$statistic)) ``` Report $w_+$, the p value, and the null hypothesis tested by the signed-rank procedure. ::: {.callout-note title="Your result"} Replace this sentence with the requested signed-rank result. ::: ## Exercise 4: Exact sign test (10 points) The [Bernoulli distribution](../../random-variables-and-distributions/bernoulli-distribution.qmd) represents whether each difference has the predicted sign. Define $$ Z_s(\omega)= \begin{cases} 1 & \text{if }D_s(\omega)<0,\\ 0 & \text{if }D_s(\omega)\ge0. \end{cases} $$ Let $\pi_D=\mathbb P(Z_s=1)$. Test $H_0:\pi_D=.5$ against $H_1:\pi_D>.5$. There are no zero differences in these data. ### 4.1 Derive the null distribution and p value Let $\mathcal K=\sum_{s=1}^{S}Z_s$ and $k_{\mathrm{obs}}=\sum_{s=1}^{S}z_s$. #### Step 1: null distribution Under $H_0$, the $Z_s$ are independent $\operatorname{Bernoulli}(.5)$ variables. Use the [binomial distribution](../../random-variables-and-distributions/binomial-distribution.qmd) for their sum. Complete one line. **Your LaTeX** $$ \mathcal K\mid H_0\sim\text{enter the null distribution}. $$ #### Step 2: upper-tail event The alternative predicts more negative differences than the null. Define the event containing counts at least as large as $k_{\mathrm{obs}}$. Complete one line. $$ B_{\mathcal K}=\{\omega\in\Omega\mid\text{enter the upper-tail condition}\}. $$ #### Step 3: exact probability Insert the binomial probability mass function and sum from the observed count through $S$. Complete one line. $$ \mathbb P_0(B_{\mathcal K})=\text{enter the binomial sum from }k_{\mathrm{obs}}\text{ to }S. $$ ### 4.2 Implement the test Calculate $k_{\mathrm{obs}}$, $k_{\mathrm{obs}}/S$, and the upper-tail probability with `pbinom(..., lower.tail = FALSE)`. Verify the calculation with `binom.test(..., alternative = "greater")` and report the one-sided 95% [exact confidence bound](../../statistical-inference/exact-binomial-intervals.qmd). ```{r} negative_indicator <- as.integer(paired_f1$difference_hz < 0) k_observed <- sum(negative_indicator) speaker_count <- length(negative_indicator) observed_proportion <- k_observed / speaker_count # pbinom(q, lower.tail = FALSE) returns P(K > q). # Setting q = k_observed - 1 therefore returns P(K >= k_observed). p_value <- pbinom( k_observed - 1, size = speaker_count, prob = .5, lower.tail = FALSE ) sign_result <- binom.test( k_observed, speaker_count, p = .5, alternative = "greater" ) c( negative_count = k_observed, observed_proportion = observed_proportion, lower_95 = sign_result$conf.int[1], p_value = p_value ) ``` Report the sign count, the observed proportion, the confidence bound, and the p value. ::: {.callout-note title="Your result"} Replace this sentence with the requested sign-test result. ::: ## Exercise 5: Fisher's exact test (25 points) We'll use Fisher's exact test to analyze **German genitive allomorphy**, a case of morphological variation in which the genitive singular ending on masculine and neuter nouns alternates between *-s* and *-es*. For instance, the genitive of *Tag* “day” may be realized as *Tages* or *Tags*. The traditional pedagogical rule states that stems ending in sibilants require *-es*. Syllable count has also been proposed to matter, with monosyllabic nouns tending to favor *-es* and polysyllabic nouns tending to favor *-s*. The data come from the German [UniMorph](https://unimorph.github.io/) paradigm inventory described by [McCarthy et al. (2020)](https://aclanthology.org/2020.lrec-1.483/). The preparation scripts reduce the genitive records to one row per lemma. They record syllable class, whether a lemma ends in one of the strings `s`, `z`, `ß`, `sch`, or `x`, and whether a genitive form ending in *-es* is attested. These are variables produced by the scripts. UniMorph does not supply them as linguistic analyses. The raw paradigm may contain both *-s* and *-es* forms for one lemma. Thus `es_attested` records attestation in the extracted paradigm rather than an exclusive suffix choice. The teaching table contains a reproducible simple random sample of 120 lexemes. Its relevant columns are `lemma`, `syllable_class`, `ends_in_sibilant`, and `es_attested`. ::: {.callout-note title="Required background"} Review [the hypergeometric distribution](../../random-variables-and-distributions/hypergeometric-distribution.qmd) and [Fisher's exact test](../../statistical-inference/fishers-exact-test.qmd) before beginning the manual calculation below. These notes develop the 2×2 contingency table, fixed margins, and the probability-ordering definition of a two-sided p value. ::: ### 5.1 Load and inspect the lexeme table Load `data/unimorph_genitive_lexemes.csv`. Verify that `lemma` is unique. Examine the levels and marginal distributions of `syllable_class`, `ends_in_sibilant`, and `es_attested`. ```{r} genitive_lexemes <- read.csv("data/unimorph_genitive_lexemes.csv") str(genitive_lexemes) stopifnot(!anyDuplicated(genitive_lexemes$lemma)) count(genitive_lexemes, syllable_class) count(genitive_lexemes, ends_in_sibilant) count(genitive_lexemes, es_attested) ``` ### 5.2 Test syllable count and genitive *-es* attestation In German genitive allomorphy, is there an association between syllable class and the attestation of an *-es* form? Construct the two-by-two table for `syllable_class` by `es_attested`, with `monosyllabic` before `polysyllabic` and `TRUE` before `FALSE`. Plot the table with a mosaic plot or grouped bar chart. Run `fisher.test()` and calculate the sample odds ratio. With this order, the odds ratio places the odds of *-es* attestation among monosyllabic lexemes in the numerator and the corresponding odds among polysyllabic lexemes in the denominator. ```{r} syllable_table <- with( genitive_lexemes, table( factor( syllable_class, levels = c("monosyllabic", "polysyllabic") ), factor(es_attested, levels = c(TRUE, FALSE)) ) ) syllable_table mosaicplot(syllable_table, color = TRUE) syllable_fisher <- fisher.test(syllable_table) a_syllable <- syllable_table[1, 1] b_syllable <- syllable_table[1, 2] c_syllable <- syllable_table[2, 1] d_syllable <- syllable_table[2, 2] syllable_sample_or <- ___ syllable_fisher syllable_sample_or ``` Report the observed table, sample odds ratio, conditional estimate and interval returned by `fisher.test()`, and p value. ::: {.callout-note title="Your result"} Replace this sentence with the requested syllable-class result. ::: ### 5.3 Construct the sibilant table Is there an association between the string-based sibilant classification and the attestation of an *-es* form? Construct the table for `ends_in_sibilant` by `es_attested`, with the character value `"True"` before `"False"` in the rows and the logical value `TRUE` before `FALSE` in the columns. With this order, the sample odds ratio places the odds of *-es* attestation among sibilant-final lexemes in the numerator and the corresponding odds among other lexemes in the denominator. ```{r} sibilant_table <- with( genitive_lexemes, table( factor(ends_in_sibilant, levels = c("True", "False")), factor(es_attested, levels = c(TRUE, FALSE)) ) ) sibilant_table rowSums(sibilant_table) colSums(sibilant_table) ``` Report the four cells and both sets of margins. ::: {.callout-note title="Your result"} Replace this sentence with the observed sibilant table and its margins. ::: ### 5.4 Derive the conditional null probability [Fisher's exact test](../../statistical-inference/fishers-exact-test.qmd) tests independence in a two-by-two table while holding both observed margins fixed. Write the observed cells in row-major order as $a,b,c,d$. Let $A$ be the random upper-left cell count among all tables with those margins. Under the null hypothesis of independence, $A$ has a [hypergeometric distribution](../../random-variables-and-distributions/hypergeometric-distribution.qmd). Derive the probability of the observed table by counting the ways to select its first-column total from the two rows and dividing by the number of selections from all observations. #### Step 1: selections within each row The first row contains $a+b$ observations and contributes $a$ observations to the first column. The second row contains $c+d$ observations and contributes $c$. Write the product of the two counts. Complete one line. **Your LaTeX** $$ N_{\mathrm{observed}}=\text{enter the product of two binomial coefficients}. $$ #### Step 2: all selections with the fixed margins There are $a+b+c+d$ observations in total, and the first column contains $a+c$ of them. Write the number of possible first-column selections. Complete one line. $$ N_{\mathrm{all}}=\text{enter one binomial coefficient}. $$ #### Step 3: conditional probability Divide the count in Step 1 by the count in Step 2. Complete one line. $$ \mathbb P_0(A=a)=\text{enter }N_{\mathrm{observed}}/N_{\mathrm{all}}\text{ in terms of }a,b,c,d. $$ ### 5.5 Implement the sibilant test Calculate the sample odds ratio $ad/(bc)$ and the observed-table probability with `dhyper()`. Run `fisher.test()` to obtain the two-sided p value and conditional confidence interval. Fisher's two-sided p value sums the probabilities of all margin-compatible tables whose null probabilities are no greater than the observed table's probability. ```{r} a <- sibilant_table[1, 1] b <- sibilant_table[1, 2] c_cell <- sibilant_table[2, 1] d_cell <- sibilant_table[2, 2] sibilant_sample_or <- ___ # dhyper(x, m, n, k) uses the upper-left cell, the two row totals, # and the first-column total in that order. observed_table_probability <- dhyper( x = a, m = a + b, n = c_cell + d_cell, k = a + c_cell ) sibilant_fisher <- fisher.test(sibilant_table) tibble( predictor = c("syllable class", "sibilant ending"), sample_odds_ratio = c(syllable_sample_or, sibilant_sample_or) ) observed_table_probability sibilant_fisher ``` Report the odds-ratio orientation, the sample odds ratio, the conditional estimate returned by `fisher.test()`, its confidence interval, and the p value. Put the sibilant and syllable-class sample odds ratios in the same table. ::: {.callout-note title="Your result"} Replace this sentence with the requested sibilant result and odds-ratio comparison. ::: ### 5.6 Estimate sampling distributions across lexicon sizes A [sampling distribution](../../statistical-inference/properties-of-estimators.qmd) is the distribution of a statistic across repeated samples from a fixed population. For this exercise, treat `data/unimorph_genitive_lexeme_population.csv` as the finite population. For each $n\in\{12,20,40,80,120\}$, draw 500 simple random samples without replacement. Use one recorded seed. For every sample, record the odds ratio, Fisher p value, confidence-interval endpoints, whether any cell is zero, and whether either margin is zero. Use explicit factor levels when constructing every table. A **zero cell** is a table entry equal to zero. A **zero margin** occurs when an entire row or column is absent. `fisher.test()` should not be called when a margin is zero. Return `NA_real_` for its p value and interval in that case. A sample odds ratio may be zero or infinite when a cell is zero. Retain that value and exclude it only from summaries that require a finite logarithm. ```{r} genitive_population <- read.csv( "data/unimorph_genitive_lexeme_population.csv" ) analyze_lexeme_sample <- function(sample_size, replicate_id) { sampled_lexemes <- genitive_population |> slice_sample(n = sample_size, replace = FALSE) sampled_table <- with( sampled_lexemes, table( factor(ends_in_sibilant, levels = c("True", "False")), factor(es_attested, levels = c(TRUE, FALSE)) ) ) zero_cell <- any(sampled_table == 0) zero_margin <- any(rowSums(sampled_table) == 0) || any(colSums(sampled_table) == 0) a <- sampled_table[1, 1] b <- sampled_table[1, 2] c_cell <- sampled_table[2, 1] d_cell <- sampled_table[2, 2] # Preserve zero and infinite sample odds ratios. sample_or <- if (b * c_cell == 0) { if (a * d_cell > 0) Inf else NA_real_ } else { ___ } if (zero_margin) { fisher_p <- NA_real_ ci_lower <- NA_real_ ci_upper <- NA_real_ } else { sampled_fisher <- fisher.test(sampled_table) fisher_p <- sampled_fisher$p.value ci_lower <- sampled_fisher$conf.int[1] ci_upper <- sampled_fisher$conf.int[2] } tibble( sample_size = sample_size, replicate_id = replicate_id, zero_cell = zero_cell, zero_margin = zero_margin, odds_ratio = sample_or, fisher_p = fisher_p, ci_lower = ci_lower, ci_upper = ci_upper ) } set.seed(___) sampling_plan <- expand_grid( sample_size = c(12, 20, 40, 80, 120), replicate_id = seq_len(500) ) sampling_results <- bind_rows(Map( analyze_lexeme_sample, sampling_plan$sample_size, sampling_plan$replicate_id )) ``` For each $n$, calculate the proportion of samples with a zero cell, the proportion with a zero margin, the proportion with a nonfinite odds ratio, the median finite log odds ratio, and the median finite confidence-interval width on the log scale. Plot the sampling distributions of the finite log odds ratios by $n$. ```{r} sampling_summary <- sampling_results |> mutate( finite_log_or = is.finite(odds_ratio) & odds_ratio > 0, finite_log_interval = is.finite(ci_lower) & ci_lower > 0 & is.finite(ci_upper) & ci_upper > 0 ) |> group_by(sample_size) |> summarize( zero_cell_rate = mean(zero_cell), zero_margin_rate = mean(zero_margin), nonfinite_or_rate = mean(!finite_log_or), median_log_or = median(log(odds_ratio[finite_log_or])), median_log_interval_width = median( log(ci_upper[finite_log_interval]) - log(ci_lower[finite_log_interval]) ), .groups = "drop" ) sampling_summary sampling_results |> filter(is.finite(odds_ratio), odds_ratio > 0) |> ggplot(aes(x = factor(sample_size), y = log(odds_ratio))) + geom_boxplot() + labs(x = "Lexemes sampled", y = "Finite log sample odds ratio") ``` Report the summary table. ::: {.callout-note title="Your result"} Replace this sentence with the sample-size summary table. ::: ## Exercise 6: Chi-squared test of independence (25 points) We'll use the chi-squared test to analyze **German plural allomorphy**. German plural marking combines several suffix classes with possible stem changes. The classes represented in this teaching table include *-e* as in *Tag* → *Tage*, *-en* as in *Frau* → *Frauen*, *-er* as in *Kind* → *Kinder*, *-s* as in *Auto* → *Autos*, and forms for which the preparation rule finds none of these overt suffixes. The second UniMorph table has one row per noun lemma. The preparation rule retains one plural form and classifies its ending as `e`, `en`, `er`, `s`, or `ø`. Here `ø` means that the rule found none of the listed overt suffixes. It is a label produced by the script, not a claim that the plural has no exponence or that it contains umlaut. The linguistic question is whether grammatical gender and the derived plural-suffix class are independent in the lexeme inventory represented by the table. This creates a contingency table larger than 2×2. ::: {.callout-note title="Required background"} Review [the chi-squared test of independence](../../statistical-inference/chi-squared-test.qmd) and [the chi-squared distribution](../../random-variables-and-distributions/chi-squared-distribution.qmd) before beginning the derivations below. These notes develop the independence factorization, expected counts, Pearson's statistic, degrees of freedom, and the large-sample reference distribution. ::: Load `data/unimorph_plural_lexemes.csv`. The relevant columns are `lemma`, `gender`, and `plural_suffix`. Construct the observed three-by-five table with gender ordered `FEM`, `MASC`, `NEUT` and suffix ordered `e`, `en`, `er`, `s`, `ø`. ```{r} plural_lexemes <- read.csv("data/unimorph_plural_lexemes.csv") stopifnot(!anyDuplicated(plural_lexemes$lemma)) plural_table <- with( plural_lexemes, table( factor(gender, levels = c("FEM", "MASC", "NEUT")), factor(plural_suffix, levels = c("e", "en", "er", "s", "ø")) ) ) plural_table ``` ### 6.1 Derive the expected counts Under independence, the expected count in cell $(i,j)$ is determined by its row total, column total, and the grand total. Follow the derivation in the [chi-squared-test notes](../../statistical-inference/chi-squared-test.qmd). #### Step 1: independence probability Under the null hypothesis, the joint cell probability is the product of its row and column probabilities. Replace both probabilities with the corresponding observed marginal proportions. Complete one line. **Your LaTeX** $$ \widehat p_{ij} =\widehat p_{i+}\widehat p_{+j} =\text{enter the product of observed marginal proportions}. $$ #### Step 2: expected count Multiply the probability in Step 1 by the grand total $O_{++}$ and cancel one factor of $O_{++}$. Complete one line. $$ E_{ij}=\text{enter the expected-count formula}. $$ #### Step 3: feminine-by-*en* cell Read the feminine row total, the *-en* column total, and the grand total from `plural_table`. Substitute those three numbers. Complete one line. $$ E_{\mathrm{FEM},\mathrm{en}} =\text{enter the numerical substitution and result}. $$ ### 6.2 Derive the statistic and reference distribution Pearson's statistic sums the standardized squared difference between every observed and expected count. #### Step 1: statistic Sum the cell contribution $(O_{ij}-E_{ij})^2/E_{ij}$ over all rows and columns. Complete one line. **Your LaTeX** $$ \mathcal X^2=\text{enter Pearson's statistic}. $$ #### Step 2: degrees of freedom An $R\times C$ table has $RC$ cells. Fixed row and column margins leave $(R-1)(C-1)$ free cells. Substitute $R=3$ and $C=5$. Complete two lines. $$ \begin{aligned} \nu&=(R-1)(C-1),\\ \nu&=\text{enter the numerical degrees of freedom}. \end{aligned} $$ #### Step 3: upper-tail event Larger values of the statistic represent larger departures from independence. Complete one line. $$ B_{\mathcal X^2}=\{\omega\in\Omega\mid\text{enter the upper-tail condition}\}. $$ #### Step 4: reference probability Write the upper-tail probability using the cumulative distribution function $F_{\chi^2_\nu}$. Complete one line. $$ \mathbb P_0(B_{\mathcal X^2})=\text{enter the chi-squared tail probability}. $$ ### 6.3 Implement the test Compute the chi-squared statistic manually to verify R's output. Without calling `chisq.test()`, calculate the complete expected table, every cell contribution, $x^2_{\mathrm{obs}}$, the degrees of freedom, and the p value. Verify that every expected count is at least five. ```{r} row_totals <- rowSums(plural_table) column_totals <- colSums(plural_table) grand_total <- sum(plural_table) expected_counts <- outer(row_totals, column_totals) / ___ cell_contributions <- ___ x_squared_observed <- sum(cell_contributions) degrees_of_freedom <- ___ p_value <- pchisq( x_squared_observed, df = degrees_of_freedom, lower.tail = FALSE ) expected_counts cell_contributions c( x_squared = x_squared_observed, df = degrees_of_freedom, p_value = p_value, minimum_expected = min(expected_counts) ) ``` Verify the statistic, degrees of freedom, p value, and expected table with `chisq.test(..., correct = FALSE)`. Calculate the Pearson residuals $(O_{ij}-E_{ij})/\sqrt{E_{ij}}$ and identify the two cells with the largest absolute residuals. ```{r} software_result <- chisq.test(plural_table, correct = FALSE) software_result stopifnot(all(expected_counts >= 5)) pearson_residuals <- ___ residual_table <- as.data.frame(as.table(pearson_residuals)) |> mutate(absolute_residual = abs(Freq)) |> arrange(desc(absolute_residual)) slice_head(residual_table, n = 2) ``` Report $x^2_{\mathrm{obs}}$, the degrees of freedom, the p value, the minimum expected count, and the two largest absolute residuals with their signs. ::: {.callout-note title="Your result"} Replace this sentence with the requested chi-squared-test result. ::: ## Submission checklist - The notebook runs from a fresh R session after the data-preparation steps. - Every derivation is entered in its LaTeX block. - Every computation is entered in the R block immediately following its instruction. - Every requested result is entered in its result block. - The paired table is joined by identifiers rather than aligned by row position. - Every hypothesis names its parameter or independence relation. - Every p value is defined relative to its null reference distribution. - Every random calculation has a recorded seed.