[fisheries_lingcod] New lecture: Bayesian fisheries stock assessment — the lingcod application

[jstac]

Goal

Add a new lecture, fisheries_lingcod.md — Bayesian Fisheries Stock Assessment: The Lingcod Application — applying the SSPM estimation method validated in the companion methodology lecture (#919) to a real stock: Pacific Coast lingcod. It produces posteriors over the management reference points (MSY, B_MSY, F_MSY), the latent biomass trajectory, and current stock-status probabilities, then sketches the bridge to risk and harvest policy.

This is the second of two fisheries lectures. The methodology and code are built and validated in #919; this lecture recaps briefly and applies. It answers "what does the method say about a real stock?"

⚠️ Critical blocker (must be resolved before drafting)

Data. The SSPM needs an annual catch series C_t and a CPUE / biomass index I_t — different columns from the intro-series MSY lecture's bundled datasets/lingcod_msy_recovery.csv (which holds B/B_MSY, F/F_MSY ratios, not raw catch + index). Before any drafting we must:

• Identify a redistributable lingcod catch + index series (RAM Legacy v4.3 and/or the PFMC 2021 assessment), check terms of use, and vendor it as lectures/datasets/lingcod_sspm.csv so the GPU build downloads nothing at execution time.
• Confirm the RAM series IDs cited in the .tex (LINGPACcatch, LINGPACcpue) actually exist and contain a usable CPUE/biomass index — they look illustrative and need verification. RAM is primarily catch/biomass; a clean CPUE series may not be present.
• Decide the catch-series start date (run from ~1940 with a long catch-only depletion period, or from the first index year) and what the φ (initial depletion) prior then refers to.
• Decide which index to use if several exist (survey index preferred over commercial CPUE — more stable q); optionally run both as a sensitivity check.

This task gates this lecture only. The methodology lecture (#919) is fully executable on synthetic data and is not blocked by it.

Position in the sequence

Fifth in the Bayesian state-space arc: unemployment (#910) → population SSM (#911) → SSPM methodology (#919) → SSPM application (this lecture). A later lecture on harvest policy (POMDP / risk-sensitive control) is future work (see "Looking ahead" below).

The connection to the MSY lecture

The intro-series MSY lecture, msy_fishery.md, builds the deterministic Schaefer model in exactly our notation, tells the same lingcod story (its B/B_MSY, F/F_MSY recovery plot), and introduces risk via stochastic shocks, the constant-effort-vs-constant-quota comparison, the knife-edge collapse analysis, and a Monte-Carlo collapse probability. This lecture is its Bayesian completion:

• The MSY lecture takes r, K, q as known; we infer them and turn each reference point into a posterior rather than a number.
• Its closing message is "MSY is a deterministic point — steering at it ignores risk." Our dual message: the reference point itself is uncertain — because of the r–K ridge and q–K confounding, MSY = rK/4 is a distribution. A manager plugging in point estimates is overconfident twice: about the future (process noise — the MSY lecture's point) and about the model (parameter uncertainty — ours).
• We re-estimate the very stock the MSY lecture plots and overlay our posterior latent-biomass trajectory on its B/B_MSY story.

Section plan

1. Recap & setup — brief; cross-ref [fisheries_sspm] New lecture: Bayesian estimation of the Schaefer surplus production model — methodology & validation #919 (method) and the MSY lecture (model + lingcod story).
2. Lingcod data — the vendored datasets/lingcod_sspm.csv; alignment of catch and index series; log-transform; scaling.
3. Fit + parameter posteriors — marginals; the (r, K) ridge with the tight MSY = rK/4 posterior overlaid (well-identified combination along a poorly-identified ridge).
4. Latent biomass trajectory — credible bands; decline → low → post-2000 recovery; set against the MSY lecture's B/B_MSY recovery narrative.
5. Current stock status — P(B_T < 0.25K | data) (overfished) and P(B_T < B_MSY | data) (below target) directly from posterior samples.
6. Posterior predictive checks — residual distribution/autocorrelation; model criticism on real data; honest discussion of structural differences vs the age-structured SS3 assessment (a qualitative external benchmark, not a like-for-like comparison).
7. Looking ahead: risk and policy.
   • Executable now (stays in posterior-summary territory, no DP): the one-period-ahead predictive of B_{T+1} under a proposed catch C_T, and p_collapse(C_T) = P(B_{T+1} < B_lim | data, C_T) — the honest version of the MSY lecture's collapse_fraction, now marginalising over both process noise and parameter uncertainty.
   • Future work (flagged, not built): optimal harvest policy as a POMDP with the Bayesian posterior as belief state (vs the MSY lecture's two fixed policies); and risk-sensitive objectives (CVaR / entropic utility / chance constraints) versus risk-neutral expected discounted reward. Worth stating prominently: even risk-neutral expected-value maximisation already induces caution here — the transition is nonlinear and collapse is near-absorbing, so driving the stock low raises the chance of falling into the low-growth region and depresses all future expected catch. Risk-sensitive objectives make that downside aversion explicit and tunable.

Correctness notes carried from #919

The model, code, and most correctness fixes live in #919 and are inherited here. Items that surface specifically in the application:

• Filtering vs smoothing: NUTS on the full model returns the smoothing posterior p(B_t | I_{1:T}) (with θ marginalised), not sequential filtering. It coincides with the filtering distribution only at the terminal time t = T — which is exactly the belief state used in the "Looking ahead" section, so the use is fine, but state this precisely rather than calling the full-sample smoother "the filtering distribution."
• Don't oversell differentiability/optimisation: the posterior samples are not differentiable w.r.t. policy; only the one-period transition mean B_{T+1}(C_T) is. And the single-period "max E[C_T] s.t. CVaR ≤ δ" is a scalar monotone constraint (a 1-D root-find), not a JAX-gradient optimisation. Keep the framing modest.

Cross-referencing & style guide (note for eventual alignment)

Same setup as #919: {ref}/{doc} within this series; plain full URLs across series (the MSY lecture at https://intro.quantecon.org/msy_fishery.html). The style-guide repo is the successor to the older QuantEcon.manual (which documents cross-referencing conventions); the whole Bayesian sequence should eventually be brought into line with it — flagging here.

Tasks

• Resolve the data blocker above (vendor datasets/lingcod_sspm.csv; verify series; decide start date and index)
• Fit the validated [fisheries_sspm] New lecture: Bayesian estimation of the Schaefer surplus production model — methodology & validation #919 model to the real series; full diagnostics
• (r, K) ridge plot with the tight MSY = rK/4 posterior overlaid
• Latent-biomass credible bands; overlay against the MSY lecture's B/B_MSY story
• Stock-status probabilities P(overfished), P(below target)
• Posterior predictive checks; qualitative SS3 sanity comparison
• One-period predictive + p_collapse(C_T) (the executable risk piece)
• "Looking ahead" section: optimal policy + risk-sensitivity as future work
• Draft lectures/fisheries_lingcod.md in MyST (GPU admonition include); add to _toc.yml after [fisheries_sspm] New lecture: Bayesian estimation of the Schaefer surplus production model — methodology & validation #919
• Cross-link [fisheries_sspm] New lecture: Bayesian estimation of the Schaefer surplus production model — methodology & validation #919, the MSY lecture (full URL), and the rest of the Bayesian sequence
• Verify execution (jupytext → py); check GPU-build runtime budget

Gating: depends on #919 (method + validated code) and on the data blocker. Original draft source (fisheries_sspm_lecture.tex) is embedded in #919.