Spread the stress scores that piled up at 40

scripts/generate-datasets.mjs

@@ -35,7 +35,10 @@ function writeWellbeingPopulation() {
   for (let i = 1; i <= N; i += 1) {
     const programme = rng() < 0.5 ? 'Psychology' : 'Business';
     // Stress is skewed: most students low, a long tail of very stressed ones.
-    const stress = Math.min(40, -12 * Math.log(Math.max(rng(), Number.EPSILON)));
+    // The scale runs 0 to 40, like the Perceived Stress Scale. Wrapping a draw
+    // past 40 back round to the start keeps the exponential shape inside the
+    // scale (it is memoryless), where clipping piled 190 students onto 40.
+    const stress = (-12 * Math.log(Math.max(rng(), Number.EPSILON))) % 40;
     const sleep = Math.min(11, Math.max(3, normal(rng, 7.2 - stress * 0.03, 1.1)));
     const exam = Math.min(100, Math.max(0, normal(rng, 62 + sleep * 2.4 - stress * 0.45, 11)));
     rows.push([i, programme, stress.toFixed(2), sleep.toFixed(2), exam.toFixed(1)].join(','));

src/content/exercises/module-05.ts

@@ -38,10 +38,10 @@ export const module05: ExerciseDef[] = [
           population <- read.csv("data/wellbeing-population.csv")
           mu <- mean(population$exam_score)
           sigma <- sd(population$exam_score)
-          expected <- pnorm(60, mean = mu, sd = sigma)   # 0.140774
+          expected <- pnorm(60, mean = mu, sd = sigma)   # 0.125920
           # tolerance 1e-3, not 1e-6: a student who computes the SD with denominator N
-          # rather than sd()'s n - 1 lands on 0.140750 - the same answer to every digit
-          # anyone would report. The three wrong answers are 0.859, 0.018 and 1.000,
+          # rather than sd()'s n - 1 lands on 0.125896 - the same answer to every digit
+          # anyone would report. The three wrong answers are 0.874, 0.017 and 1.000,
           # all far outside this band.
           if (isTRUE(all.equal(value, expected, tolerance = 1e-3, check.attributes = FALSE))) {
             list(pass = TRUE, message = paste0("Correct: ", round(expected, 4), ". Under this normal model about ", round(100 * expected), " out of every 100 students score below 60, and the actual count in the population is ", sum(population$exam_score < 60), " out of 5000 - the model is close, not exact."))
@@ -99,10 +99,10 @@ export const module05: ExerciseDef[] = [
           population <- read.csv("data/wellbeing-population.csv")
           mu <- mean(population$exam_score)
           sigma <- sd(population$exam_score)
-          expected <- (85 - mu) / sigma    # 0.946517
+          expected <- (85 - mu) / sigma    # 0.916170
           # tolerance 1e-3 for the same reason as m5-1-a: an SD computed with
-          # denominator N gives 0.946612, which rounds identically at two decimals.
-          # The wrong answers are 11.69, 6.88 and -0.95.
+          # denominator N gives 0.916262, which rounds identically at two decimals.
+          # The wrong answers are 11.11, 7.01 and -0.92.
           if (isTRUE(all.equal(value, expected, tolerance = 1e-3, check.attributes = FALSE))) {
             list(pass = TRUE, message = paste0("Correct: z = ", round(expected, 2), ". That student is just under one standard deviation above the mean, which puts them above about ", round(100 * pnorm(expected)), " per cent of students."))
           } else if (isTRUE(all.equal(value, -expected, tolerance = 1e-3, check.attributes = FALSE))) {
@@ -155,8 +155,8 @@ export const module05: ExerciseDef[] = [
           list(pass = FALSE, message = "Both z_exam and z_sleep should be single numbers.")
         } else {
           population <- read.csv("data/wellbeing-population.csv")
-          want_exam <- (85 - mean(population$exam_score)) / sd(population$exam_score)      # 0.946517
-          want_sleep <- (9.1 - mean(population$sleep_hours)) / sd(population$sleep_hours)  # 1.996967
+          want_exam <- (85 - mean(population$exam_score)) / sd(population$exam_score)      # 0.916170
+          want_sleep <- (9.1 - mean(population$sleep_hours)) / sd(population$sleep_hours)  # 1.986021
           # tolerance 1e-3 throughout, as in m5-2-a: an SD computed with denominator N
           # is within 1e-4 of sd(), and no wrong answer here is closer than 0.25.
           ok_exam <- isTRUE(all.equal(z_exam, want_exam, tolerance = 1e-3, check.attributes = FALSE))
@@ -219,9 +219,9 @@ export const module05: ExerciseDef[] = [
           population <- read.csv("data/wellbeing-population.csv")
           mu <- mean(population$exam_score)
           sigma <- sd(population$exam_score)
-          expected <- qnorm(0.90, mean = mu, sd = sigma)   # 89.13958
-          # tolerance 1e-4: an SD computed with denominator N gives 89.13800, a
-          # relative difference of 1.8e-5. The wrong answers are 57.5, 1.28 and
+          expected <- qnorm(0.90, mean = mu, sd = sigma)   # 89.42980
+          # tolerance 1e-4: an SD computed with denominator N gives 89.42825, a
+          # relative difference of 1.7e-5. The wrong answers are 58.4, 1.28 and
           # essentially 0, so nothing plausible sneaks through this band.
           if (isTRUE(all.equal(value, expected, tolerance = 1e-4, check.attributes = FALSE))) {
             list(pass = TRUE, message = paste0("Correct: ", round(expected, 2), ". In the actual population ", sum(population$exam_score > expected), " of 5000 students score above it - close to the 500 the model predicts."))
@@ -281,8 +281,8 @@ export const module05: ExerciseDef[] = [
           sigma <- sd(population$exam_score)
           lower <- pnorm(65, mean = mu, sd = sigma)
           upper <- pnorm(85, mean = mu, sd = sigma)
-          expected <- upper - lower   # 0.577321
-          # tolerance 1e-3, as in m5-1-a: the denominator-N SD gives 0.577290.
+          expected <- upper - lower   # 0.588578
+          # tolerance 1e-3, as in m5-1-a: the denominator-N SD gives 0.588624.
           if (isTRUE(all.equal(value, expected, tolerance = 1e-3, check.attributes = FALSE))) {
             list(pass = TRUE, message = paste0("Correct: ", round(expected, 4), ". Just under three students in five fall in that twenty-point band, and the actual population count is ", sum(population$exam_score > 65 & population$exam_score < 85), " of 5000."))
           } else if (value > 1) {

src/content/exercises/module-06.ts

@@ -187,7 +187,7 @@ export const module06: ExerciseDef[] = [
           expected <- sigma / sqrt(40)
           # 1e-3 rather than 1e-6: admits a population SD computed with denominator N
           # instead of sd()'s n - 1, which differs by ~1e-4 here. Both wrong answers
-          # below (sigma ~10.6, sigma/40 ~0.27) are far outside this band.
+          # below (sigma ~9.2, sigma/40 ~0.23) are far outside this band.
           # check.attributes = FALSE on each comparison: a second guard, so names never pick the wrong branch.
           if (isTRUE(all.equal(value, expected, tolerance = 1e-3, check.attributes = FALSE))) {
             list(pass = TRUE, message = paste0("Correct: ", round(expected, 3), ". Individual students vary by about ", round(sigma, 2), ", but sample means of 40 vary by only ", round(expected, 3), "."))

src/content/exercises/module-07.ts

@@ -122,7 +122,7 @@ export const module07: ExerciseDef[] = [
           # sort(): an interval written upper-first is the same interval.
           value <- sort(value)
           if (isTRUE(all.equal(value, expected, tolerance = 1e-6, check.attributes = FALSE))) {
-            list(pass = TRUE, message = paste0("Correct: [", round(expected[1], 2), ", ", round(expected[2], 2), "]. The population mean is 73.31, so this interval does contain it - which about 19 intervals in 20 will."))
+            list(pass = TRUE, message = paste0("Correct: [", round(expected[1], 2), ", ", round(expected[2], 2), "]. The population mean is 73.89, so this interval does contain it - which about 19 intervals in 20 will."))
           } else if (isTRUE(all.equal(value, sort(m + c(-1, 1) * 1.96 * se), tolerance = 1e-6, check.attributes = FALSE))) {
             list(pass = FALSE, message = "You used 1.96, which comes from the normal distribution and assumes you KNOW the SD. You estimated it from 50 students, so the multiplier is qt(0.975, df = 49) = 2.01 - a slightly wider interval, and that extra width is the price of not knowing sigma.")
           } else if (isTRUE(all.equal(value, sort(m + c(-1, 1) * se), tolerance = 1e-6, check.attributes = FALSE))) {

src/content/exercises/module-08.ts

@@ -271,17 +271,18 @@ export const module08: ExerciseDef[] = [
         } else {
           # A band, not a value: power here is estimated by simulation, and every
           # legitimate route draws its own random numbers. The true value for
-          # n = 100 against mu = 70 is about .75, and 1000 replicates give a
-          # standard error of about .013, so .67 to .85 is about six standard errors
-          # either way. The wrong answers land at about .24 (the complement),
-          # .05 (alpha), .06 (the mean p-value) and .12 (n = 10).
-          if (value >= 0.67 && value <= 0.85) {
+          # n = 100 against mu = 70 is about .89, and 1000 replicates give a
+          # standard error of about .010, so .83 to .95 is about six standard errors
+          # either way. The wrong answers land at about .11 (the complement),
+          # .14 (n = 10), .03 (alpha) and .02 (the mean p-value). The first two are
+          # too close to tell apart by value, so they share a message.
+          if (value >= 0.83 && value <= 0.95) {
             list(pass = TRUE, message = paste0("About ", round(100 * value), " studies in 100 would detect this effect. That is the power of the design - and it means roughly ", round(100 * (1 - value)), " in 100 would miss it and report a non-significant result, even though the effect is genuinely there."))
-          } else if (value >= 0.15 && value <= 0.33) {
-            list(pass = FALSE, message = "That is the proportion of studies that FAILED to reject - the Type II error rate. Power is its complement: count p_values below alpha, not above it.")
-          } else if (value <= 0.10) {
+          } else if (value >= 0.09 && value <= 0.20) {
+            list(pass = FALSE, message = "That is far below this design's power. Two mistakes land here. You may have counted the studies that FAILED to reject, which is the Type II error rate: power is its complement, so count p_values below alpha, not above it. Or each study has far fewer than the 100 students the exercise specifies.")
+          } else if (value < 0.09) {
             list(pass = FALSE, message = "That is close to alpha, the rejection rate you get when the null is TRUE. Power is the rejection rate when the null is FALSE, so the test has to be against mu = 70 - a value the population mean is genuinely not equal to. (If you averaged the p-values instead of counting rejections, that also lands here.)")
-          } else if (value < 0.67) {
+          } else if (value < 0.83) {
             list(pass = FALSE, message = paste0("A power of ", round(value, 2), " is too low for this design. Check the sample size inside replicate(): the exercise specifies 100 students per study, and power climbs steeply with n."))
           } else {
             list(pass = FALSE, message = paste0("A power of ", round(value, 2), " is higher than this design can deliver. Check that alpha is 0.05 and that each study samples 100 students."))
@@ -340,7 +341,7 @@ export const module08: ExerciseDef[] = [
         } else if (abs(abs(got) - abs(target)) > 1e-6) {
           list(pass = FALSE, message = paste0("d is ", round(got, 3), ", but it should be ", round(target, 3), ". Check the pooled SD: each group's variance times its n - 1, added, divided by n1 + n2 - 2, then the square root."))
         } else {
-          list(pass = TRUE, message = paste0("Correct: d = ", format(round(abs(target), 2), nsmall = 2), ". The whole population's d is 0.50, so this significant small study overestimates the effect, as significant small studies tend to."))
+          list(pass = TRUE, message = paste0("Correct: d = ", format(round(abs(target), 2), nsmall = 2), ". The whole population's d is 0.47, so this significant small study overestimates the effect, as significant small studies tend to."))
         }
       }
     `,

src/content/lessons/05-1-density-and-area.mdx

@@ -12,9 +12,9 @@ population <- read.csv("data/wellbeing-population.csv", stringsAsFactors = TRUE)
 
 population %>% summarise(mu = mean(exam_score), sigma = sd(exam_score), n = n())`} />
 
-Two numbers: a mean of about 73.3 and a standard deviation of about 12.4. (sd()
-divides by n − 1; the textbook's population σ divides by N and gives 12.35
-instead of 12.36, a difference that changes nothing in this module.) The
+Two numbers: a mean of about 73.9 and a standard deviation of about 12.1. (sd()
+divides by n − 1; the textbook's population σ divides by N and gives 12.122
+instead of 12.124, a difference that changes nothing in this module.) The
 claim of this lesson is that those two numbers are almost the whole story.
 
 <CodeBlock id="hist" code={`library(ggplot2)
@@ -53,11 +53,11 @@ rest of this module treats it as if it were the truth.
 
 <Predict
   id="p-height"
-  question="dnorm(73.3, mean = 73.3, sd = 12.4) returns 0.0322. What does that number tell you?"
+  question="dnorm(73.9, mean = 73.9, sd = 12.1) returns 0.0330. What does that number tell you?"
   choices={[
-    { text: 'About 3.2% of students score exactly 73.3', response: 'No student scores *exactly* 73.3 in any meaningful sense, and a density is not a percentage. Read on — this is the confusion the lesson exists to clear up.' },
-    { text: 'The height of the curve at 73.3, which only becomes a probability once you multiply it by a width', correct: true, response: 'Exactly. A density is a height. Probability is height times width, which is to say area.' },
-    { text: 'About 3.2% of students score below 73.3', response: 'That is what pnorm() would tell you, and the answer would be close to 50%, not 3.2%.' },
+    { text: 'About 3.3% of students score exactly 73.9', response: 'No student scores *exactly* 73.9 in any meaningful sense, and a density is not a percentage. Read on — this is the confusion the lesson exists to clear up.' },
+    { text: 'The height of the curve at 73.9, which only becomes a probability once you multiply it by a width', correct: true, response: 'Exactly. A density is a height. Probability is height times width, which is to say area.' },
+    { text: 'About 3.3% of students score below 73.9', response: 'That is what pnorm() would tell you, and the answer would be close to 50%, not 3.3%.' },
   ]}
 />
 
@@ -80,13 +80,13 @@ pnorm(85, mean = mu, sd = sigma) - pnorm(65, mean = mu, sd = sigma)
 # And the model can be checked against the population it describes.
 mean(population$exam_score < 60)`} />
 
-Compare the first line with the last. The model says about 14 per cent of students score
-below 60; the population itself says about 15 per cent. Close, not identical —
+Compare the first line with the last. The model says about 12.6 per cent of students score
+below 60; the population itself says 13.3 per cent. Close, not identical —
 which is what "the data are approximately normal" means in practice.
 
 Notice also what the area of a single point must be. A point has no width, so
-`pnorm(73.3, mu, sigma) - pnorm(73.3, mu, sigma)` is exactly zero: the probability of scoring
-*exactly* 73.3 is zero under a continuous model. That is not a paradox. It is
+`pnorm(73.9, mu, sigma) - pnorm(73.9, mu, sigma)` is exactly zero: the probability of scoring
+*exactly* 73.9 is zero under a continuous model. That is not a paradox. It is
 the price of modelling a score as a real number rather than as a count.
 
 ## See it move
@@ -114,9 +114,9 @@ at 1 however wide the curve spreads.
   id="i-5-1"
   question="You are writing up the distribution of exam scores for a methods report. Which sentence reports what you have actually computed?"
   choices={[
-    { text: 'Exam scores were approximately normally distributed (μ = 73.31, σ = 12.36). Under this model, 14% of students are expected to score below 60.', correct: true, response: 'Correct. It names the model, gives the two parameters in APA form, and reports the tail area as an expected proportion rather than as a claim about any one student.' },
-    { text: 'The density of exam scores at 60 was .018, so 1.8% of students scored 60.', response: 'That is dnorm(60), a height. Heights are not probabilities, and no student scores exactly 60 under a continuous model.' },
-    { text: 'The probability that a student scored exactly 60 was .14.', response: '.14 is the area to the LEFT of 60, covering every score below it. The probability of any exact score under a continuous model is zero.' },
-    { text: 'Exam scores were normally distributed (μ = 73.31, σ = 12.36), so 14% of students scored below 60.', response: 'Two problems. "Were normally distributed" overstates it (the fit is close, not exact), and the 14% is what the model expects, while the population itself holds about 15%.' },
+    { text: 'Exam scores were approximately normally distributed (μ = 73.89, σ = 12.12). Under this model, 13% of students are expected to score below 60.', correct: true, response: 'Correct. It names the model, gives the two parameters in APA form, and reports the tail area as an expected proportion rather than as a claim about any one student.' },
+    { text: 'The density of exam scores at 60 was .017, so 1.7% of students scored 60.', response: 'That is dnorm(60), a height. Heights are not probabilities, and no student scores exactly 60 under a continuous model.' },
+    { text: 'The probability that a student scored exactly 60 was .13.', response: '.13 is the area to the LEFT of 60, covering every score below it. The probability of any exact score under a continuous model is zero.' },
+    { text: 'Exam scores were normally distributed (μ = 73.89, σ = 12.12), so 13% of students scored below 60.', response: 'Two problems. "Were normally distributed" overstates it (the fit is close, not exact), and the 13% is what the model expects (12.6%), while the population itself holds 13.3%.' },
   ]}
 />