Lab 2 — Instructor answer key

Companion file: solution.R (in this folder) — the lab chunk code with every TODO completed, in chunk order with the same seeds. Run it from the lab root (Rscript solutions/solution.R) and it prints every checkpoint, including the optional Cox and design-dial rows, and writes data/case_cohort_sample.csv. Figure chunks are omitted (they contain no TODOs); their numeric readouts are reproduced.

Contents: TODO code answers · expected checkpoints · model answers Q1–Q9 · understanding notes

TODO code answers

TODO	Answer
1a	`filter(!is.na(hospital_death))` — drops the 18 blank-status stays (all in the no-early-initiation group). This makes the lab a complete-case analysis, acknowledged in the handout as not ideal (missing-data methods are out of scope for this lab). Note that the drop happens before the subcohort draw, so the seed-13 checkpoints depend on it.
1b	`R1 <- a1 / (a1 + b1)`; `R0 <- a0 / (a0 + b0)`; `rr_cohort <- R1 / R0`; `rd_cohort <- R1 - R0`
1c	`early_vasopressor == 1` in the first line, `early_vasopressor == 0` in the second
2a	`samp <- coh \|> filter(subc \| hospital_death == 1)` — the union of the two design pieces: in the subcohort or a case. (`&` in place of `\|` keeps only the 29-person overlap — a giveaway `n_sample` of 29.)
3a	`stop = followup_days - eps` — a Type-2 (subcohort case) comparison row stops just before the event
3b	`start = pmax(0, followup_days - eps)` and `w = 1` — every case’s event enters exactly once, as the final sliver of their follow-up, at weight 1 (a Type-3 case contributes only this row)
5a	`slice_max(stop, n = 1, with_ties = FALSE)` — each person’s final person-time segment. This collapses the time-varying weights to one per person: 1 for events (Types 2–3), 1/q = 5 for sampled non-events (Type 1); equivalently, `ifelse(hospital_death == 1, 1, 1/q)`. (`slice_max(start, ...)` happens to select the same row; `slice_max(w, ...)` does not — it picks the weight-5 comparison row for the 29 subcohort cases, and `sum(w)` jumps from 1,291 to 1,407.)

Expected checkpoints

Verified against the canonical Lab 1 dataset (1,292 rows). All robust GLM standard errors use vcovHC(type = "HC3"), as in the lab chunks; the Cox robust SEs come from coxph(..., cluster = id).

setup:   nrow(lab1) = 1292
Step 1:  n_coh 1274 | 2x2: a1 6, b1 21, a0 130, b0 1117
         R1 0.222, R0 0.104, RR 2.132, RD 0.118
         person-days 168.3 (exposed), 7310.1 (unexposed), 7478.4 total
         rates 3.564 / 1.778 per 100 person-days; rate ratio 2.004
Step 2:  (q = 0.20, seed 13) subcohort 260; cases inside 29, outside 107;
         n_sample 367 distinct people
Step 3:  nrow(cc) 396 = 231 (Type 1) + 58 (Type 2: 29 people x 2 rows)
           + 107 (Type 3)
         weighted person-time 7664.6 vs cohort 7478.4 (+2.5%)
         weighted events exactly 136; time-zero weighted risk set 1300
Step 4:  KM risks at 5/10/30 days: full cohort .099/.154/.335;
           weighted case-cohort .096/.151/.354
         Cox log HR: full cohort 0.722 (SE 0.419); case-cohort 0.725,
           naive SE 0.419 (= the full-cohort SE), robust SE 0.655
Step 5:  exit weights: events 1, sampled non-events 5
         weighted 2x2: 6 / 20 (exposed), 130 / 1135 (unexposed);
           sum(w) = 1291
         weighted Poisson: log RR 0.8090 (RR 2.246), naive SE 0.4176
           (= the full-cohort model SE), robust SE 0.5858
         weighted linear: RD 0.1280, robust SE 0.1336 (robust only)
         full cohort: RR 2.132, robust SE 0.3830; RD 0.1180, robust 0.0835
         side-by-side tables:
           RR  full 2.132 [naive .418 | robust .383 | CI 1.006-4.516]
               cc   2.246 [naive .418 | robust .586 | CI 0.712-7.079]
           RD  full 0.118 [robust .084 | CI -0.046-0.282]
               cc   0.128 [robust .134 | CI -0.134-0.390]
Step 6:  bad person-time 8397.1 = 7664.6 + 732.5, and 732.5 is exactly
           the 136 cases' total follow-up; weighted events still 136
         KM .089/.139/.324 — about 0.92x the correct curve, everywhere
         rates per 100 pd: cohort 1.819 | cc 1.774 | double-counted 1.620
         Cox log HR: full 0.722/0.419 | cc 0.725/0.655 |
           double-counted 0.716/0.596
         simulation callout (see Step 6 of the handout): HR truth 2.00 |
           full 1.97 | correct split 2.04 | double-counted 1.79;
           person-time check overshoots by +33%
Step 7:  (seed 42, R = 2000) full_cohort | pooled | middle half:
         risk_ratio          2.132 | 2.125 | 1.727-2.791
         risk_difference     0.118 | 0.117 | 0.081-0.183
         rate_per100pd       1.819 | 1.819 | 1.720-1.935
         rate_double_counted 1.819 | 1.657 | 1.574-1.752
Step 8:  data/case_cohort_sample.csv, 367 rows x 15 columns
Dial (optional; seed 13 per q):
         q .10: n 261, cost $39,150, naive SE .418, robust SE .949 (RR 2.41)
         q .20: n 367, cost $55,050, naive SE .418, robust SE .586 (RR 2.25)
         q .40: n 584, cost $87,600, naive SE .418, robust SE .448 (RR 1.82)

If your Step-2 numbers differ, check that you kept set.seed(13), that the complete-case filter (1a) ran before the draw, and that you did not re-sort coh between the filter and the rbinom() call — the draw depends on row order.

Model answers, Q1–Q9

Q1. Why must the subcohort be drawn blind to the outcome, and why is the 29-person overlap expected? Drawn blind, the subcohort is a miniature of the baseline cohort, so scaled by 1/q it represents everyone’s denominators — person-time and risk denominators alike, for any outcome. Drawn with knowledge of the outcome it would represent nothing but itself. The overlap is the design working: a random 20% of a cohort containing 136 future cases should catch about q × 136 ≈ 27 of them (we got 29), and the analysis — not deletion — handles their dual role via the weight switch in Step 3. (Lecture’s lineage: the subcohort is exactly a case-base control series with the event times retained.)

Q2. What would be wrong with an unweighted risk of death computed on the 367 sampled people? Unweighted, the sample’s death risk is 136/367 = 0.371 against the cohort’s 0.107 — three and a half times too high, and not by chance: cases are 37% of the sample because we sampled all of them and only 20% of everyone else, so the proportion of cases is a design constant, not an estimate of anything. The 907 never-sampled non-cases (1,274 − 367) are exactly what the weights put back: until they do, every absolute quantity — risks, rates, the risk curve — is an artifact of the sampling rule rather than a feature of the cohort.

Q3. Why must a Type-2 person’s comparison row stop just before their event? If the comparison row ran to the event at weight 5, the Type-2 person would stand in their own event’s risk set twice — once as five people of comparison time, once as the weight-1 case — and, in lecture’s terms, the subcohort’s 1/q-scaled person-time already represents all of the cohort’s cases in expectation, so any additional appearance of sampled-case material double-counts. The weight must switch from 1/q to 1 at the event so that every case’s event is counted exactly once, at weight 1, like every other case’s.

Q4. Why must the naive SE of the weighted Cox fit reproduce the full-cohort SE, and why is that precision unearned? The weighted risk sets are scaled back up to full-cohort size (the time-zero weighted risk set is 1,300 ≈ 1,274), so the model-based information of the weighted partial likelihood is (essentially) the full cohort’s, and its inverse must reproduce the full-cohort SE — 0.419, to three decimals. The claim is unearned because only 367 exposures were measured: the model-based formula is blind to the fact that a different 20% draw would have delivered different risk sets and a different estimate. The robust, person-clustered SE (0.655) adds exactly that between-possible-subcohorts variability. The “missing information” sits in the 907 never-abstracted charts.

Q5. The methods-section sentence. “In a case-cohort sample of the 1,274-patient cohort — a 20% random subcohort (n = 260) plus all 136 in-hospital deaths, giving 367 patients with abstracted exposure — early vasopressor initiation was associated with a hazard ratio of 2.06 (95% CI 0.57–7.46) for in-hospital death, from a weighted Cox model with time-varying case-cohort weights (1/q over subcohort person-time, 1 at each event) and robust standard errors clustered on patient.” The essentials a reader needs: the subcohort fraction and size, all-cases ascertainment, the weighting scheme, and the variance estimator.

Q6. Which cells of the weighted 2×2 are exact, and which are estimates? The case cells (6 and 130) are exact: cases enter with probability 1 at weight 1, so the sample’s case cells are the cohort’s. The non-case cells are Horvitz–Thompson estimates: each sampled non-case stands for 5, so 4 exposed non-cases become 20 (truth 21) and 227 unexposed become 1,135 (truth 1,117). One piece of the design per cell type: “all the cases” delivers the numerators exactly; the subcohort estimates the denominators.

Q7. Whose time was double-counted in Step 6, and by which two routes? The 136 cases’ pre-event follow-up (732.5 person-days total) was double-counted. Route one: the subcohort’s comparison rows, scaled by 1/q, already represent the whole cohort’s person-time — including the follow-up of people who go on to become cases. Route two: the double-counted case rows add that same follow-up again, at weight 1, from time zero. The event check passes because the event rows are untouched: one weight-1 event per case, 136 in total.

Q8. The two-sentence rebuttal to “the ratios were basically right, so the shortcut is fine.” “The ratios survived because both arms’ denominators were inflated by nearly the same fraction — exposed person-time by 8.3%, unexposed by 9.8% — a numerical coincidence of how case follow-up happens to distribute across arms in this cohort, not a property of the error; the simulation aside shows the identical mistake cutting a true hazard ratio of 2.00 to 1.79 when the outcome is common. Meanwhile every absolute quantity — risks, rates, the entire risk curve — is about 10% wrong everywhere, and delivering valid absolute risks is precisely what the case-cohort design promises over case-control; a design whose selling point is broken is not vindicated by ratios that got lucky.”

Q9. What drives the single-draw spread, why can’t we oversample the exposed, and what could we legitimately dial? The driver is k1, the number of exposed non-cases the subcohort catches: k1 ~ Binomial(21, 0.2), expectation ≈ 4. It is so influential because the weighted exposed denominator is 6 + 5·k1 — a handful of people carrying a fifth of the estimate — a direct consequence of this cohort’s 2.1% exposure prevalence. Oversampling exposed people into the subcohort cannot be implemented as stated because membership would then depend on the exposure, the very variable the design exists to avoid measuring on everyone. Legitimate dials: a larger q (the dial table), a larger cohort, or stratifying the subcohort on a cheap, already-recorded correlate of exposure with stratum-specific sampling fractions carried into the weights (lecture’s stratified extension).

Understanding notes

The key questions are Q3 and Q7 (the same double-counting mechanism, approached from two directions — connect them), Q4 (where the missing information physically sits), and Q2/Q6 (which numbers are design constants, which are estimates). In every case what matters is the mechanism — who is represented by what, at which weight — not the vocabulary.
The lab’s core distinction, worth stating once cleanly: time-to-event analyses (KM, Cox) are built on time-indexed risk sets and need the person-time split, because a sampled case’s weight changes at their event (1/q while serving as comparison time, 1 at the event). The binary-outcome regressions (log-Poisson risk ratio, linear risk difference) have no time axis, so each person carries one weight — their final segment’s. The split is not about time-varying exposures; the exposure is fixed at baseline. It exists because the sampling probability — hence the weight — changes over time.
The two ways this analysis dies quietly, and their catches: a model-based standard error that reports full-cohort precision the abstraction budget never bought (caught by the robust, person- clustered variance — note the naive SEs 0.419/0.4176 equal the full-cohort model SEs to three or four decimals, in both the Cox and the Poisson); and outcome-selected case person-time smuggled into the denominators (caught by the person-time design check, which overshoots by exactly the cases’ total follow-up, 732.5 person-days).
One nuance worth not over-teaching: “naive is always too small” is false as a slogan — the full-cohort Poisson’s model-based SE (0.4176) is larger than its robust SE (0.3830), because Poisson variance is conservative for a binary outcome. The correct statement is that the weighted fit’s model-based SE reproduces the full-cohort model-based SE, i.e., it prices 1,274 abstracted charts when only 367 were bought.