如何识别R语言ABM新冠模型函数的性能瓶颈以提升运行速度？

架构老林

2026-4-30

如何分析并优化R语言基于主体的模型（ABM）的性能瓶颈？

我跟着Ed Boone的YouTube视频实现了一个新冠疫情的ABM模型，代码在1000个样本时运行正常，但扩大样本量后速度骤降。我试过手动加Sys.time()来计时，但代码里有多层嵌套for循环，这种方法太复杂了——还要考虑每行代码的执行次数，根本不方便。想问问有没有更高效的方法定位性能瓶颈，然后优化这些函数的运行速度？

以下是我的代码：

# ABM_Covid - II
start.time <- Sys.time()
Data_Generator <- function(nPop1, E0, I0) {
  # Create a population of susceptibles
  Data <- data.frame(
    AgentNo=1:nPop1,
    State="Susceptible",
    Mixing= runif(nPop1,0,1),
    TimeE = 0,
    TimeI = 0,
    stringsAsFactors = FALSE)
  Data$State[1:E0] <- "Exposed"
  # This just says that the first person is exposed, since the mixing is random anyway, this is not an issue (because the mixing of Exposed is random)
  Data$Time[1:E0] <- rbinom(E0, 13, 0.5) + 1 # Exposure up to 14 days
  Data$State[(E0+1):(E0+I0)] <- "Infected"
  Data$Time[(E0+1):(E0+I0)] <- rbinom(I0, 12, 0.5)
  return(Data)
}
ABM_Covid <- function(Data, parameters, runtime){
  nPop1 <- nrow(Data)
  # runtime <- 15
  Results <- data.frame(
    Susceptible = rep (0, runtime),
    Exposed = rep (0, runtime),
    Infected = rep (0, runtime),
    Recovered = rep (0, runtime),
    Deaths = rep (0, runtime))
  # Move people through time
  for (k in 1:runtime){
    # Moving people through time
    StateSusceptible <- (1:nPop1)[Data$State == "Susceptible"]
    StateSusceptible_or_Exposed <- (1:nPop1)[Data$State == "Susceptible" | Data$State == "Exposed"]
    for (i in StateSusceptible) {
      # Determine if they like to meet others
      Mix1 <- Data$Mixing[i]
      # How many agents will they meet? The plus one meets everybody meets somebody
      Meetings <- round(Mix1*parameters$MaxMix,0) + 1
      # Grab the agents they will meet
      People_met <- sample(StateSusceptible_or_Exposed, Meetings, replace=TRUE, prob = Data$Mixing[StateSusceptible_or_Exposed])
      for (j in 1:length(People_met)) {
        # Grab who they will meet
        Meetingsa <- Data[People_met[j], ]
        # If exposed change State
        if(Meetingsa$State== "Exposed") {
          Urand1 <- runif(1,0,1)
          if (Urand1 < parameters$S2E){
            Data$State[i] <- "Exposed"
          }
        }
      }
    }
    # Grab those who have been exposed and increment
    StateE1 <- (1:nPop1)[Data$State== "Exposed"]
    Data$TimeE[StateE1] = Data$TimeE[StateE1] + 1
    StateE2 <- (1:nPop1)[Data$State== "Exposed" & Data$TimeE > 14]
    Data$State[StateE2] <- "Recovered"
    # Grab those who could become sick
    StateE3 <- (1:nPop1)[Data$State== "Exposed" & Data$TimeE > 3]
    for (i in StateE3){
      Urand1 <- runif(1,0,1)
      # randomly assign whether they get sick or not
      if ( Urand1 < parameters$E2I ) {
        Data$State[i] <- "Infected"
      }
    }
    # Update how long they have been sick
    StateI1 <- (1:nPop1)[Data$State== "Infected"]
    Data$TimeI[StateI1] = Data$TimeI[StateI1] + 1
    # Recovered bin
    StateI2 <- (1:nPop1)[Data$State== "Infected" & Data$TimeI > 14]
    Data$State[StateI2] <- "R"
    # Not recovered could potentially die
    StateI3 <- (1:nPop1)[Data$State== "Infected" & Data$TimeI < 15]
    Data$State[StateI3] <- ifelse(runif(length(StateI3), 0, 1 ) > parameters$I2D, "Infected", "Deaths")
    Results$Susceptible[k] <- length(Data$State[Data$State=="Susceptible"])
    Results$Exposed[k] <- length(Data$State[Data$State=="Exposed"])
    Results$Infected[k] <- length(Data$State[Data$State=="Infected"])
    Results$Recovered[k] <- length(Data$State[Data$State=="Recovered"])
    Results$Deaths[k] <- length(Data$State[Data$State=="Deaths"])
  }
  return(Results)
}
Data <- Data_Generator(1000, E0=5, I0=2)
parameters <- data.frame(
  MaxMix = 10,
  S2E = 0.25,
  E2I = 0.1,
  I2D = 0.1)
Model1 <- ABM_Covid(Data, parameters, runtime=25)
plot(1:25, Model1$Susceptible, type="l", col="purple", ylim = c(0,1000))
lines(1:25, Model1$Exposed, type="l", col="orange")
lines(1:25, Model1$Infected, type="l", col="red")
lines(1:25, Model1$Recovered, type="l", col="seagreen")
lines(1:25, Model1$Deaths, type="l", col="black")
end.time <- Sys.time()
time.taken <- end.time - start.time

第一步：用专业工具定位性能瓶颈

手动加计时语句确实麻烦，R有专门的性能分析工具，能帮你快速找到最慢的代码块：

profvis包：这是可视化性能分析的首选工具，能直观展示每行代码的耗时和调用次数。使用方法很简单：
1. 安装并加载包：install.packages("profvis"); library(profvis)
2. 把你的模型运行代码包在profvis()里：
```
profvis({
  Data <- Data_Generator(5000, E0=5, I0=2) # 用大样本测试
  parameters <- data.frame(MaxMix=10, S2E=0.25, E2I=0.1, I2D=0.1)
  Model1 <- ABM_Covid(Data, parameters, runtime=25)
})
```
  运行后会弹出一个交互式界面，能看到哪个循环、哪行代码占了最多的运行时间——大概率是你那三层嵌套的for循环。
microbenchmark包：如果想对比某几段代码的运行速度，这个包很合适，能给出多次运行的统计数据，帮你确认优化效果。

第二步：针对性优化你的ABM代码

从你的代码来看，性能问题主要来自嵌套for循环和逐行操作数据框——这在R里是效率很低的做法，因为R更擅长向量化操作。这里给几个具体的优化方向：

1. 把内层的人际接触循环改成向量化操作

你现在的代码里，每个易感者要循环遍历自己遇到的所有人，逐个判断是否被感染。可以改成一次性处理所有接触事件：

# 替换原来的StateSusceptible循环部分
if(length(StateSusceptible) > 0) {
  # 批量生成每个易感者的接触人数
  n_meetings <- round(Data$Mixing[StateSusceptible] * parameters$MaxMix, 0) + 1
  # 生成所有接触事件：每个易感者对应n_meetings个接触对象
  contact_df <- data.frame(
    susceptible_id = rep(StateSusceptible, n_meetings),
    contact_id = sample(StateSusceptible_or_Exposed, sum(n_meetings), replace=TRUE, prob=Data$Mixing[StateSusceptible_or_Exposed])
  )
  # 筛选接触对象是Exposed的行
  exposed_contacts <- contact_df[Data$State[contact_df$contact_id] == "Exposed", ]
  # 批量生成随机数判断是否感染
  exposed_contacts$infected <- runif(nrow(exposed_contacts)) < parameters$S2E
  # 去重：每个易感者只要被感染一次就够了
  newly_exposed <- unique(exposed_contacts$susceptible_id[exposed_contacts$infected])
  # 一次性更新状态
  Data$State[newly_exposed] <- "Exposed"
}

这样就把三层循环改成了向量化操作，速度会快很多。

2. 用`data.table`代替普通数据框

普通数据框的逐行修改效率很低，换成data.table（专门用于高效处理大数据集的包），能大幅提升状态更新的速度。比如把Data_Generator改成：

library(data.table)
Data_Generator <- function(nPop1, E0, I0) {
  Data <- data.table(
    AgentNo=1:nPop1,
    State="Susceptible",
    Mixing= runif(nPop1,0,1),
    TimeE = 0,
    TimeI = 0
  )
  Data[1:E0, State := "Exposed"]
  Data[1:E0, Time := rbinom(E0, 13, 0.5) + 1]
  Data[(E0+1):(E0+I0), State := "Infected"]
  Data[(E0+1):(E0+I0), Time := rbinom(I0, 12, 0.5)]
  return(Data)
}

data.table的:=操作符是原地修改，比普通数据框的$赋值快得多。

3. 减少重复的子集筛选

你的代码里每次循环都要多次筛选状态子集（比如StateSusceptible、StateE1等），可以把这些筛选结果存起来重复使用，避免重复计算：

for (k in 1:runtime){
  # 一次性筛选所有状态索引
  idx_S <- Data$State == "Susceptible"
  idx_E <- Data$State == "Exposed"
  idx_I <- Data$State == "Infected"
  idx_S_or_E <- idx_S | idx_E
  
  # 后面直接用这些逻辑向量，比如StateSusceptible就是which(idx_S)
  # ... 其余代码修改为使用这些索引 ...
}

这样能减少重复的子集计算开销。

4. 用预生成的随机数代替循环内生成

你现在在循环里多次调用runif(1)，可以提前生成所有需要的随机数，然后在循环里取用——生成随机数本身也有开销，批量生成更高效。

第三步：验证优化效果

每次优化后，用profvis或者microbenchmark对比优化前后的运行时间，确认瓶颈已经被解决。比如用microbenchmark测试：

library(microbenchmark)
microbenchmark(
  original = {
    Data <- Data_Generator(5000, E0=5, I0=2)
    Model1 <- ABM_Covid(Data, parameters, runtime=25)
  },
  optimized = {
    Data_opt <- Data_Generator(5000, E0=5, I0=2) # 优化后的生成函数
    Model_opt <- ABM_Covid_opt(Data_opt, parameters, runtime=25) # 优化后的模型函数
  },
  times = 5
)

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

如何识别R语言ABM新冠模型函数的性能瓶颈以提升运行速度？

第一步：用专业工具定位性能瓶颈 #

第二步：针对性优化你的ABM代码 #