Skip to content

Commit f4567b1

Browse files
jstacclaude
andauthored
Use bar for steady states and asterisk for maxima in MSY lecture (#765)
* Use bar for steady states and asterisk for maxima in MSY lecture Separate the two roles the asterisk was previously overloaded with: * steady states now use a bar --- the steady-state stock and sustainable yield become \bar{x}(e) and \bar{y}(e) * the asterisk is reserved for the maximizer --- the MSY effort becomes e^* The stock at MSY (= K/2) is the steady state at the optimal effort and is written \bar{x}(e^*) throughout. Code variables are renamed to match (x_star -> x_bar, e_msy -> e_star, x_msy -> x_bar_star). The RAM-database reference points B_{MSY} and F_{MSY} in the lingcod section are left as is. Co-Authored-By: Claude Opus 4.8 (1M context) <noreply@anthropic.com> * Fix backspace escape in yield-effort ylabel The ylabel string was a non-raw literal, so the \b in \bar was interpreted as a backspace before matplotlib parsed the mathtext, breaking notebook execution. Make it a raw string. Co-Authored-By: Claude Opus 4.8 (1M context) <noreply@anthropic.com> --------- Co-authored-by: Claude Opus 4.8 (1M context) <noreply@anthropic.com>
1 parent 6c511da commit f4567b1

1 file changed

Lines changed: 50 additions & 50 deletions

File tree

lectures/msy_fishery.md

Lines changed: 50 additions & 50 deletions
Original file line numberDiff line numberDiff line change
@@ -233,13 +233,13 @@ mystnb:
233233
---
234234
e_demo = 20.0 # an illustrative fixed effort level
235235
236-
def x_star(e):
236+
def x_bar(e):
237237
"Sustainable (steady-state) stock at effort e."
238238
return K * (1 - q * e / r)
239239
240240
fig, ax = plt.subplots(figsize=(4.95, 4.95))
241241
plot_45(ax, lambda x: update(x, e_demo), x0=80, x_max=1100,
242-
steady_state=x_star(e_demo), ss_label=r'$x^*(e)$',
242+
steady_state=x_bar(e_demo), ss_label=r'$\bar{x}(e)$',
243243
map_label=r'$x_{t+1}=x_t+G(x_t)-qex_t$')
244244
plt.tight_layout()
245245
plt.show()
@@ -248,7 +248,7 @@ plt.show()
248248
On the 45-degree diagram, fishing pulls the update curve down by the harvest
249249
term, so it now crosses the $45^\circ$ line at a lower steady state.
250250
251-
Since this steady state is a function of $e$ now, we denote it by $x^*(e)$.
251+
Since this steady state is a function of $e$ now, we denote it by $\bar{x}(e)$.
252252
253253
Here's the dynamics for two different levels of $e$, with the staircases omitted.
254254
@@ -265,7 +265,7 @@ fig, ax = plt.subplots(figsize=(4.95, 4.95))
265265
ax.plot(grid, grid, color='0.6', lw=1, ls='--', label=r'$45^\circ$ line')
266266
for e in (10.0, 30.0):
267267
line, = ax.plot(grid, update(grid, e), lw=2, label=f'$e={e:.0f}$')
268-
xs = x_star(e)
268+
xs = x_bar(e)
269269
ax.plot([xs], [xs], 'o', color=line.get_color(), ms=8, zorder=5)
270270
271271
ax.set_xlabel('stock this year $x_t$')
@@ -278,7 +278,7 @@ plt.tight_layout()
278278
plt.show()
279279
```
280280
281-
Not surprisingly, the steady state $x^*(e)$ is decreasing in fishing effort
281+
Not surprisingly, the steady state $\bar{x}(e)$ is decreasing in fishing effort
282282
(more fishing leads to less fish).
283283
284284
@@ -298,31 +298,31 @@ Writing it out and ignoring the empty-ocean case $x = 0$, we get
298298
$$
299299
r\left(1 - \frac{x}{K}\right) = qe
300300
\quad\Longrightarrow\quad
301-
x^*(e) = K\left(1 - \frac{qe}{r}\right).
301+
\bar{x}(e) = K\left(1 - \frac{qe}{r}\right).
302302
$$ (eq:xstar)
303303
304304
305-
Provided $e < r/q$, we see that each effort level $e$ pins down one steady-state stock $x^*(e)$.
305+
Provided $e < r/q$, we see that each effort level $e$ pins down one steady-state stock $\bar{x}(e)$.
306306
307307
The catch this steady state delivers is called the **sustainable yield**, and defined as
308308
309309
$$
310-
y^*(e) := q\,e\,x^*(e).
310+
\bar{y}(e) := q\,e\,\bar{x}(e).
311311
$$ (eq:yield)
312312
313313
```{code-cell} ipython3
314314
def sustainable_yield(e):
315315
"Sustainable catch delivered each year at effort e."
316-
return q * e * x_star(e)
316+
return q * e * x_bar(e)
317317
```
318318
319319
To visualize the sustainable yield, we plot the growth curve $G(x)$ together
320320
with the harvest line $q e x$ for a single effort level $e$.
321321
322322
The two cross where growth exactly replaces the catch: that crossing sits at the
323-
steady-state stock $x^*(e)$.
323+
steady-state stock $\bar{x}(e)$.
324324
325-
The height of the line there is the sustainable catch $y^*(e)$.
325+
The height of the line there is the sustainable catch $\bar{y}(e)$.
326326
327327
```{code-cell} ipython3
328328
---
@@ -337,13 +337,13 @@ fig, ax = plt.subplots()
337337
ax.plot(x, G(x), lw=2, label=r'growth $G(x)$')
338338
ax.plot(x, q * e_demo * x, lw=2, label=r'harvest $q e x$')
339339
340-
xs = x_star(e_demo)
340+
xs = x_bar(e_demo)
341341
ys = sustainable_yield(e_demo)
342342
ax.plot([xs], [ys], 'o', color='black', ms=8, zorder=5)
343343
ax.vlines(xs, 0, ys, ls='--', color='black', lw=1)
344344
ax.hlines(ys, 0, xs, ls='--', color='black', lw=1)
345-
ax.annotate(r'$x^*(e)$', xy=(xs, 0), xytext=(xs + 12, 8), fontsize=12)
346-
ax.annotate(r'$y^*(e)$', xy=(0, ys), xytext=(10, ys + 5), fontsize=12)
345+
ax.annotate(r'$\bar{x}(e)$', xy=(xs, 0), xytext=(xs + 12, 8), fontsize=12)
346+
ax.annotate(r'$\bar{y}(e)$', xy=(0, ys), xytext=(10, ys + 5), fontsize=12)
347347
348348
ax.set_xlabel('stock biomass $x$')
349349
ax.set_ylabel('catch')
@@ -354,7 +354,7 @@ plt.tight_layout()
354354
plt.show()
355355
```
356356
357-
Different effort levels give different sustainable catches $y^*(e)$.
357+
Different effort levels give different sustainable catches $\bar{y}(e)$.
358358
359359
In the next figure we plot the growth curve $G(x)$ together with the harvest
360360
lines $q e x$ for several values of $e$.
@@ -374,7 +374,7 @@ labels = [r'low $e$', r'moderate $e$', r'high $e$']
374374
375375
for e, lab in zip(efforts, labels):
376376
line, = ax.plot(x, q * e * x, lw=2, label=lab)
377-
ax.plot([x_star(e)], [q * e * x_star(e)], 'o', color=line.get_color(), ms=7, zorder=5)
377+
ax.plot([x_bar(e)], [q * e * x_bar(e)], 'o', color=line.get_color(), ms=7, zorder=5)
378378
379379
ax.set_xlabel('stock $x$')
380380
ax.set_ylabel('catch')
@@ -394,7 +394,7 @@ falls --- exactly the trade-off that the maximum sustainable yield captures.
394394
395395
### The yield-effort curve
396396
397-
To visualize the MSY, we plot the sustainable catch $y^*(e)$ against effort $e$.
397+
To visualize the MSY, we plot the sustainable catch $\bar{y}(e)$ against effort $e$.
398398
399399
This gives the classic dome-shaped Schaefer curve that fisheries managers use.
400400
@@ -409,9 +409,9 @@ e_grid = np.linspace(0, r / q, 400)
409409
y_grid = sustainable_yield(e_grid)
410410
411411
fig, ax = plt.subplots()
412-
ax.plot(e_grid, y_grid, lw=2, label=r'$y^*(e)=qeK\,(1-qe/r)$')
412+
ax.plot(e_grid, y_grid, lw=2, label=r'$\bar{y}(e)=qeK\,(1-qe/r)$')
413413
ax.set_xlabel('fishing effort $e$')
414-
ax.set_ylabel('sustainable yield $y^*(e)$')
414+
ax.set_ylabel(r'sustainable yield $\bar{y}(e)$')
415415
ax.set_xlim(0, r / q)
416416
ax.set_ylim(0, y_grid.max() * 1.25)
417417
ax.legend(loc='upper right', frameon=False)
@@ -431,12 +431,12 @@ It is the largest steady state catch attainable, assuming a constant effort rate
431431
The maximum sustainable yield can be defined more mathematically as
432432
433433
$$
434-
\text{MSY} \;=\; \max_{0 \le e < r/q} \; y^*(e).
434+
\text{MSY} \;=\; \max_{0 \le e < r/q} \; \bar{y}(e).
435435
$$
436436
437-
The effort level that solves this maximization problem is represented by $e_{MSY}$.
437+
The effort level that solves this maximization problem is represented by $e^*$.
438438
439-
We can compute it either by calculus or numerically, by evaluating $y^*(e)$ on a fine grid of effort
439+
We can compute it either by calculus or numerically, by evaluating $\bar{y}(e)$ on a fine grid of effort
440440
levels and picking the largest.
441441
442442
Below we look at the solution using calculus.
@@ -447,22 +447,22 @@ For now we will use the numerical approach:
447447
e_search = np.linspace(0, r / q, 100_001)
448448
i = np.argmax(sustainable_yield(e_search))
449449
450-
e_msy = e_search[i]
451-
MSY = sustainable_yield(e_msy)
452-
x_msy = x_star(e_msy)
450+
e_star = e_search[i]
451+
MSY = sustainable_yield(e_star)
452+
x_bar_star = x_bar(e_star)
453453
454-
print(f"effort at MSY e_MSY = {e_msy:.2f}")
455-
print(f"stock at MSY x* = {x_msy:.1f} tonnes (= K/2)")
454+
print(f"effort at MSY e* = {e_star:.2f}")
455+
print(f"stock at MSY = {x_bar_star:.1f} tonnes (= K/2)")
456456
print(f"maximum yield MSY = {MSY:.1f} tonnes/year (= rK/4)")
457457
```
458458
459-
The search lands on the round numbers $x^* = K/2$ and $\text{MSY} = rK/4$; the
459+
The search lands on the round numbers $\bar{x}(e^*) = K/2$ and $\text{MSY} = rK/4$; the
460460
optional calculus section below shows why.
461461
462462
463463
### Dynamics
464464
465-
What happens to biomass when we set $e = e_{MSY}$ from some arbitrary initial
465+
What happens to biomass when we set $e = e^*$ from some arbitrary initial
466466
$x_0$?
467467
468468
Below we run the yearly recursion {eq}`eq:update` forward from several
@@ -485,11 +485,11 @@ def simulate(x0, e, years=40):
485485
486486
fig, ax = plt.subplots()
487487
for x0 in [100, 300, 700, 950]:
488-
t, xt = simulate(x0, e_msy)
488+
t, xt = simulate(x0, e_star)
489489
ax.plot(t, xt, 'o-', ms=3, lw=2, label=f'$x_0={x0}$')
490490
491-
ax.axhline(x_msy, ls='--', color='black', lw=1)
492-
ax.annotate(r'$x^*=K/2$', xy=(0, x_msy), xytext=(1, x_msy + 25), color='black')
491+
ax.axhline(x_bar_star, ls='--', color='black', lw=1)
492+
ax.annotate(r'$\bar{x}(e^*)=K/2$', xy=(0, x_bar_star), xytext=(1, x_bar_star + 25), color='black')
493493
ax.set_xlabel('year $t$')
494494
ax.set_ylabel('stock biomass $x_t$')
495495
ax.set_xlim(0, 40)
@@ -498,7 +498,7 @@ plt.tight_layout()
498498
plt.show()
499499
```
500500
501-
We see that, under the MSY effort, every path climbs or falls toward $x^* = K/2$ and stays.
501+
We see that, under the MSY effort, every path climbs or falls toward $\bar{x}(e^*) = K/2$ and stays.
502502
503503
As a result, the catch settles down to the MSY.
504504
@@ -532,19 +532,19 @@ $$
532532
$$
533533
534534
The sustainability condition is unchanged --- a steady state still needs
535-
$G(x) = qex$ --- so $x^*(e) = K(1 - qe/r)$ and $y^*(e) = qeK(1 - qe/r)$ exactly as
535+
$G(x) = qex$ --- so $\bar{x}(e) = K(1 - qe/r)$ and $\bar{y}(e) = qeK(1 - qe/r)$ exactly as
536536
before.
537537
538-
Now $y^*$ is a smooth, concave parabola in $e$, and its peak is where the slope
538+
Now $\bar{y}$ is a smooth, concave parabola in $e$, and its peak is where the slope
539539
vanishes:
540540
541541
$$
542-
\frac{d y^*}{d e} = qK\left(1 - \frac{2qe}{r}\right) = 0
542+
\frac{d \bar{y}}{d e} = qK\left(1 - \frac{2qe}{r}\right) = 0
543543
\quad\Longrightarrow\quad
544-
e_{MSY} = \frac{r}{2q},
544+
e^* = \frac{r}{2q},
545545
$$
546546
547-
which gives $x^* = K/2$ and $\text{MSY} = rK/4$.
547+
which gives $\bar{x}(e^*) = K/2$ and $\text{MSY} = rK/4$.
548548
549549
550550
```{note}
@@ -577,12 +577,12 @@ To follow the fishery we use two ratios.
577577
The first is $B / B_{MSY}$, the stock biomass relative to the biomass $B_{MSY}$
578578
that supports the MSY.
579579
580-
(In our model this reference stock is $x^* = K/2$.)
580+
(In our model this reference stock is $\bar{x}(e^*) = K/2$.)
581581
582582
The second is $F / F_{MSY}$, fishing pressure relative to the pressure
583583
$F_{MSY}$ that achieves the MSY.
584584
585-
(In our model this is the MSY effort $e_{MSY}$.)
585+
(In our model this is the MSY effort $e^*$.)
586586
587587
The data come from the RAM Legacy Stock Assessment Database {cite:t}`ricard2012`.
588588
@@ -716,7 +716,7 @@ How should we set the catch $h_t$?
716716
717717
We compare two policies, both designed to take the MSY *on average*:
718718
719-
* **constant effort**: set $h_t = q\, e_{MSY}\, x_t$, so the catch scales with the
719+
* **constant effort**: set $h_t = q\, e^*\, x_t$, so the catch scales with the
720720
current stock, and
721721
* **constant quota**: set $h_t = \text{MSY}$ every year, regardless of the stock.
722722
@@ -725,18 +725,18 @@ yield each year".
725725
726726
It is also closer to how catch limits have often been set in practice.
727727
728-
Let's simulate both, starting every fishery from the MSY stock $x^* = K/2$.
728+
Let's simulate both, starting every fishery from the MSY stock $\bar{x}(e^*) = K/2$.
729729
730730
```{code-cell} ipython3
731-
def simulate_stochastic(policy, σ, years, rng, x0=x_msy):
731+
def simulate_stochastic(policy, σ, years, rng, x0=x_bar_star):
732732
"Simulate the stochastic fishery under a given harvest policy."
733733
x = np.empty(years + 1)
734734
x[0] = x0
735735
ξ = shocks(σ, years, rng)
736736
for t in range(years):
737737
growth = ξ[t] * G(x[t])
738738
if policy == 'effort':
739-
harvest = q * e_msy * x[t]
739+
harvest = q * e_star * x[t]
740740
else: # constant quota
741741
harvest = MSY
742742
x[t + 1] = max(0.0, x[t] + growth - harvest)
@@ -762,12 +762,12 @@ fig, axes = plt.subplots(2, 1, figsize=(6, 8), sharey=True)
762762
763763
for ax, policy, label in zip(
764764
axes, ['effort', 'quota'],
765-
['constant effort $h_t = q e_{MSY} x_t$',
765+
['constant effort $h_t = q e^* x_t$',
766766
'constant quota $h_t = \\mathrm{MSY}$']):
767767
for k in range(8):
768768
path = simulate_stochastic(policy, σ, years, rng)
769769
ax.plot(path, lw=2, alpha=0.8)
770-
ax.axhline(x_msy, ls='--', color='black', lw=1)
770+
ax.axhline(x_bar_star, ls='--', color='black', lw=1)
771771
ax.text(0.5, 0.92, label, transform=ax.transAxes, ha='center', va='top')
772772
ax.set_xlabel('year $t$')
773773
ax.set_ylabel('stock biomass $x_t$')
@@ -780,7 +780,7 @@ plt.show()
780780
781781
The difference is stark.
782782
783-
Under **constant effort**, the paths jiggle around $x^* = K/2$ but always recover:
783+
Under **constant effort**, the paths jiggle around $\bar{x}(e^*) = K/2$ but always recover:
784784
a bad year cuts the catch (because the catch is proportional to the stock), giving
785785
the population room to bounce back.
786786
@@ -912,7 +912,7 @@ def collapse_fraction_quota(quota, σ=0.15, years=100, n_paths=2000, seed=0):
912912
rng = np.random.default_rng(seed)
913913
collapses = 0
914914
for _ in range(n_paths):
915-
x = x_msy
915+
x = x_bar_star
916916
ξ = shocks(σ, years, rng)
917917
for t in range(years):
918918
x = max(0.0, x + ξ[t] * G(x) - quota)
@@ -1009,9 +1009,9 @@ The key formulas of the deterministic Schaefer model are collected below.
10091009
10101010
| Quantity | Formula | Value |
10111011
|---|---|---|
1012-
| Stock at MSY | $x_{MSY} = K/2$ | 500 t |
1012+
| Stock at MSY | $\bar{x}(e^*) = K/2$ | 500 t |
10131013
| Maximum sustainable yield | $\text{MSY} = rK/4$ | 125 t/yr |
1014-
| Effort at MSY | $e_{MSY} = r/2q$ | 25 |
1014+
| Effort at MSY | $e^* = r/2q$ | 25 |
10151015
10161016
The MSY is a deterministic, single-species, steady-state concept.
10171017

0 commit comments

Comments
 (0)