企业官网建设流程全解析

简介：hesim 是一套面向卫生经济评价与决策分析场景的 R 语言建模工具包，适合卫生经济学研究者、药物经济学评价人员及具备一定 R 基础的数据分析者使用。它提供将统计分析与经济评估相集成的通用框架，支持队列离散时间状态转换模型、N 状态分区生存模型以及个体级连续时间状态转换模型，覆盖 Markov 与半 Markov 过程，并借助 Rcpp 与 data.table 实现高效仿真、概率敏感性分析与患者异质性建模。资源包共 697 个文件，约 9.72MB，包含 60 个 R 源码、20 个 C++ 实现、102 个 rd 帮助文档、258 个 html 页面及若干 rmd 教程与 pdf 说明，便于查阅函数用法与模型示例。目前已有 471 人学习下载。读者可据此掌握从参数分布或自举采样、PSA 不确定性传播到成本效益分析与决策不确定性呈现的完整流程，并参考源码理解模型构建与计算加速思路。

1. 从一张成本效果散点图说起：hesim 到底解决什么问题

如果你做过卫生经济学评价，大概率经历过这样的场景：模型跑完，ICER 算出来是每 QALY 三万块，决策者问「这个结论稳不稳」，你只能硬着头皮说「我们做了单因素敏感性分析」。问题在于，单因素分析一次只动一个参数，而真实世界里成本、效用、转移概率是同时波动的。hesim 这个 R 包就是冲着这个痛点来的——它把微观模拟、离散事件模拟和分区生存模型统一到一套框架里，让你能跑概率敏感性分析、算成本效果可接受曲线，还能做个体层面的异质性建模。

这个包适合两类人：一类是做卫生技术评估的从业者，手里有 Markova 或 TreeAge 的模型想迁移到开源环境；另一类是流行病学或公共卫生背景的研究者，需要把干预措施的成本效果算清楚但不想被商业软件锁死。它的核心价值在于「概率化」——不是给你一个点估计，而是给你一个分布。你拿到的不只是「这个方案更划算」，而是「在 95% 的模拟里，这个方案有 78% 的概率更划算」。

hesim 的底层逻辑是「先定义模型结构，再生成模拟数据，最后做决策分析」。它不要求你一开始就写对全部参数，而是允许你用define_model()搭骨架，用define_parameters()填血肉，最后用run_model()跑出结果。这种分层设计的好处是，你可以在参数层反复调整而不动模型结构，也可以在模型层切换微观模拟和队列模拟而不重写参数。接下来我会按「装环境 → 搭模型 → 跑模拟 → 排错 → 进阶」的顺序，把每个环节的代码和参数含义拆开讲。

2. 环境准备与核心对象：把 hesim 装进工作流

2.1 安装与依赖检查

hesim 在 CRAN 上有稳定版本，但如果你要用最新的微观模拟功能，建议从 GitHub 装开发版。我一般会先检查 R 版本，低于 4.0 的话data.table和R6的某些特性会报错。

# 检查 R 版本，低于 4.0 建议先升级 getRversion() # 从 CRAN 安装稳定版 install.packages("hesim") # 如果需要开发版，用 devtools 从 GitHub 装 # install.packages("devtools") # devtools::install_github("hesim-dev/hesim") # 加载核心包 library(hesim) library(data.table)

装完之后别急着跑模型，先确认几个关键依赖：data.table用于高效处理模拟数据，R6用于面向对象建模，flexsurv用于参数生存模型拟合。如果flexsurv没装，后面做生存分析会直接卡住。

# 检查关键依赖是否可用 pkgs <- c("data.table", "R6", "flexsurv", "ggplot2") sapply(pkgs, requireNamespace, quietly = TRUE)

返回全是TRUE才算环境干净。有一个FALSE就去补装，别跳过。

2.2 理解 hesim 的三个核心对象

hesim 的建模流程围绕三个对象展开：hesim_data管数据，params管参数，model管结构。很多人第一次用觉得绕，是因为把这三个混在一起写了。我习惯先画一张表把职责分清：

hesim_data()是最先要建的。它要求你提供三个要素：strategies（干预策略名称）、patients（患者层面的协变量）、states（健康状态）。这三个要素的维度必须对齐，否则后面run_model()会报「维度不匹配」——这是新手最常见的翻车点。

# 定义三个干预策略 strategies <- data.table( strategy_id = 1:3, strategy_name = c("标准治疗", "新药A", "新药B") ) # 定义患者特征，这里模拟 1000 个患者 set.seed(42) patients <- data.table( patient_id = 1:1000, age = rnorm(1000, mean = 65, sd = 10), sex = rbinom(1000, 1, 0.5), risk_score = runif(1000, 0, 1) ) # 定义健康状态 states <- data.table( state_id = 1:3, state_name = c("稳定", "进展", "死亡") ) # 组装成 hesim_data 对象 hesim_dat <- hesim_data( strategies = strategies, patients = patients, states = states ) # 打印结构确认 print(hesim_dat)

这段代码的逻辑是：先把三个维度的数据分别建好，再用hesim_data()打包。set.seed(42)是为了让患者特征可复现，你做概率敏感性分析时如果每次患者样本都变，结果就没法对比。age用正态分布、sex用二项分布、risk_score用均匀分布，这是卫生经济评价里最常见的协变量生成方式。参数说明：strategies的strategy_id必须是整数且从 1 开始，patients的patient_id同理，states的state_id也要连续。这三个 ID 是 hesim 内部索引的锚点，断了就会报错。

3. 模型搭建：从转移概率到微观模拟

3.1 用define_model()定义分区生存结构

分区生存模型的核心是转移矩阵。hesim 允许你用两种方式定义：一种是直接给转移概率矩阵，另一种是用公式指定每个转移的回归关系。我一般先用矩阵方式跑通流程，再换成公式方式做异质性分析。

# 定义三个状态之间的转移矩阵 # 行是起点状态，列是终点状态 # 这里假设不能从死亡回到其他状态 transition_matrix <- matrix( c(0, 0.3, 0.1, # 从稳定出发：留在稳定0，去进展0.3，去死亡0.1 0, 0, 0.2, # 从进展出发：去稳定0，留在进展0，去死亡0.2 0, 0, 1), # 从死亡出发：只能留在死亡 nrow = 3, byrow = TRUE ) # 用 define_model 定义模型结构 mod <- define_model( transition_matrix = transition_matrix, n_cycles = 20, # 模拟 20 个周期 cycle_length = 1, # 每个周期 1 年 starting_state = 1 # 所有患者从稳定状态开始 ) # 查看模型摘要 summary(mod)

transition_matrix的每一行代表一个起点状态，每一列代表终点状态。第 1 行第 2 列的 0.3 表示从稳定到进展的概率是 0.3。注意第 3 行第 3 列的 1，表示死亡状态是吸收态，进去就出不来。n_cycles = 20意味着模拟 20 年，cycle_length = 1表示每年一个周期。如果你的研究时限是 10 年但每半年一个周期，就把n_cycles设成 20、cycle_length设成 0.5。

提示：转移矩阵的每一行概率之和必须等于 1，否则run_model()会报「概率不合法」。我见过有人把 0.3 和 0.1 填进去忘了补留在原状态的概率，结果跑出来患者凭空消失。

3.2 用define_parameters()绑定成本和效用

模型结构搭好后，下一步是把成本和效用参数绑上去。hesim 的参数对象支持三种形式：固定值、分布、回归方程。做概率敏感性分析时，你至少要把关键参数设成分布。

# 定义成本参数：每个状态每年的成本 cost_params <- params_lm( formula = cost ~ 1, data = data.frame( state_id = c(1, 2, 3), cost = c(5000, 15000, 0) # 稳定5000，进展15000，死亡0 ) ) # 定义效用参数：每个状态的健康效用值 utility_params <- params_lm( formula = utility ~ 1, data = data.frame( state_id = c(1, 2, 3), utility = c(0.8, 0.5, 0) # 稳定0.8，进展0.5，死亡0 ) ) # 把参数绑定到模型 mod <- set_parameters( mod, cost = cost_params, utility = utility_params )

params_lm()是线性模型参数容器，这里虽然只填了固定值，但它的结构支持后续替换成回归系数。formula = cost ~ 1表示截距模型，data里的state_id对应健康状态。成本参数里死亡状态设为 0 是惯例，因为死者不产生医疗费用。效用参数里死亡状态也设为 0，这是 QALY 计算的标准做法。

如果你要做概率敏感性分析，把固定值换成分布：

# 用正态分布替换固定值 cost_params_psa <- params_lm( formula = cost ~ 1, data = data.frame( state_id = c(1, 2, 3), mean = c(5000, 15000, 0), sd = c(500, 1500, 0) # 标准差设为均值的 10% ), uncertainty = "normal" # 指定分布类型 )

uncertainty = "normal"告诉 hesim 在每次模拟中从正态分布抽样。sd设成均值的 10% 是经验做法，如果你有文献来源的置信区间，用(上限 - 下限) / 3.92反推标准差更严谨。

3.3 跑微观模拟与结果提取

参数绑好后就可以跑模型了。hesim 支持两种模拟方式：run_model()跑队列模拟，run_microsim()跑个体模拟。做异质性分析必须用后者。

# 跑微观模拟，1000 个患者，100 次概率抽样 results <- run_microsim( mod, n_patients = 1000, n_sims = 100, seed = 123 ) # 提取成本效果结果 ce_results <- summarize_ce( results, type = "icer", comparator = "标准治疗" ) # 查看前几行 head(ce_results)

n_patients = 1000控制个体模拟的样本量，n_sims = 100控制概率敏感性分析的抽样次数。这两个参数直接决定运行时间——1000 患者乘 100 次模拟等于 10 万次个体轨迹，普通笔记本大概跑 2 到 5 分钟。seed = 123保证结果可复现。summarize_ce()的comparator参数指定基线策略，ICER 是相对于它算的。

跑完之后你会拿到一个包含strategy_id、cost、qaly、icer的 data.table。如果icer列出现NaN，说明该策略在所有模拟中都不比基线更有效，这是正常的，不是 bug。

4. 避坑与排查：那些让我重跑一整天的问题

4.1 转移概率矩阵行和不为 1

现象：run_model()报错Error: transition probabilities must sum to 1。原因：手动填矩阵时漏了留在原状态的概率。比如从稳定出发，你填了去进展 0.3、去死亡 0.1，但忘了留在稳定应该是 0.6。解决：用rowSums()检查每一行，或者用normalize()函数自动归一化。

# 检查行和 rowSums(transition_matrix) # 自动归一化（慎用，会掩盖逻辑错误） transition_matrix_norm <- transition_matrix / rowSums(transition_matrix)

我一般不建议自动归一化，因为它会把你的逻辑错误悄悄修掉。手动检查每一行，确认概率来源有依据。

4.2 患者 ID 与策略 ID 不连续

现象：hesim_data()创建成功，但run_microsim()报Error: patient_id must be consecutive integers。原因：你从外部 CSV 读入患者数据时，ID 可能跳号或者从 0 开始。解决：用setorder()排序后重新生成连续 ID。

# 重新生成连续 ID patients[, patient_id := .I] setorder(patients, patient_id)

.I是 data.table 的内置变量，返回当前行的序号。这招对策略 ID 和状态 ID 同样适用。

4.3 效用值超过 1 或为负

现象：QALY 算出来大于生命周期年限，或者出现负值。原因：效用参数没有做边界约束，正态分布抽样时抽到了大于 1 或小于 0 的值。解决：在params_lm()里加lower和upper约束，或者用 Beta 分布替代正态分布。

# 用 Beta 分布约束效用值在 0 到 1 之间 utility_params_beta <- params_lm( formula = utility ~ 1, data = data.frame( state_id = c(1, 2, 3), mean = c(0.8, 0.5, 0), sd = c(0.05, 0.08, 0) ), uncertainty = "beta" )

Beta 分布的参数化方式与正态不同，hesim 内部会用矩估计法转换。如果你的效用均值接近 0 或 1，Beta 分布比正态更合适。

4.4 模拟次数太少导致 ICER 不稳定

现象：每次跑run_microsim()得到的 ICER 差异很大，置信区间宽得没法用。原因：n_sims设得太小，概率空间没覆盖全。解决：先跑 10 次看 ICER 的变异系数，如果超过 20%，加到 100 次以上。

# 快速检查 ICER 稳定性 test_results <- run_microsim(mod, n_patients = 500, n_sims = 10, seed = 1) test_ce <- summarize_ce(test_results, type = "icer", comparator = "标准治疗") sd(test_ce$icer) / mean(test_ce$icer) # 变异系数

变异系数低于 0.1 说明模拟次数够了，高于 0.2 就继续加。别一上来就设 1000 次，先跑小样本探路。

4.5 生存模型外推超出数据范围

现象：用flexsurv拟合生存曲线后，外推到 20 年时生存概率变成负数。原因：参数生存模型在尾部可能不稳定，尤其是 Gompertz 或 Weibull 分布。解决：限制外推年限，或者换用限制性立方样条。

# 检查外推边界 fit <- flexsurvreg(Surv(time, status) ~ 1, data = survival_data, dist = "weibull") summary(fit, t = c(1, 5, 10, 20)) # 看 20 年时的生存概率

如果 20 年概率为负，把模型时限缩到 15 年，或者在define_model()里加max_cycles限制。

5. 进阶技巧：用概率敏感性分析做决策曲线

5.1 算成本效果可接受曲线

概率敏感性分析的最终产出不是 ICER 点估计，而是成本效果可接受曲线（CEAC）。它回答的是「在 willingness-to-pay 阈值下，每个策略成为最优的概率是多少」。

# 定义支付意愿阈值范围 wtp_range <- seq(0, 100000, by = 5000) # 算 CEAC ceac <- calc_ceac( results, wtp = wtp_range, comparator = "标准治疗" ) # 画图 library(ggplot2) ggplot(ceac, aes(x = wtp, y = prob, color = strategy_name)) + geom_line(size = 1) + labs(x = "支付意愿阈值 (元/QALY)", y = "成为最优策略的概率") + theme_minimal()

calc_ceac()的内部逻辑是：对每个 WTP 阈值，算每个策略的净货币收益（NMB = WTP × QALY - Cost），然后统计每个策略在 100 次模拟中 NMB 最高的频率。wtp_range从 0 到 10 万，步长 5000，覆盖了大多数中高收入地区的阈值范围。如果你的研究针对特定地区，把上限调到当地人均 GDP 的 3 倍即可。

5.2 用情景分析替代单因素敏感性分析

单因素敏感性分析一次只动一个参数，但决策者往往关心「如果新药价格降 20% 且效用提高 0.05 同时发生」会怎样。hesim 的情景分析可以同时改多个参数。

# 定义情景：新药价格降 20%，效用提高 0.05 scenario_params <- params_lm( formula = cost ~ 1, data = data.frame( state_id = c(1, 2, 3), cost = c(5000, 12000, 0) # 进展状态成本从 15000 降到 12000 ) ) # 把情景参数绑到模型 mod_scenario <- set_parameters(mod, cost = scenario_params) # 跑情景模拟 scenario_results <- run_microsim(mod_scenario, n_patients = 1000, n_sims = 100, seed = 123) scenario_ce <- summarize_ce(scenario_results, type = "icer", comparator = "标准治疗")

对比ce_results和scenario_ce的 ICER 变化，你就能量化「降价 + 提效」的联合影响。这比单因素分析更接近真实决策场景。

5.3 一个我踩过的坑：模拟种子与并行计算

最后说一个血泪经验。hesim 支持并行计算加速，但如果你在run_microsim()里设了seed又开了并行，不同核心上的随机数序列会冲突，导致结果不可复现。我当时的解决方案是：要么关并行，要么用future包显式管理随机数种子。

# 方案一：关并行，老老实实跑 results <- run_microsim(mod, n_patients = 1000, n_sims = 100, seed = 123) # 方案二：用 future 管理并行种子 library(future) plan(multisession, workers = 4) results <- run_microsim(mod, n_patients = 1000, n_sims = 100, seed = 123, parallel = TRUE) plan(sequential) # 跑完恢复默认

从那以后我每次跑概率敏感性分析都强制走一遍「小样本探路 → 检查变异系数 → 确认种子 → 再放大样本」的流程，再也没出现过结果对不上的情况。希望帮到你。

本文还有配套的精品资源，点击获取