阿坝州大骨节病区病情动态仿真

四川省自然资源科学研究院生物生态所(610041) 杜 川

为了启动“阿坝州扶贫开发和综合防治大骨节病试点工作”(以下简称“试点工作”)，州政府于2007年在大骨节病区开展了前所未有的大规模入户详查，收集到占全州人口半数的44.3万余人的病情相关基本信息数据，囊括大骨节病患者38322人，占全州总数的93%。经上报四川省和国务院，于2008年确认后，启动“试点工作”。经统计，截止2007年底，州内有11个县的95个乡、镇流行大骨节病，379个村被列入病区村[1]，占行政村总数的28%，患者总数41184人。

利用这批数据来探讨大骨节病区未来病情和综合防治策略，既可提升这批数据的价值，又有利于总结过去、预见未来。

KBD病情发展的随机性

1.KBD病情转变概率

大骨节病是一种地方性、多发性、变形性骨关节病，医学界称本病为Kaschin-Beck Disease(以下简称KBD)[2]。因KBD病区存在众多因素，使得一个人患不患KBD或者病人病情发展情况都具有不确定性。将前述病区数据所记录的人群看作样本空间Ω[3]，根据KBD病情与年龄的分类统计(表1)，定义二维随机向量(x,y),其中，x∈X={0，1，2，3，4},代表KBD病情，y∈Y={1,2,…,112},代表年龄,

表1 KBD病区按病情和年龄分类统计的人数

(1)

这里，“正常人”仅指非KBD患者。显然，x、y相互独立。x、y的边际分布律见表2。

表中，Px(x=0)>Px(x=2)>Px(x=3)>Px(x=4)>Px(x=1),正常人在KBD病区占大多数，之后依次是Ⅰ、Ⅱ和Ⅲ度病人，而早期病人最少。

表2 随机向量(x、y)的分布

KBD病区内个体病情发展一般按0-1-2-3-4的顺序进行，即病情较重的人都来自病情较轻的人。记KBD病情由k转变为k+1的概率为Pc(k,k+1)，根据条件概率定义和性质，有

(2)

从表2和(2)式的计算结果(表3)看，虽然正常人患KBD的可能性并不大，但早期病人转变成Ⅰ度病人的可能性较大。其中，Pc(0,1)=1.41%，与某些专家的抽样调查结果接近，但也有一些地区的抽样调查显示1.6%或更高[4-6]。

表3 KBD病情转变概率

2.KBD病情转变延迟

KBD病人的病情发展会经历多年。对k∈X且k≠4，记KBD病区人群的病情从k转变为k+1所需年数为Dk，称Dk为“KBD病情转变延迟”。

对k∈X且k≠0，设A(k)为病情为k的人数达峰值时的年龄，求得：A(1)=20,A(2)=43,A(3)=64,A(4)=79。又设a为早期发病年龄，虽然调查漏掉了这项，但由KBD发病规律知1≤a≤15。故不难证明

(3)

根据(3)式可以估算Dk(k=1,2,3,4)所在区间(表4)，例如20-15≤D1≤20-1，即5≤D1≤19。

表4 KDB病情转变延迟

综合上述，D1平均12年、max{P(x,y=1)}=P(x=2,y=1)=1.01%、max {Pc(k,k+1)}=Pc(1,2)=82.32%，表明早期病人少的原因除了早期发病集中于15岁以内人群[4-5]且阿坝州儿童中的患病率较低有关外[6]，还与当地自然条件恶劣、生活艰苦，使得早期KBD病人病情较不稳定而容易转变成Ⅰ度的实际情况有密切关系。所以，应实施“异地育人”，将病区儿童尽早转移到条件较好的非病区上学、生活。

可以认为Dk具有均匀分布。在动态仿真模型中Pc(k,k+1)、Dk都用作可调参数。

KBD病区病情动态仿真模型

1.系统动力学建模

系统动力学[7](System Dynamics,简称SD)出现于1956年，创始人为美国麻省理工学院(MIT)的福瑞斯特(J.W.Forrester)教授。SD是系统科学理论与计算机仿真紧密结合、研究系统反馈结构与行为的一门科学，是系统科学与管理科学的一个重要分支。动态仿真指以获取对象随时间变化的行为特征信息为目的，在电脑上运行对象的数学模型。SD建模的重点是用代表各类变量的“累积”、“速率”和“参数”三要素，以及刻画要素相互关系的有向弧来绘制SD流图。以SD为理论基础催生的一批类似本文所使用的Vensim这样的计算机建模系统[8-9]，既能将数学的复杂性隐蔽在SD流图内，又能依流图自动生成动态仿真程序，带来极大的方便。

有了KBD病情转变概率和延迟，可以为SD建模做这样的假设，即KBD病情的发展是分4个阶段有序进行的。本文运用SD所建立的动态仿真模型由一个四阶段KBD病情子系统、一个两阶段人口子系统[10-12]和历史检验模块组成，有6个累积、12个速率和16个参数以及其他变量。“早期病人增”速率使人口子系统产生的KBD病人流入病情子系统。“综合措施效果”等系来自其他子模型的参数(图1)。

表1给出6个累积的初始值，而两个人口累积初始值应扣除病人数。参数中，人口出生率和死亡率均取自阿坝州统计年鉴[13-16]，表3、4 给出“KBD病情转变概率”、“KBD病情转变延迟”，Ⅲ度病人死亡率应略大于正常人口死亡率0.0055，其他死亡率都预设为10-4。所有参数都有4个数值：参照、最大、最小和步长，供动态调整用。累积和速率的表达式也要逐一设置。

2.模型有效性检验

SD模型的有效性检验是为了验证模型运行所产生信息、行为与实际对象的特征、变化规律之间能否很好的吻合，能否正确理解所要考虑的问题。主要检验方法有：灵敏度分析、历史检验。

参数灵敏度检验模型行为模式是否会因某些参数的微小变动而产生与实际情况迥异的现象。Vensim内置的蒙特卡洛[17-18]灵敏度分析表明，在KBD病情转变概率和延迟共8个参数的随机扰动下，KBD病人数的置信范围仅随时间平稳扩大，这与KBD病区病人数相对稳定的实际情况相符，因而通过了灵敏度检验。

历史检验就是验证历史数据与模型输出数据的拟合程度。根据试点工作总结，2007-2012年KBD病区病情稳定且无新增病人，因此2008年所固定的病人数据(表1)就是历史数据，存入变量01-04，其与仿真所输出的对应5年各度病人数之间的相对误差存于变量d1-d4；0t为历史病人总数，dt为其相对误差。仿真输出的总体平均相对误差最大值为2.75%(表5)，可见5年历史检验的总体拟合度应在97%以上。

3.动态仿真结果

SD建模搭建起动态仿真平台后，还要针对问题选择部分参数组成情景模式，只需在Vensim窗体上调整变量数值，就可以切换情景，动态改变输出图表，预测不同情境的KBD病情趋势，实现“牵一发而动全身”且无实物消耗的策略实验，寻求解决问题方案。其中，“综合措施力度”变量可以控制所有KBD病情转变和延迟参数。仿真年限在1～200之间可调。对4种代表性情景(表6)做50年仿真，其结果简述如下。

情景1：在不综合治理但执行阿坝州现行计划生育政策的情况下，病区病人总数不会减少，病人占总人口的比例呈对数曲线，逐渐趋稳，即病区人群将持续贫病交加(图2)。

情景2：不综合治理但强调计划生育的情况下，病人总数与人口数成正比，但病人比例最终将在8%上下变化。可见单靠计划生育政策是不能扭转病区面貌的。

情景3：进行综合治理并杜绝新增早期病人的前提下，50年内，病区每年新增早期发病人将由566人下降到0人，早期病人数由现有645人下降到23人，病人总数由现有38322人下降到36350人，病人比例由现8.64%下降到6.58%，。这意味着农村劳动力恢复，政府供养负担减轻，病区贫病交加的局面得以扭转(图3)。

图1 运行中的阿坝州KBD病情动态仿真模型的SD流图

表5 模型历史检验相对误差

情景4：综合治理但不能杜绝新增早期病人的情况下，如果能较好地解决更换粮食和饮水安全等问题，病人病情将得到缓解，Ⅱ、Ⅲ度病人数趋于减少，病人数和比例下降不明显。

表6 情景模式与情景表

讨 论

“试点工作”自2008年5月6日启动到2012年结束，一直受到中央和四川省的高度重视，从九个方面采取综合措施(图1)开展治理，调动了大量资源，五年累计投入资金55亿元。虽然病区群众的生产、生活和医疗卫生条件都得到很大改善，对缓解病情、斩断病根作用明显，但要继续巩固提高、防止回潮，仍有大量工作。

统计学方法 采用统计学软件SPSS 20.0处理数据，计数资料以百分率 （%）表示，计量资料以均数 ±标准差 （x±s）表示。

通过多情景多方案动态仿真所生成的预测数据，可以为当地政府决策提供大量有关KBD病区病情未来趋势的信息，在KBD防治机制与政策研究方面也是有意义的。可以将情景3看作最有效的策略，这有利于扭转病区贫病交加的恶性循环又有利于减少政府用于供养和治疗的开支。因此，继续做好规划，优化推进模式，将有助于病区以最小代价换取最好效果。

图3 情景3的输出图表

1.GB16397-1996,大骨节病防制效果判定.北京:中国标准出版社,1996.

2.孙殿军.地方病学.北京:人民卫生出版社,2011.

3.盛骤,谢式千,潘承毅.概率论与数理统计(第四版).北京:高等教育出版社,2008:6.

4.黄慧,李富忠,邓佳云,等.2005-2009年四川省大骨节病病情调查结果分析.预防医学情报杂志,2010,26(12):966-968.

5.赵生成,白玛次旺,扎西桑珠,等.西藏藏族自治区大骨节病病情调查结果分析.中国地方病学杂志,2011,30(5):524-526.

6.樊效鸿,王鑫灵,黄勇,等.四川阿坝地区儿童大骨节病病情现状及影响因素调查分析.医药前沿,2012,2(2):377-378.

8.钟永光,贾晓菁,李旭.系统动力学.北京:科学出版社,2009.

9.Lambert M,Surhone,Miriam T.Timpledon,Susan F.Marseken.Vensim.VDM Publishing House,2010.

10.宋健,田雪原,于景元.人口预测和人口控制.北京:人民出版社,1982.

11.陈强.人口系统模型及人口状态分析.南京理工大学,2004.

12.丹尼斯·米都斯.增长的极限.李宝恒,译.长春:吉林人民出版社,1997.

13.四川省阿坝藏族羌族自治州统计局.2009年阿坝统计年鉴.阿坝州:四川省阿坝藏族羌族自治州统计局,2009.

14.四川省阿坝藏族羌族自治州统计局.2010年阿坝统计年鉴.阿坝州:四川省阿坝藏族羌族自治州统计局,2010.

15.四川省阿坝藏族羌族自治州统计局.2011年阿坝统计年鉴.阿坝州:四川省阿坝藏族羌族自治州统计局,2011.

16.四川省阿坝藏族羌族自治州统计局.2012年阿坝统计年鉴.阿坝州:四川省阿坝藏族羌族自治州统计局,2012.

17.杜比.蒙特卡洛方法在系统工程中的应用.西安:西安交通大学出版社,2010.

18.曾平,赵晋芳.非嵌套模型选择的Vuong检验及其在ZIP模型中的应用.中国卫生统计,2011，28(1):8.