Chapter 5 — Main study and independent validation

Chapter 5 — Main study and independent validation

Lesson map

Lesson 05 published. Map, example, EX028.

flowchart TD
  N0["Learning outcomes and connection"]
  N1["5.1 Designing the simulation experiment"]
  N2["5.2 The confirmatory simulation study"]
  N3["5.3 The operator tabletop study"]
  N4["5.4 Validating the model"]
  N5["5.5 Independent validation in a second context"]
  N6["5.6 Reporting the main study"]
  N7["Chapter synthesis"]
  N0 --> N1 --> N2 --> N3 --> N4 --> N5 --> N6 --> N7
  N7 --> EX["Worked example"]
  EX --> Q["Exercise EX028"]

Logical drawing

flowchart LR
  N0["Learning outcomes and connection"]
  N1["5.1 Designing the simulation experiment"]
  N2["5.2 The confirmatory simulation study"]
  N3["5.3 The operator tabletop study"]
  N4["5.4 Validating the model"]
  N5["5.5 Independent validation in a second context"]
  N6["5.6 Reporting the main study"]
  N7["Chapter synthesis"]
  N0 --> N1 --> N2 --> N3 --> N4 --> N5 --> N6 --> N7 --> Q["EX028"]

[Learning outcomes and connection]
    |
    v
[5.1 Designing the simulation experiment]
    |
    v
[5.2 The confirmatory simulation study]
    |
    v
[5.3 The operator tabletop study]
    |
    v
[5.4 Validating the model]
    |
    v
[5.5 Independent validation in a second context]
    |
    v
[5.6 Reporting the main study]
    |
    v
[Chapter synthesis]
    |
    v
[Worked example] --> [Exercise EX028]

Learning outcomes and connection

By the end of this chapter the learner will be able to:

1. Design a paired simulation experiment with common random numbers, a justified number of draws and an oracle bound.
2. Distinguish Monte Carlo uncertainty from parameter uncertainty and structural uncertainty, and say which one a reported interval describes.
3. Interpret a mixed pattern of supported and unsupported hypotheses without rescuing the favorable ones.
4. Design and analyse a blinded operator tabletop study whose outcome is the consequence of the choices made, not stated preference.
5. Validate a model using face, historical, extreme-condition and operational validity evidence, and state what each establishes.
6. Run an independent validation in a second context selected on predeclared criteria and reason about transportability.

Chapter 4 verified the simulator, reproduced the baselines and piloted the procedures. This chapter runs the confirmatory study exactly as preregistered and then tests whether its result travels.

5.1 Designing the simulation experiment

5.1.1 Paired comparisons and common random numbers

Each simulation draw samples one scenario: restoration durations for every feeder and buffers for every service. Every policy is evaluated on the same draw. This common random numbers design removes shared between-scenario variance from the comparison: the quantity of interest is the per-draw difference between two policies, whose variance is Var(A) + Var(B) − 2 Cov(A, B). Pairing helps exactly as much as the two policies' losses are positively correlated across draws. When two policies trade off — one protects water on a draw where the other protects telecom — the correlation is weak and pairing helps little; EX028 shows this happening in Riverbend. Analyse paired differences anyway, because they are the correct estimand, but do not assume that pairing guarantees precision.

5.1.2 How many draws?

The Monte Carlo standard error of a mean difference is s_d / √n, where s_d is the standard deviation of the per-draw differences. In the pilot (Chapter 4) the standard deviation of DEP − CRITFLAG differences in Riverbend was about 7 weighted service-hours, so 500 draws give a 95% half-width of about 1.96 × 7 / √500 ≈ 0.6. The smallest effect of interest for H2 is 10% of CRITFLAG's loss; because CRITFLAG's loss in Riverbend is small (about 4), that threshold is about 0.4 — smaller than the half-width. The preregistration therefore fixed a three-way decision rule: H2 is supported if the interval lies wholly below −10% of CRITFLAG's mean loss, rejected if it lies wholly above that threshold (a reduction of at least 10% is ruled out), and inconclusive otherwise, in which case the planned remedy was 2,000 additional draws without the oracle. Record the calculation and the rule; do not increase n after seeing a "nearly significant" result outside the rule you wrote.

5.1.3 What a Monte Carlo interval means

A bootstrap or normal interval over simulation draws describes Monte Carlo uncertainty — how much the estimate would move with more draws from this model with these parameter distributions. It says nothing about parameter uncertainty (whether the elicited distributions are right) or structural uncertainty (whether the model's form is right). With enough draws, any difference becomes "statistically significant" inside a model. Report Monte Carlo intervals for completeness, but base claims on effect sizes relative to the smallest effect of interest, and address parameter and structural uncertainty with sensitivity analyses and validation (§5.4 and Chapter 6).

5.2 The confirmatory simulation study

5.2.1 The simulator

The event-based simulator below implements the preregistered design. A service is unavailable while its feeder is de-energized beyond its buffer, or while a service it depends on is unavailable beyond its on-site buffer; a one-hour restart lag follows restoration. Crews take the next feeder in the policy's order as soon as they are free. Customer-hours are computed per feeder. ORACLE enumerates every ordering for each draw. Save the code as lab5_mc.py; Chapter 6 imports it.

# Lab 5 - paired Monte Carlo comparison of restoration orderings (offline, standard library only)
# Event-based: a service is unavailable while its feeder is down beyond its buffer, or while a
# service it depends on is unavailable beyond its on-site buffer; restart lag applies after restoration.
import itertools, random, statistics

RIVERBEND = {
    "customers": {"F1": 9000, "F2": 1200, "F3": 2500, "F4": 6000},
    "crews": 1, "dur": (2, 6, 4),          # triangular(low, high, mode) hours per feeder
    "services": {  # name: (feeder, buffer triangular, weight, (depends_on, on-site hours) or None)
        "WATER": ("F3", (4, 14, 8), 4.0, None),
        "TEL":   ("F2", (4, 12, 8), 3.0, None),
        "HOSP":  ("F2", (8, 48, 24), 5.0, ("WATER", 4.0)),
    },
}
PORT_HALVARD = {
    "customers": {"H1": 7000, "H2": 3000, "H3": 1500, "H4": 4000, "H5": 800, "H6": 5000},
    "crews": 2, "dur": (1.5, 5, 3),
    "services": {
        "WATER2":   ("H4", (10, 30, 20), 4.0, None),
        "HOSP2":    ("H2", (48, 96, 72), 5.0, ("WATER2", 6.0)),
        "TEL2":     ("H3", (2, 6, 4), 3.0, None),
        "DISPATCH": ("H5", (1, 3, 2), 4.0, None),
    },
}
H, LAG = 72.0, 1.0

def draw(region, rng, buffer_scale=1.0):
    dur = {f: rng.triangular(*region["dur"]) for f in sorted(region["customers"])}
    buf = {s: rng.triangular(*v[1]) * buffer_scale for s, v in sorted(region["services"].items())}
    return dur, buf

def restore_times(order, dur, crews):
    free, done = [0.0] * crews, {}
    for f in order:
        i = min(range(crews), key=lambda k: free[k])
        free[i] += dur[f]
        done[f] = free[i]
    return done

def outage(region, done, buf):
    """Return {service: (start, end)} of unavailability, clipped to [0, H]."""
    out = {}
    def solve(s):
        if s in out:
            return out[s]
        feeder, _, _, dep = region["services"][s]
        ivs = []
        if buf[s] < done[feeder]:
            ivs.append((buf[s], done[feeder] + LAG))
        if dep:
            d = solve(dep[0])
            if d and d[0] + dep[1] < d[1]:
                ivs.append((d[0] + dep[1], d[1]))
        out[s] = (min(a for a, _ in ivs), max(b for _, b in ivs)) if ivs else None
        return out[s]
    for s in region["services"]:
        solve(s)
    return {s: (max(0.0, iv[0]), min(H, iv[1])) for s, iv in out.items() if iv}

def evaluate(region, order, dur, buf):
    done = restore_times(order, dur, region["crews"])
    loss = sum(region["services"][s][2] * (e - b) for s, (b, e) in outage(region, done, buf).items() if e > b)
    cust_h = sum(c * min(H, done[f]) for f, c in region["customers"].items()) / 1000.0
    return loss, cust_h

def silo(region):
    return sorted(region["customers"], key=lambda f: -region["customers"][f])

def critflag(region, believed_feeder=None):
    feeder_of = {s: v[0] for s, v in region["services"].items()}
    feeder_of.update(believed_feeder or {})   # the critical-customer list may also be out of date
    crit = set(feeder_of.values())
    return sorted(region["customers"], key=lambda f: (f not in crit, -region["customers"][f]))

def dep(region, believed_feeder=None):
    svc = dict(region["services"])
    if believed_feeder:                       # planner's (possibly wrong) dependency map
        for s, f in believed_feeder.items():
            svc[s] = (f,) + svc[s][1:]
    mean = {s: sum(v[1]) / 3 for s, v in svc.items()}
    tti = {}
    def t(s):
        if s not in tti:
            dep_ = svc[s][3]
            tti[s] = min(mean[s], t(dep_[0]) + dep_[1]) if dep_ else mean[s]
        return tti[s]
    key = {f: min([t(s) for s, v in svc.items() if v[0] == f], default=float("inf")) for f in region["customers"]}
    return sorted(region["customers"], key=lambda f: (key[f], -region["customers"][f]))

def experiment(region, n=500, seed=7, buffer_scale=1.0, believed_feeder=None, oracle=True):
    rng = random.Random(seed)
    policies = {"SILO": silo(region), "CRITFLAG": critflag(region, believed_feeder), "DEP": dep(region, believed_feeder)}
    res = {k: [] for k in list(policies) + (["ORACLE"] if oracle else [])}
    cust = {k: [] for k in policies}
    perms = list(itertools.permutations(region["customers"]))
    for _ in range(n):
        dur, buf = draw(region, rng, buffer_scale)
        for k, order in policies.items():
            l, c = evaluate(region, order, dur, buf)
            res[k].append(l)
            cust[k].append(c)
        if oracle:
            res["ORACLE"].append(min(evaluate(region, list(p), dur, buf)[0] for p in perms))
    return policies, res, cust

def boot_ci(x, seed=11, b=1000):
    rng = random.Random(seed)
    means = sorted(statistics.fmean(rng.choices(x, k=len(x))) for _ in range(b))
    return round(means[int(0.025 * b)], 2), round(means[int(0.975 * b)], 2)

if __name__ == "__main__":
    for name, region in (("Riverbend", RIVERBEND), ("Port Halvard", PORT_HALVARD)):
        policies, res, cust = experiment(region)
        print(f"== {name}: orders {policies}")
        for k, v in res.items():
            extra = f"  customer-hours(k) {statistics.fmean(cust[k]):7.1f}" if k in cust else ""
            print(f"  {k:8} weighted critical loss {statistics.fmean(v):6.2f}{extra}")
        d = [a - b for a, b in zip(res["DEP"], res["CRITFLAG"])]
        print("  DEP - CRITFLAG loss: mean", round(statistics.fmean(d), 2), "95% CI", boot_ci(d))
        r = statistics.fmean(cust["DEP"]) / statistics.fmean(cust["SILO"]) - 1
        r2 = statistics.fmean(cust["DEP"]) / statistics.fmean(cust["CRITFLAG"]) - 1
        print("  DEP customer-hours vs SILO:", f"{r:+.1%}", "| vs CRITFLAG:", f"{r2:+.1%}")

Expected output:

== Riverbend: orders {'SILO': ['F1', 'F4', 'F3', 'F2'], 'CRITFLAG': ['F3', 'F2', 'F1', 'F4'], 'DEP': ['F2', 'F3', 'F1', 'F4']}
  SILO     weighted critical loss  52.75  customer-hours(k)   133.0
  CRITFLAG weighted critical loss   3.97  customer-hours(k)   222.5
  DEP      weighted critical loss   4.32  customer-hours(k)   227.6
  ORACLE   weighted critical loss   1.79
  DEP - CRITFLAG loss: mean 0.35 95% CI (-0.23, 0.97)
  DEP customer-hours vs SILO: +71.2% | vs CRITFLAG: +2.3%
== Port Halvard: orders {'SILO': ['H1', 'H6', 'H4', 'H2', 'H3', 'H5'], 'CRITFLAG': ['H4', 'H2', 'H3', 'H5', 'H1', 'H6'], 'DEP': ['H5', 'H3', 'H4', 'H2', 'H1', 'H6']}
  SILO     weighted critical loss  53.93  customer-hours(k)   103.2
  CRITFLAG weighted critical loss  31.63  customer-hours(k)   148.5
  DEP      weighted critical loss   9.61  customer-hours(k)   163.2
  ORACLE   weighted critical loss   9.61
  DEP - CRITFLAG loss: mean -22.01 95% CI (-22.4, -21.59)
  DEP customer-hours vs SILO: +58.2% | vs CRITFLAG: +9.9%

5.2.2 Reading the results against the hypotheses

The honest summary is mixed. Criticality-aware ordering of any kind protects critical services dramatically compared with customer-first ordering, at a large cost in customer-hours that scheduling theory predicted (Chapter 3). Whether timing adds value over a critical-customer list depends on the region: none in Riverbend, a great deal in Port Halvard. The preregistered failure condition "H2 fails in one region and holds in the other" applies: narrow the claim to the conditions that distinguish the regions and test the mechanism.

5.2.3 Explaining the Riverbend anomaly without rescuing the hypothesis

In Riverbend, DEP restores F2 (telecom and hospital) before F3 (water), because telecom's expected time-to-impact (8 h) is slightly earlier than water's (8.7 h). CRITFLAG restores F3 first because it has more customers. But F3 protects two weighted services — water directly and the hospital through its water dependency — while restoring F2 first protects telecom and only the hospital's power path, which has a long buffer. An earliest-time-to-impact rule ignores weights and fan-out; Chapter 3 named this as the heuristic's known weakness. ORACLE's much lower loss confirms that neither rule is close to optimal in Riverbend.

This explanation is post hoc. It is legitimate to report it as a hypothesis generated by the confirmatory study, to be tested in new data (Chapter 6 tests it across many synthetic regions); it is not legitimate to re-run the main study with a revised rule and report that as confirmatory.

5.3 The operator tabletop study

5.3.1 Why people, not only models

A restoration rule that operators will not follow has no consequence. The tabletop study asks whether, given the effective time-to-impact display piloted in Chapter 4, experienced operators choose orderings that the model rates as lower-consequence, and why or why not.

5.3.2 Design (synthetic)

• Participants. Twelve operators and emergency managers (P01–P12) from three synthetic utilities, recruited through partners with individual consent, none from the pilot. Participation was voluntary and unpaid beyond expenses; managers did not learn who took part.
• Within-subject crossover. Each participant ordered restoration in two isomorphic scenarios — Port Halvard-like topology with renamed feeders and shuffled positions — once with a conventional display (customer counts and critical-customer flags) and once with the effective time-to-impact display added. Order of displays and scenario versions was counterbalanced.
• Blinding. Scenarios contained no outcome information; injects were checked in the pilot for premature cues. Analysts scoring the free-text rationales did not know which display condition produced them.
• Outcome. The chosen ordering, evaluated in the simulator on the same 500 draws as the main study — so the outcome is the model-estimated consequence of the participant's decision, not their stated preference. Secondary: decision time; rationale codes.
• Analysis. Paired comparison of whether each participant chose the DEP ordering under each display; exact McNemar test on discordant pairs; descriptive consequence differences; thematic analysis of rationales.

5.3.3 Results (synthetic)

All seven discordant pairs moved toward DEP. Under the null hypothesis that a switch in either direction is equally likely, the exact two-sided probability of seven or more switches in one direction out of seven is 2 × 0.5⁷ = 0.0156. With twelve participants from three utilities, this is evidence that the display can change choices in a tabletop, not an estimate of how often operators would change practice. Among the five who did not adopt DEP, four kept a CRITFLAG-like order. Rationale codes (two coders, disagreements resolved and logged) included: institutional expectation ("the board measures us on customers restored"; 4 participants), distrust of buffer figures ("I've seen hospital fuel logs — they're never what's on paper"; 3) and field logistics ("crews are near the big feeder already"; 3). The last theme identifies a variable the model omits: travel time between feeders, which makes restoration durations order-dependent.

5.3.4 What the tabletop can and cannot support

It supports: the display is comprehensible; it shifts tabletop choices toward lower modelled consequence; operators raise specific objections that the model must address. It cannot support: real-world adoption; effects under stress, fatigue and incomplete information; or effects in utilities unlike the three sampled. Tabletops are also performances — participants know they are observed — which pushes toward "correct-looking" answers.

5.4 Validating the model

Following Sargent (2013), validation is evidence that a model is adequate for its intended purpose — here, ranking restoration orderings by critical-service consequence — not that it is "true".

The historical validation produced an important correction: the model treats buffers as fixed, but real operators extend buffers — trucked water, fuel deliveries, load shedding at the hospital. A model without these human buffers is conservative about consequence and may overstate the benefit of reordering, because it ignores alternatives that operators would use anyway. The thesis must state this.

5.5 Independent validation in a second context

5.5.1 Selection on predeclared criteria

Port Halvard was chosen before simulation, on criteria derived from the mechanism and recorded in the preregistration amendment of Chapter 3: more critical feeders (four), divergence between buffer order and customer order among critical feeders, two crews, and long hospital buffers backed by fuel contracts. A different modeller built the Port Halvard configuration from partner documentation without access to Riverbend results. These steps reduce, but do not remove, the risk of a validation context that favours the method: the criteria themselves were chosen to give the mechanism a chance to operate. That is acceptable for a test of the mechanism, not for an estimate of how often the method helps across real regions. Chapter 6 addresses the second question with many randomly generated regions.

5.5.2 Transportability

A result transports from one context to another when the causal structure relevant to it is shared and the differences are measured. Pearl and Bareinboim (2014) formalize when experimental findings can be transported between populations. Informally, for the Cascade Study: the mechanism (buffers, sequential restoration, weighted consequence) is shared; the moderators that differ — number of critical feeders, order divergence, crew count, buffer lengths — must be measured in any target region before the result is applied. A utility that knows these four things about itself can estimate whether DEP is likely to beat its current list; one that does not, cannot.

5.6 Reporting the main study

Report, at minimum: the preregistration and all amendments with dates; the deviation log; every hypothesis with its verdict, including failures; effect sizes relative to smallest effects of interest; Monte Carlo intervals labelled as such; ORACLE gaps; tabletop design, consent, attrition and results with exact tests and qualitative themes; each validation type with what it establishes; and the two regions' configuration at the level of detail allowed by the information-protection tiers (Chapter 7). Follow a recognized reporting structure for simulation studies where one exists in the target venue, and for the tabletop follow the transparency expectations of human-subjects research in your field.

5.7 Common errors and limits

• Analysing two policies' means separately instead of paired differences.
• Increasing the number of draws until a difference is "significant".
• Presenting a Monte Carlo interval as uncertainty about the world.
• Re-running the confirmatory study with a rule revised after seeing results and calling it confirmatory.
• Choosing a validation context because the method will do well there, without saying so.
• Measuring operators' stated preferences instead of the consequences of their choices.

• Treating face validity as validation of parameters.

• Reading the Riverbend failure of H2 as a reason to revise DEP and rerun the main study under the same preregistration.

• Reporting a statistically significant tabletop shift as evidence of adoption, when twelve participants from three utilities only show that a display can change choices.
• Treating one historical event, such as the 2019 ice storm, as validation for cyber-induced scenarios.
• Omitting the ORACLE gap, which shows how far any rule is from the best ordering.

Exercises — chapter 5

EX027 — Run and interpret the main study. Objective: reproduce and interpret Lab 5. Inputs: the code in §5.2.1. Constraints: offline; record exact output. Steps: (1) run the lab; (2) fill the hypothesis table in §5.2.2 from your output; (3) compute the relative change of DEP versus CRITFLAG in each region; (4) write a 150-word results paragraph that a hostile reviewer would accept. Expected evidence: output, table, paragraph. Rubric: 3 = exact output, H2 reported as failed in Riverbend and supported in Port Halvard, H4a failure reported, no rescue; 2 = correct numbers, paragraph overstates; 1 = errors; 0 = not run.

EX028 — Paired versus unpaired. Objective: see why common random numbers matter. Inputs: Lab 5. Constraints: offline. Steps: (1) for Riverbend, compute the standard deviation of DEP losses, of CRITFLAG losses and of the per-draw differences; (2) compute an unpaired 95% interval for the difference of means using the normal approximation; (3) compare with the paired interval; (4) explain the difference in 100 words. Expected evidence: three standard deviations, two intervals, explanation. Rubric: 3 = SDs about 6.36 (DEP), 4.57 (CRITFLAG) and 6.79 (differences); unpaired interval about −0.34 to 1.04 versus paired −0.24 to 0.95; explains that pairing helped only modestly because the policies' losses are weakly correlated — they trade off which service fails on a given draw — and that the difference SD can exceed each policy's SD while still being below √(Var A + Var B); 2 = numbers right, explanation assumes pairing always shrinks variance; 1 = error; 0 = none.

EX029 — Uncertainty taxonomy. Objective: classify uncertainty. Inputs: five statements: (a) "with 5,000 draws the interval would narrow"; (b) "the hospital's fuel may be lower than the facilities director said"; (c) "travel time between feeders is not modelled"; (d) "the weights reflect one panel's values"; (e) "buffers may not be triangular". Constraints: classify each as Monte Carlo, parameter, structural or value uncertainty and name the method that addresses it. Steps: (1) classify each statement; (2) justify the class in one sentence; (3) name the method that reduces or exposes that uncertainty. Expected evidence: table. Rubric: 3 = (a) Monte Carlo, (b) parameter, (c) structural, (d) value, (e) parameter or structural with justification, with appropriate methods; 2 = one error; 1 = two or more; 0 = none.

EX030 — Tabletop design critique. Objective: strengthen a human-subjects design. Inputs: §5.3.2. Constraints: propose at least five improvements within ethical and practical limits. Steps: (1) identify the threats — demand characteristics, small sample, three utilities, the absence of time pressure, isomorphic-scenario recognition; (2) propose a design change for each; (3) state what each change costs. Expected evidence: table. Rubric: 3 = five or more real threats with feasible changes and honest costs; 2 = three or four; 1 = generic; 0 = none.

EX031 — Validation evidence matrix. Objective: plan validation for your own model. Inputs: a model from your research or the Cascade Study. Constraints: at least four validation types; for each, the evidence, what it establishes and what it does not. Steps: (1) choose at least four validation types from §5.4; (2) for each, state the evidence you would collect; (3) state what that evidence establishes and what it does not; (4) identify which validation would be most damaging to your claim if it failed. Expected evidence: matrix and paragraph. Rubric: 3 = types well chosen, limits honest, most damaging test identified and planned early; 2 = limits thin; 1 = face validity only; 0 = none.

EX032 — Transportability memo. Objective: advise a third utility. Inputs: a synthetic utility with five feeders, three carrying critical services whose buffer order matches their customer order, one crew. Constraints: 300 words; use the moderators in §5.5.2. Steps: (1) measure the four moderators for this utility; (2) predict DEP's gain over CRITFLAG; (3) state what data you would need to change your prediction. Expected evidence: memo. Rubric: 3 = predicts no gain because the orders coincide, explains the moderators, and notes that weights or fan-out could still matter; 2 = prediction right, reasoning thin; 1 = predicts a gain; 0 = none.

Chapter synthesis

The main study compared four orderings on paired draws with an oracle bound. Criticality-aware ordering of any kind protected critical services dramatically compared with customer-first restoration and cost far more customer-hours than the predeclared margin — a trade-off for policy, not a free gain. Timing information added nothing over a critical-customer list in Riverbend and a great deal in Port Halvard; the post hoc explanation (weights and fan-out ignored by an earliest-time rule) is a new hypothesis, not a rescue. A blinded tabletop showed that a time-to-impact display can shift operators' choices toward lower modelled consequence and surfaced omitted variables such as crew travel time. Validation established that the model is adequate for ranking orderings in familiar cases and conservative because it omits human buffers. Monte Carlo intervals described simulation noise only; parameter, structural and value uncertainty remain for Chapter 6.

Continue studying — chapter 05 of 08

← Chapter 04 · Course chapter index · Chapter 06 →

Continue to the next authored lesson.