七鰓鰻族群動力學建模

這道數學建模題研究七鰓鰻的性別比如何隨環境資源改變，以及這種適應能力對族群與整個生態系統有何影響。

問題背景

很多物種出生時接近 1:1 性別比，亦有物種會因環境而偏離，稱為適應性性別比分配。例如孵化溫度可影響美洲短吻鱷性別。

七鰓鰻在不同地區的生態與文化角色複雜：某些湖泊把牠們視為影響魚類的寄生生物，斯堪的納維亞、波羅的海及北美太平洋西北部分原住民族群亦會把牠們作食物。

海七鰓鰻幼體生長速度受食物供應影響，並與性別分化相關。資源低時，雄性比例可接近 78%；資源較充足時，觀測雄性比例約 56%。題目要求分析：

可變性別比如何影響較大生態系統？
對七鰓鰻族群有何利弊？
性別比改變如何影響系統穩定性？
其他生物是否可能從這種改變獲益？

由捕食者—獵物模型開始

經典 Lotka–Volterra 系統為

\frac{dx}{dt}=\alpha x-\beta xy, \qquad \frac{dy}{dt}=\delta xy-\gamma y,

其中 $x$ 是獵物、 $y$ 是捕食者， $\alpha$ 是獵物內在增長率， $\beta$ 是捕食強度， $\delta$ 把攝食轉成捕食者增長， $\gamma$ 是捕食者死亡率。

這個模型適合建立互動直覺，卻未包含資源上限、生命階段、性別限制繁殖、季節、遷徙或密度依賴。若研究問題核心是性別比，模型必須令性別比真正改變出生數或配偶限制；只把性別比畫成一條輸出曲線並不足夠。

合成 Python 模型

先匯入函式庫並固定種子：

import numpy as np
import matplotlib.pyplot as plt

np.random.seed(32)

基本族群類別為：

class Organism:
    def __init__(self, initial_population):
        self.population = initial_population

    def update_population(self, birth_rate, mortality_rate):
        self.population += int(
            self.population * (birth_rate - mortality_rate)
        )
        self.population = max(0, self.population)

七鰓鰻類別保存食物與性別比，並以 logistic 增長和捕食項更新族群：

class LampreyPopulation(Organism):
    def __init__(self, initial_population, initial_food):
        super().__init__(initial_population)
        self.food = initial_food
        self.male_ratio = 0.5
        self.female_ratio = 0.5

    def update_sex_ratio(self, temperature):
        temp_factor = np.clip((temperature - 10) / 20, 0, 1)
        mean_ratio = 0.78 - 0.22 * (self.food / 1000)
        noisy_ratio = (
            mean_ratio - 0.1 * temp_factor
        ) * np.random.normal(1, 0.1)
        self.male_ratio = np.clip(noisy_ratio, 0.56, 0.78)
        self.female_ratio = 1 - self.male_ratio

    def reproduce(self, birth_rate):
        return int(
            self.population * self.female_ratio * birth_rate
        )

    def update_population(self, predator_population, r, K, a):
        delta = (
            r * self.population * (1 - self.population / K)
            - a * predator_population * self.population / K
        )
        self.population = max(0, self.population + int(delta))

    def update_food(self, algae_population, consumption_rate):
        self.food = max(
            0,
            algae_population * 0.1
            - self.population * consumption_rate
        )

相對原稿，上面把噪音後的雄性比例再次限制在 $[0.56,0.78]$ ；原寫法在 clip 後乘隨機數，實際仍可越界。不過加入溫度影響並沒有由題目提供的七鰓鰻證據支持，只可視為待檢驗假設。

捕食者與藻類類別為：

class Predator(Organism):
    def update_population(self, lamprey_population, a, b, d, K):
        delta = (
            b * a * lamprey_population * self.population / K
            - d * self.population
        )
        self.population = max(0, self.population + int(delta))


class Algae(Organism):
    def __init__(self, initial_population, carrying_capacity):
        super().__init__(initial_population)
        self.carrying_capacity = carrying_capacity

    def grow(self, growth_rate, temperature, nutrient_level):
        temp_factor = np.clip((temperature - 5) / 25, 0, 1)
        nutrient_factor = min(1, nutrient_level / 100)
        effective_rate = (
            growth_rate * temp_factor * nutrient_factor
        )
        delta = (
            self.population * effective_rate
            * (1 - self.population / self.carrying_capacity)
        )
        self.population = np.clip(
            self.population + int(delta),
            0,
            self.carrying_capacity
        )

模擬器更新溫度、營養、藻類、食物、性別比、七鰓鰻與捕食者：

def simulate_ecosystem(years, lamprey, predator, algae, params):
    history = {
        'lamprey_pop': [], 'predator_pop': [],
        'algae_pop': [], 'male_ratio': [],
        'temperature': [], 'nutrient_level': []
    }

    for year in range(years):
        phase = np.random.uniform(0, 1)
        temperature = 20 + 10 * np.sin(
            2 * phase * np.pi * year / 4
        )

        params['nutrient_level'] = np.clip(
            params['nutrient_level']
            + np.random.normal(0, 5),
            0, 100
        )

        algae.grow(
            params['algae_growth_rate'],
            temperature,
            params['nutrient_level']
        )
        lamprey.update_food(
            algae.population,
            params['food_consumption_rate']
        )
        lamprey.update_sex_ratio(temperature)

        lamprey.update_population(
            predator.population,
            params['lamprey_growth_rate'],
            params['lamprey_carrying_capacity'],
            params['predation_rate']
        )
        predator.update_population(
            lamprey.population,
            params['predation_rate'],
            params['predator_reproduction_rate'],
            params['predator_death_rate'],
            params['lamprey_carrying_capacity']
        )

        if year % 5 == 0:
            lamprey.population = int(
                lamprey.population
                * (1 - params['fishing_rate'])
            )
            params['nutrient_level'] += (
                params['pollution_rate'] * 10
            )

        history['lamprey_pop'].append(lamprey.population)
        history['predator_pop'].append(predator.population)
        history['algae_pop'].append(algae.population)
        history['male_ratio'].append(lamprey.male_ratio)
        history['temperature'].append(temperature)
        history['nutrient_level'].append(
            params['nutrient_level']
        )

    return history

參數及初始值：

params = {
    'lamprey_growth_rate': 0.3,
    'lamprey_carrying_capacity': 10000,
    'predation_rate': 0.2,
    'predator_reproduction_rate': 0.35,
    'predator_death_rate': 0.025,
    'algae_growth_rate': 0.1,
    'food_consumption_rate': 0.1,
    'nutrient_level': 50,
    'fishing_rate': 0.3,
    'pollution_rate': 0.2
}

lamprey = LampreyPopulation(1000, 2000)
predator = Predator(100)
algae = Algae(10000, 50000)

history = simulate_ecosystem(
    400, lamprey, predator, algae, params
)

輸出可畫成族群、藻類、雄性比例、溫度、營養及七鰓鰻變化率：

fig, axes = plt.subplots(3, 2, figsize=(16, 16))

axes[0, 0].plot(history['lamprey_pop'], label='Lamprey')
axes[0, 0].plot(history['predator_pop'], label='Predator')
axes[0, 0].legend()

axes[0, 1].plot(history['algae_pop'])
axes[1, 0].plot(history['male_ratio'])
axes[1, 1].plot(history['temperature'])
axes[2, 0].plot(history['nutrient_level'])

lamprey_values = np.asarray(history['lamprey_pop'])
rate = np.divide(
    np.diff(lamprey_values),
    lamprey_values[:-1],
    out=np.full(len(lamprey_values)-1, np.nan),
    where=lamprey_values[:-1] != 0
)
axes[2, 1].plot(rate)

plt.tight_layout()
plt.show()

原專案圖為：

模型實際回答了甚麼

目前程式可以產生帶資源、捕食、捕撈與環境波動的合成族群軌跡，但仍未回答題目核心。reproduce() 雖然定義了雌性比例對出生數的作用，卻從未在 simulate_ecosystem() 中呼叫；七鰓鰻總增長只由固定 lamprey_growth_rate 決定。因此，改變雄性比例不會回饋到七鰓鰻出生或族群更新。

換言之，圖中的性別比與族群同時變化，不等於性別比造成族群改變。要研究 adaptive sex ratio，至少要把出生項改寫成

B_t=b\,L_t\,g(F_t,M_t),

其中 $F_t=(1-r_t)L_t$ 、 $M_t=r_tL_t$ ， $g$ 表示配偶限制，例如 $\min(F_t,\eta M_t)$ 或其他有生物根據的 mating function。其後應比較：

固定 1:1 性別比；
固定觀測平均比例；
隨資源改變的性別比；
在相同環境隨機序列下的族群持續性、變異、滅絕風險與生態影響。

還要校準時間單位、生命階段、單次產卵後死亡、幼體期長度與遷徙。現有年度整數更新及任意參數只適合作為程式示範。

結語

捕食者—獵物模型為七鰓鰻生態問題提供一個起點，但真正研究焦點不是畫出幾條振盪曲線，而是讓資源、性別分化、繁殖限制與物種互動形成可檢驗的因果鏈。

目前合成程式展示如何組織多族群模擬，也揭示一項重要建模教訓：若一個變量沒有進入決定研究結果的方程，它即使被計算和繪圖，也不能支持關於其影響的結論。