如何用R语言基于给定函数生成500代种群数据并绘制周期图？

我希望用以下公式生成500代的种群数量（y值）：

#3a 野兔初始种群变化计算
delta_H <- function(r1,Hn,Ln,d1)
  return(r1*Hn*(1-(Ln/d1)))
r1 <- 0.05
Hn <- 300000 
Ln <- 1500
d1 <- 1.5e+07
delta_H(r1,Hn,Ln,d1)

#3b 山猫初始种群变化计算
delta_L <- function(r2,Hn,Ln,d2)
  return(r2*Ln*(1-(Hn/d2)))
r2 <- -0.06
Hn <- 300000 
Ln <- 1500
d2 <- -1.8e+07
delta_L(r2,Hn,Ln,d2)

我尝试了以下代码，但没得到预期的周期图：

initial_hare_population <- 300000
Hn <- initial_hare_population
Hn <- 300000
initial_lynx_population <- 1500
Ln <- initial_lynx_population
Ln <- 1500
pop_Hn <- c(300000)
pop_Ln <- c(1500)
generations <- seq(0,500)
for (generation in generations) {
  delta_Hn <- delta_H(r1,Hn,Ln,d1)
  delta_Ln <- delta_L(r2,Hn,Ln,d2)
  Hn <- Hn+delta_Hn
  Ln <- Ln+delta_Ln
  pop_Hn <- append(pop_Hn,c(Hn))
  pop_Ln <- append(pop_Ln,c(Ln))
}
df <- data.frame(pop_Hn, pop_Ln, generations)
plot(df$generations,pop_Hn, xlab="野兔世代", ylab="野兔种群数量")
plot(df$generations, pop_Ln, xlab="山猫世代", ylab="山猫种群数量")

还尝试了这段代码：

curve(delta_H(r1,Hn,Ln,d1), from = 0, to = 500, n = 500, 
      xlab = "野兔世代", ylab = "野兔种群数量", 
      main = "500代野兔种群数量变化")

但R返回错误：

Error in curve(delta_H(r1, Hn, Ln, d1), from = 0, to = 500, n = 500, xlab = "野兔世代",  : 
  'expr' must be a function, or a call or an expression containing 'x'

请问如何用上述公式在R语言中绘制展示500代种群数量变化的周期图？

原代码问题分析 #

1. 长度不匹配：generations <- seq(0,500)会生成501个值，但初始pop_Hn和pop_Ln已有1个初始值，循环执行501次后，pop_Hn/pop_Ln会变成502个值，导致后续构建数据框时列长度不一致。
2. 冗余赋值：代码中重复对Hn和Ln赋值（比如Hn <- initial_hare_population后又Hn <- 300000），属于无效代码。
3. curve函数误用：curve要求传入的表达式必须包含变量x，但你的delta_H是计算单步种群变化的函数，不是关于世代x的函数，因此无法直接用curve绘制迭代过程的曲线。

修正后的完整代码 #

# 定义种群变化函数
delta_H <- function(r1, Hn, Ln, d1) {
  return(r1 * Hn * (1 - (Ln / d1)))
}

delta_L <- function(r2, Hn, Ln, d2) {
  return(r2 * Ln * (1 - (Hn / d2)))
}

# 设置参数
r1 <- 0.05
d1 <- 1.5e+07
r2 <- -0.06
d2 <- -1.8e+07

# 初始化种群
initial_hare <- 300000
initial_lynx <- 1500
Hn <- initial_hare
Ln <- initial_lynx

# 存储每代种群数量（初始值已加入）
pop_H <- c(initial_hare)
pop_L <- c(initial_lynx)

# 循环计算500代（从第1代到第500代）
for (i in 1:500) {
  delta_Hn <- delta_H(r1, Hn, Ln, d1)
  delta_Ln <- delta_L(r2, Hn, Ln, d2)
  Hn <- Hn + delta_Hn
  Ln <- Ln + delta_Ln
  pop_H <- c(pop_H, Hn)
  pop_L <- c(pop_L, Ln)
}

# 构建数据框，世代从0到500共501个值
generations <- 0:500
df <- data.frame(generation = generations, hare = pop_H, lynx = pop_L)

# 绘制周期图（双折线，展示两种群的周期变化）
plot(df$generation, df$hare, type = "l", col = "blue", 
     xlab = "世代", ylab = "种群数量", main = "500代野兔与山猫种群数量周期变化")
lines(df$generation, df$lynx, col = "red")
legend("topright", legend = c("野兔", "山猫"), col = c("blue", "red"), lty = 1)

代码说明 #

• 调整循环次数为1:500，确保pop_H/pop_L最终有501个值（对应0到500世代），和generations长度一致。
• 用lines()函数在同一张图上绘制两种群的数量变化，更直观展示捕食者与猎物的数量波动周期关系。
• 移除冗余赋值，代码更简洁清晰。

内容的提问来源于stack exchange，提问作者Peace