在疾控和流行病学相关领域待过的人大概率都遇到过这种场景拿着一份病例数据图上明明看得出来有几个乡镇的颜色比周围深但领导或者审稿人问你“这个聚集是真实存在的还是随机波动”的时候你没法用一张色阶图去回答。这时候就需要上空间扫描统计而业内最常用、也绕不开的工具就是SaTScan。SaTScan这个名字我最早是在一篇研究某地区手足口病时空聚集的论文里看到的后来自己接手一个县级市的食源性疾病数据分析项目时正式用了它。这篇文章不打算复读官方手册而是从实际使用的角度把从数据准备、参数设置、结果解读到地图导出的完整链路捋一遍顺便把那些官方文档里没写、但你在实操中一定会撞上的坑讲清楚。无论你是公共卫生领域的研究生、疾控中心的业务人员还是做疾病监测、动物疫病甚至生态学的数据分析这篇文章都可以当一份上手参考。1. 为什么空间扫描统计不可或缺以及 SaTScan 到底做了什么1.1 传统统计方法在地理聚集检测上的局限性在接触 SaTScan 之前大多数人判断“有没有聚集”用的方法无非两种一种是直接看地图把发病率按乡镇或街道分级填色颜色深的自然觉得有问题另一种是把整个研究区域当成一个整体算一个全局发病率再看每个小地区的发病率跟全局均值的差异。这两种做法的问题都很明显。地图填色只能看出“哪里颜色不一样”但给不出一个统计检验的结论而且颜色深浅很容易受到分级方式的影响——你把发病率分成五级和分成七级视觉上可能得出完全不同的结论。全局比较的做法也有缺陷一个县有几十个乡镇逐个做假设检验的话多重比较会让假阳性的概率大幅上升更麻烦的是地理数据本身存在空间自相关紧挨着的两个乡镇发病率往往互相影响简单统计模型很难处理这种依赖关系。空间扫描统计解决的就是这个问题。它不是在分析前就指定“我觉得这几个镇有问题”而是让一个动态窗口在整个研究区域内反复扫描每到一个位置就用似然比检验来判断该窗口内的病例数是否显著高于窗口外。这样一来不用预先假设聚集的位置和大小还能同时输出统计学显著性指标也就是 P 值。SaTScan 就是这套方法最成熟的软件实现。1.2 扫描窗口、似然比检验和蒙特卡罗模拟SaTScan 的三大支柱要顺利使用 SaTScan至少得理解它背后的三个核心概念这三个概念直接决定了你在软件里怎么设参数。第一是扫描窗口。SaTScan 的默认窗口是圆形圆心你可以理解为研究区域内任意一个可能的位置点圆的半径从小到大不停变化每换一次半径就构成一个新的候选聚集区域。除了圆形新版本也支持椭圆扫描椭圆的好处是能捕捉条带状、沿道路或河流走向分布的聚集我后面会单独讲。第二是似然比检验。对每一个候选窗口SaTScan 都会计算一个似然比Log Likelihood RatioLLR。通俗地说这个值衡量的是“如果窗口内的发病风险确实比窗外更高数据出现当前情况的可能性”和“假设整个区域风险一致数据出现当前情况的可能性”之间的比值。窗口内观测病例数相对期望病例数偏离越大LLR 就越大这个窗口是真正聚集的可能性也就越高。第三是蒙特卡罗模拟。因为地理数据的统计分布很难用理论公式直接算 P 值SaTScan 的做法是保持每个位置的人口数和病例总数不变但把所有病例的位置随机打乱成千上万次每次重新扫描一遍统计模拟出的最大 LLR 有多大。如果你的真实数据中最大 LLR 超过了绝大多数模拟结果就说明这个聚集不太可能是随机产生的。模拟次数就是我们常说的蒙特卡罗重复次数这个参数后面会细说。搞懂这三个概念SaTScan 的参数面板就不再是一堆看不懂的按钮而是一套可以直接对应到具体决策的选项了。2. 跑通第一轮分析前必须先吃透的数据格式SaTScan 上手最大的门槛不是软件本身而是数据文件的组织和清洗。这个软件对输入格式有严格要求而且不像很多现代工具那样给你直观的报错提示——数据不对跑出来的结果就是错的或者干脆静默失败。我当初第一次跑的时候就因为坐标文件里多了一个表头软件直接给了一堆莫名其妙的结果。2.1 四类核心输入文件的作用与格式示例SaTScan 的纯空间和时空分析通常需要准备以下四类文本文件全部可以用普通记事本编辑扩展名没有硬性要求但内容格式必须规范。病例文件Case File是最核心的输入。它记录每个地理位置上的发病数最简单的格式就是两列第一列是位置 ID第二列是该位置的病例数。位置 ID 可以是一串编号也可以是乡镇名称但必须和坐标文件、人口文件里的 ID 完全一致。如果你做的是时空分析那在病例文件里还需要加一列时间信息告诉软件病例发生在哪个时间段。人口文件Population File记录每个位置的人口基数格式和病例文件类似第一列位置 ID第二列人口数。这里有一个容易懵的地方如果你用的是 Bernoulli 模型病例-对照那第二个文件就不叫人口文件而是控制文件Control File内容换成对照组的人数。坐标文件Coordinates File是所有文件里格式最敏感的。三列内容分别是位置 ID、X 坐标、Y 坐标。坐标可以是经纬度注意 WGS84 坐标系经度在前纬度在后也可以是投影坐标比如 UTM 的东向和北向值。SaTScan 对坐标的解读方式取决于计算距离的选项这块我在下一节要重点展开。除了以上三个必选文件还有一个地理信息文件Geo File属于可选项当你使用椭圆扫描或者需要限制扫描窗口不能跨过某些区域边界时才需要。它描述的是地图边界线上各个节点的连接关系格式相对复杂对大多数第一次用 SaTScan 的人来说可以先不管它。我以一份县级食源性疾病监测数据为例给你看一个最简化的病例文件大概长什么样PositionID Cases 101 12 102 5 103 18 104 3注意不要在这个文件里写“例”、“人”这种中文单位也不要用 Excel 保存后残留的空格或特殊符号列与列之间用空格或 Tab 分隔都可以但全文件必须统一。2.2 坐标系选择与 ID 匹配内藏的暗坑关于坐标文件我踩过两个比较典型的坑在这里有必要多说几句。第一个坑是经纬度顺序搞反。SaTScan 在默认状态下要求坐标文件中的第一列坐标是 X对应经度或投影坐标的东向第二列是 Y。如果你从 ArcGIS 或 QGIS 导出数据时习惯性地先导出了纬度再导出经度软件不会报错但所有距离计算都会变成在错误的坐标系里进行聚集位置会偏到完全荒谬的地方。我的习惯是在坐标文件里用表头明确标注 Longitude、Latitude跑完分析后第一时间把输出结果里的聚集中心点拿到地图上核对一下看是不是落在真实的市区或乡镇所在位置。验证数据永远比信任软件选项靠谱。第二个坑是位置 ID 的口径不一致。比如病例文件里用的是乡镇界编码“110101”但坐标文件里用的是乡镇名称拼音“huangcun”或者一个文件里是文本型数字另一个文件里是数值型数字SaTScan 不会提示你哪些 ID 对不上只会把对不上的位置当成数据缺失处理。我建议所有文件统一从同一个行政区划代码字典中取 ID而且不要用 Excel 默认格式打开编辑后另存为——Excel 对超长数字 ID 的自动转换很容易把事情搞坏。优先用记事本、VS Code、或者写个小脚本统一生成这三个文件。2.3 用案例数据构建一份可直接运行的输入文件集理论讲了不少我给你一个可以直接照着操作的例子。假设我们要分析某县 12 个乡镇一年的乙肝病例分布数据结构包括乡镇编码、乡镇名称、经度、纬度、年病例数、年人口数。第一步新建一个文本文件命名 case.txt内容如下ID cases 1 23 2 17 3 44 4 6 5 12 6 8 7 25 8 19 9 31 10 7 11 14 12 9第二步新建 pop.txt内容如下ID population 1 33520 2 28760 3 45210 4 21980 5 30750 6 26430 7 38560 8 34620 9 41280 10 25670 11 29340 12 22150第三步新建 coord.txt里面填乡镇所在地的经纬度格式为ID longitude latitude 1 111.231 32.452 2 111.442 32.218 ...把三个文件放到同一个目录下打开 SaTScan在 Input 标签页里分别选入这三个文件然后在 Analysis 标签页选择纯空间 Poisson 概率模型直接跑。跑完你就会看到一份完整的结果报表。至于模型选择为什么选 Poisson、参数怎么调是接下来要重点讲的内容。3. 参数配置的核心逻辑模型选择与窗口上限为何如此关键很多人第一次打开 SaTScan 的参数面板会被吓一跳选项确实多。但真正影响分析结果的其实就三类概率模型、窗口大小上限、蒙特卡罗设置。把这三类搞明白软件的其他配置基本都不会跑偏。3.1 Poisson、Bernoulli 与时空置换模型怎么选SaTScan 提供了好几种概率模型日常用得最多的是三种。第一种是 Poisson 模型适用于病例数加人口基数的数据比如我前面例子里的“乡镇乙肝发病数乡镇人口数”。它的默认假设是每个人发病概率相等期望病例数跟人口规模成正比实际观测数显著超过这个期望值的地方就是可疑聚集。绝大多数传染病空间聚集研究、肿瘤发病地理分析用的都是这个模型。第二种是 Bernoulli 模型适用于病例-对照类型的数据比如你收集了一批发病病例和一批健康对照每个个体都有明确的位置属性。这时候不需要人口数只需要病例文件和控制文件SaTScan 会判断哪些位置上病例的比例异常升高。动物疫病、某些病例对照研究的空间聚类分析常用到它。第三种是时空置换模型Space-Time Permutation这个模型比较特殊它不需要人口文件只需要病例文件而且病例文件里必须包含时间信息。它通过随机置换病例的时间和空间标签来生成期望值特别适合没有准确分母数据、但正在发生的暴发监测场景比如基于急诊就诊记录做肠道传染病异常信号探测。我的建议很直接如果你有人口或对照组数据优先用 Poisson 或 Bernoulli因为检验效力更可靠如果确实拿不到可靠分母只有病例个案数据和发病时间那就用时空置换模型兜底总比不做分析强。3.2 空间窗口与时间窗口的取值策略从 50% 的默认值往外退窗口上限是 SaTScan 里最容易被随意设置、又对结果影响巨大的一个参数。默认值是总人口的 50%意思是扫描窗口的最大半径不超过覆盖区域总人口的 50%。这个上限设得越大计算时间越长而且可能扫出一些横跨半个研究区的巨型聚集这类聚集往往没有实际公共卫生意义。我在疾控项目里比较推荐的策略是如果你对这个疾病的聚集规模一无所知先按总人口的 50% 跑一遍看看大尺度格局如果已经知道疾病大概集中在村级尺度或者乡镇尺度比如怀疑水源污染导致某几个行政村高发那么把上限压到 10% 到 30% 之间更合适。窗口上限的本质是你对最大聚集范围的事先约束——它不是一个需要到处套用的魔法数字而是一个需要结合传染病传播特点、行政边界、人口密度来认真思考的假设。如果做时空分析还要面对时间窗口上限的设置。默认也是总研究期的 50%。比如你手上有 2018 年到 2023 年总共六年的月度数数据默认允许聚集持续三年。这样处理的好处是能发现持续性聚集但对急性暴发可能造成过度平滑。我建议对急性传染病直接设成 30% 或更短的时间窗口重点捕捉短期集中暴发而慢性病分析则可以保留较长时间窗口。3.3 蒙特卡罗重复次数、随机种子与计算成本的权衡蒙特卡罗模拟次数直接影响 P 值的稳定性和计算耗时。SaTScan 的默认值是 999 次这个数值在绝大多数情况下够用。但如果你要发表论文或者结果要用于监管决策我通常会调到 9999 次。原因在于999 次模拟对 P 值的最小分辨率是 0.001而 9999 次可以分辨到 0.0001在多重比较校正和敏感性分析时这个精度差别有意义。代价是计算时间几乎线性增加。跑一次 999 次模拟可能需要几分钟调到 9999 次可能就要几十分钟。我的操作习惯是先用 999 次快速摸一遍结果确认数据没有低级错误、参数没有明显问题之后再切到 9999 次跑最终版。随机种子在 SaTScan 里默认是随机的如果你想让自己某次分析的结果可以精确复现需要手动固定随机种子这个选项在高级设置里论文的复现性要求高的时候就靠它了。4. 结果文件的正确打开姿势从 P 值到相对风险SaTScan 跑完以后会在你指定的输出目录里生成一堆文件。新手最常犯的一个错误是只盯着结果文本文件的最后几行看忽略了文件里包含的丰富信息。这一节我带你逐个拆解。4.1 输出文件清单与各自用途SaTScan 输出的核心文件是结果文本文件通常是 .txt 或 .html 格式里面包含了所有检测到的聚集簇信息每个聚集簇的编号、地理坐标或位置列表、半径、观测病例数、期望病例数、相对风险Relative RiskRR、对数似然比LLR、P 值。这是你写论文、写分析报告时最主要的信息来源。除了文本结果软件还会输出一系列地图相关文件。如果你在输出选项里勾选了地图文件输出会得到 .map、.gis、.col、.geo 等格式的文件这些是给 ArcGIS、QGIS 用的几何信息。另一个很有用的输出是 KML/KMZ 文件直接用 Google Earth 打开就能看到聚集圈的叠加显示。时空分析还会额外输出 .gif 动画文件它按时间顺序动态展示聚集圈的变化过程做汇报PPT时效果很好。4.2 读懂结果表LLR、RR、P 值的实际含义打开结果文本文件你会看到一张类似这样的表Clusters 1.Location: 4,8,11 Coordinates: (111.67,32.36) Radius: 12.3 km Population: 88740 Number of cases: 126 Expected cases: 78.5 Relative risk: 1.61 Log likelihood ratio: 15.32 P-value: 0.001我们逐个看。坐标和半径描述了这个聚集圆圈的地理位置和覆盖范围Population 是圈内的总人口Number of cases 是圈内的实际病例数Expected cases 是如果全区域发病风险一致、按人口比例算圈内应该出现的病例数。相对风险 RR 等于实际病例数与期望病例数的比值1.61 意味着这个区域内发病风险是全区平均水平的 1.61 倍。LLR 是对数似然比数值越大代表聚集证据越强。这里最容易误读的是 P 值。SaTScan 报告的 P 值不是针对单个聚集簇的常规显著性检验而是基于蒙特卡罗模拟的经验 P 值。它的含义是在全区域随机分布的前提下出现一个比当前聚集簇 LLR 更大的聚集簇的概率。因此P 值越小说明这个聚集越不可能是随机形成的。但要注意P 值和 RR 不完全等价——一个 RR 很高的小聚集可能因为病例数太少而 P 值不显著反之一个 RR 只有 1.2 的大聚集因为覆盖人口多P 值可能非常显著。解读时两个指标都要看。4.3 多个聚集簇同时出现时先看哪个数据分析中经常出现一次扫描出十多个聚集簇的情况。SaTScan 会自动把最不可能由随机产生的那个标记为主聚集簇Primary Cluster它拥有最大的 LLR。其余那些在空间上与主聚集簇有重叠的簇会标记为次要簇Secondary Cluster其中有一部分可能是主聚集簇的“附属”并不是独立的新聚集。我的处理习惯是先把主聚集簇和标记为“没有与任何更强的簇重叠”的次要簇找出来因为这些是真正独立的聚集信号然后再看那些与主聚集簇重叠的簇思考它代表的是同一个暴发事件的延伸还是确实存在两个独立风险源。这一步需要结合流行病学专业知识判断软件本身不会替你下这个结论。5. 让分析结果走出去地图可视化的导出链路跑完分析只是第一步把聚集结果可视化地呈现在地图上才是向领导、合作单位或审稿人展示的关键。SaTScan 的可视化能力不如专业 GIS 软件但它的导出链路设计得还算顺滑。5.1 从 SaTScan 到 Google Earth / QGIS / ArcGIS 的操作路径最简单粗暴又实用的方式是直接在 SaTScan 输出设置里勾选生成 KML 文件然后用 Google Earth 打开。KML 里会显示每个聚集圈的边界并且可以直接看到聚集圈半径、人口、病例数等属性信息。Google Earth 的好处是对新手零门槛交互流畅适合快速查验聚集位置是否合理。需要做更正式的地图出图时我一般走两条路。一是把 SaTScan 输出的 .gis / .map / .col 文件组合导入 ArcGIS 里叠加到底图上这在老牌 GIS 用户里很常见。二是在 QGIS 里操作QGIS 可以通过插件或者直接读取 SaTScan 的文本结果文件按坐标生成点或圆形的矢量图层。我个人现在偏向 QGIS 路线因为免费、跨平台而且对 SaTScan 结果文件的兼容性反而更好不用解析那些老旧的 .gis 格式。还有一个土办法也值得提一下如果只是做简单汇报可以直接复制结果文本中的聚集圆心坐标和半径自己在绘图工具里画圆标注。数据量小的时候这个方法最省事。5.2 时空分析结果的时间动画解读如果你做的是时空扫描SaTScan 生成的时间序列 GIF 图非常直观地图底色不变聚集圈随着时间推进一个一个出现或增大能让人一眼看出暴发是从哪个区域开始、朝哪个方向扩展的。严格说这个 GIF 反映的只是扫描窗口在时间维度的移动情况解读时需要结合结果表里的时间范围列来看不能只凭动画脑补出“传播路径”。用 GIF 做汇报时我建议配合一个时间线表格聚集簇编号、起止时间、持续时间、涉及的地区列表。这样观众既能看到动态画面又能定位到具体的统计数据和风险等级。5.3 地图坐标对齐中的常见细节地图叠加时最花时间的环节其实是坐标对齐。如果 SaTScan 用的是经纬度 WGS84而你的底图是某地方坐标系或 GCJ02 坐标别问为什么会有这种情况国内做地图对接时很常见直接叠加会出现几百米的偏移。此时需要在导入 GIS 前先做坐标转换或者确保输入 SaTScan 的坐标文件和底图用的已经是对应的投影坐标系。动手之前先在底图上验一个已知地标点再批量处理能省掉大量返工时间。6. 实操中反复踩到的坑和我的应对习惯最后这部分是我最想写的因为软件的官方文档不会提醒你这些但它们每一个都真实地消耗过我的时间。6.1 病例为 0 的区域怎么处理做空间扫描时大量乡镇的病例数是 0。这是常态尤其对于低发病率疾病。Poisson 模型可以处理零病例区域不需要把它们剔除也不需要做平滑。但要注意的是人口文件里对应位置如果填了 0那就会出问题——人口为 0 但病例也为 0 的地方相当于分母为零的未定义区域软件计算时会绕过它但如果人口为 0 的地方又有病例那数据本身就不合理需要回到源头核查。6.2 最大空间窗口上限对结果影响究竟有多大我做一个模拟数据实验给你看同样是 100 个位置点、200 个病例把最大空间窗口从 50% 改成 15% 后主聚集簇从覆盖 27 个位置的大圆变成了覆盖 8 个位置的小圆中心位置也有轻微偏移。两个结果都“统计显著”但公共卫生解释完全不同。这一点一定要牢记SaTScan 输出的聚集位置和范围在很大程度上受到你预设的窗口上限影响。所以我每次分析都会至少跑三组参数50%、30%、15%然后对比主聚集簇的稳定性。如果三组结果的核心位置一致只是范围大小不同那结论就比较可靠如果三组结果完全跳变我就要反思数据质量或模型选择有没有问题。6.3 分母数据的时间口径问题做时空扫描时人口文件一般只有一个数值但病例数据跨越多年。如果期间区域人口结构发生了明显变化比如某地因为城镇化人口大幅增长而病例总数也增长那用静态人口数当分母会扭曲聚集检测结果。解决思路是尽量获取每年的年龄别人口数把数据分层处理或者至少在校验时看看聚集检测结果是否对人口数据版本敏感。SaTScan 的高级功能可以处理分层数据但这部分配置复杂我建议一般用户先从静态人口跑通流程再把敏感性分析作为论文的补充内容。6.4 关于结果解读和论文写作的提醒SaTScan 给的是统计聚集不是因果结论。哪怕一个聚集簇的 P 值小于 0.001也只能说明这个地方的发病数在空间上偏离了随机分布至于为什么偏离——是环境污染、人口老龄化、还是诊断偏倚——需要靠后续的流行病学调查去验证。我在实际报告里都会明确写“本分析识别出需要优先关注的空间区域”而不是“本分析证实了某地存在某风险源”。这个表述上的分寸既是学术严谨的要求也是避免把统计结果过度现实化的基本操守。6.5 数据处理的小习惯备份项目文件与记录参数SaTScan 支持把当前分析配置保存成项目文件.ssi 格式但我见过太多人跑完就走从不保存配置。过了一个月要复现某个结果时完全想不起来当时用的窗口上限是 20% 还是 30%蒙特卡罗次数是 999 还是 9999。我的习惯是每次分析都复制一份项目文件命名带上日期和参数简述比如“20240115_hepB_poisson_win30pct_9999rep.ssi”同时把输入文本文件也一并归档。这样以后不管是自己复查还是回应审稿人要求都能在五分钟内找回原始分析条件。最后再说一个听起来很小、但确实坑过我的细节SaTScan 对中文路径的处理不算特别稳定。如果输入文件所在的文件夹路径里有中文或者空格某些版本的软件在导出地图文件时会出现失败或生成空文件的情况。我现在所有分析项目一律使用纯英文路径这个习惯从我第一次遇到导出失败就一直保持到了今天。虽然 SaTScan 的界面操作已经足够直观但把数据文件和项目路径都规范好整个分析流程才能真正顺畅地跑起来。以上就是我基于实际使用经验整理的全部内容希望对准备上手或在数据里挣扎的同行有一点实际的帮助。