1、中国人口增长预测模型的建立与分析摘要针对我国人口发展过程中出现的老龄化进程加快，出生人口性别比持续升高，乡村人口城镇化的新特点，我们基于LESLIE 矩阵，着重考虑城镇与乡村间的人口迁移及女性人口比例变化对我国人口增长的影响，经过两次改进建立了便于计算机求解的差分方程模型，对我国2005年以后45年的人口增长进行了预测。随后利用时间段参数设置法，对差分方程模型又进行了一次改进。然后运用等维灰色系统预测法对该差分方程模型的中短期预测进行了检验，同时根据2001年人口基本数据运用此模型对2001年2005年进行了预测，并用实际数据对预测结果进行了检验。我们将预测区间分为20062020年、2021
2、2035年、20362050年三个区间，以量化短期、中期与长期。通过调整模型中相关参数及输入条件，定量地分析了男女性别比例、老龄化和乡村人口城镇化对我国人口增长的影响。预测结果表明，从短期来看，我国的出生性别比变化不明显，将在短期内维持基本不变，老龄化进程在15年内在上升了8个百分点，人口扶养比持续升高，这将加重我国的人口压力，乡村人口城镇化水平进展缓慢；从中期来看，总人口性别比将保持在1与1.1之间 ，老龄化进程将呈线性增加趋势，乡村人口城镇化水平将持续发展；从长期来看，老龄化进程将在2035到2045年经历老龄人口高峰平台，老龄人口比重在0.3以上，育龄妇女人数持续下降，总人口数将在202
3、3年达到峰值14.05亿。 关键词：LESLIE矩阵，人口预测，性别比例，城镇化，老龄化，灰色系统预测一、 问题的重述人口问题是中国社会发展的重要问题，对中国人口的中长期预测有助于政府制定相应的政策保持中国的长治久安。现需要解决的问题如下： 1.主要根据20012005年的人口统计数据，对中国人口增长的中短期和长期趋势作出预测，特别要关注老龄化，出生人口性别比及乡村人口城镇化等因素。2.指出所建模型的优点和不足之处。二、 模型假设1.在未来50年人口生存的社会环境相对稳定（即没有战争及毁灭性灾难）。2.国际人口迁入与迁出量相等。3.在本世纪中叶前，我国计划生育政策稳定。4.题目所给抽样数据是随
4、机的，真实地反映了整体实际情况。三、 符号说明分别表示乡村、镇、市第年岁人口的死亡率；分别表示乡村、镇、市第年岁的人口数；分别表示乡村、镇、市第年岁的女性生育率；分别表示乡村、镇、市第年岁人口的女性比；分别表示乡村、镇、市第年的婴儿死亡率；分别表示乡村、镇、市第年的出生人数；分别表示乡村、镇、市第年岁女性的生育模式；分别表示乡村、镇、市第年的总和生育率；分别表示乡村、镇、市第年的总人数； 分别表示乡村、镇、市第年岁女性的总人数；分别表示乡村、镇、市第年岁女性的死亡率；分别表示乡村、镇、市第年岁男性的总人数； 分别表示乡村、镇、市第年岁男性的死亡率；表示为迁移人口中女性所占比例；分别表示乡村、镇
5、、市出生人口中女性所占的比例；四、 问题的分析人口发展过程的定量预测，需要预测出未来的人口发展趋势，包括人口总数、人口的性别、年龄和城乡构成，人口出生、死亡和自然增长率的变化以及在未来的人口构成中劳动力和抚养水平及老龄化水平等各项人口指数全部测算出来。人口增长的决定因素为出生率、死亡率和人口基数，但人口分布，人口素质，宏观政策和人口结构（如：年龄结构，性别比例等）等众多因素能够影响出生率与死亡率的波动，从而从根本上影响我国人口的增长。鉴于我国人口问题已有多方面的研究，我们针对近年来我国的人口发展出现的一些新特点，忽略国际人口流动，故可以认为我国人口为一个封闭的系统。对于封闭的系统来说 ，某时刻
6、人口总量=人口基数+新生人口数死亡人口数。为了提供更多关于市、镇、乡的人口增长分布趋势，我们对三者分别进行研究。但由于市、镇、乡人口在我国城镇化进程中迁移流动量大，故必须把人口迁移考虑进去。在本论文中，我们希望对我国人口进行中短期及长期预测。现在给出三个区间，分别表征从2005年算起的短期、中期与长期。建立LESLIE矩阵差分方程模型，采用人口普查的相关数据、运用灰色系统预测方法对模型的合理性进行检验，从而可以预测我国人口发展的中短期和长期趋势。五、 模型建立与求解5.1 基于LESLIE矩阵的模型的建立考虑到市，镇，乡村的人口性别比例，妇女生育率以及人口的死亡率都有所差别，我们分别建立市，镇
7、，乡村的差分方程模型。市，镇，乡村合起来即可得到全国人口增长的差分方程模型。（1）首先，在不考虑人口迁移的情况下（以人口的户籍变动为准），以乡村为例，建立人口增长模型。记人口最大年龄为岁（由所给数据取为90）。为乡村第年岁的人口数（用上标1表示乡村，上标2表示镇，上标3表示市），为乡村第年岁人口的死亡率，由人口死亡率数据为所占该类年龄段人口的千分比（），我们可以得到乡村第年岁人口的死亡率为： 1000 ；于是第年岁人口数为： （1）记为乡村第年岁女性生育率，该生育率表示所有岁女性平均生育婴儿数，生育率为该年龄妇女生育子女与该类年龄妇女的千分比，为育龄区间（根据题目附录二数据取15到49岁为育龄
8、区间），为乡村第年岁人口的女性比率，则乡村第年的出生人口数为 （2）考虑到婴儿并不是全部都能活到年统计时刻（如2005年统计时刻为11月1日零时零分）其中有些婴儿由于疾病等原因死亡，能够活到年统计时刻的婴儿数是，因此， 就是年到年的婴儿死亡数，记婴儿死亡率为， 则： 1000 于是第年的婴儿数 （3）由（1）、（2）、（3）式可得在年的一岁人口数： （4）将分解为 （5）其中：是生育模式，用以调整育龄妇女在不同年龄时的生育率高低，且 ：；是乡村妇女的总和生育率，。 记 则由（4）、（5）式可得 (6)于是我们可以得到乡村人口差分方程模型为： 将其表示为矩阵形式 ，记 此处，即为LESLIE矩阵
9、；同时，记 则由上述方程得：运用同样方法可分别得出镇与城的人口差分方程模型为：其中，与，的推导方法相同。（2）由于市，镇，乡村之间并不是相互独立的，他们之间有着频繁的人口流动，在实际问题中不能被忽视，下面我们在考虑人口迁移的情况下对上述模型进行改进。考虑到在实际发生的人口迁移中多数由于贫富差距引起，我们在对模型进行改进时仅考虑乡村、镇、市的人口净迁移人口量，可以看到镇、市人口净迁移量都为正。我们假设每年乡村到城镇的人口迁移数为上年人口总数的倍，注意到一个地区人口数量与经济发达程度有很大联系，我们以市，镇总人口的比例来分配乡村到其人口迁移的数量。由于我们以每年的各个年龄段为预测变量，必须考虑各个
10、年龄段的迁移数量，为简化起见我们以2001至2005年各个年龄段五年人口数量之和占五年总人口数量的比例分配各个年龄段的迁移数量，并记为年龄所占比例。记 ，分别为乡村，镇，市在第年的总人口数，显然满足： ，由以上分析，第年乡村1岁人口数量为：岁人口的数量应为： 第年镇的1岁人口数量为： 岁人口的数量应为： 第年市的1岁人口数量应为： 岁的人口数量应为： 于是我们可以得到在考虑人口迁移的情况下市，镇，乡村的差分方程模型为： （3）从理论上，男女出生比例应该为1：1，但是根据题目附录二材料，我国男女出生比例的实际情况是有差异的，生男孩比例普遍比生女孩高。我们应该对每一个孕龄妇女所生子女的性别加以区分
11、，而不能笼统地概括为1：1，这样可以更加精确地得到每个年龄男女比例的预测值。同时，由于城镇化进程中，乡村人口的迁移带来了城镇人口男女比例的变化。在这个过程中，我们也应该着重对性别因素进行考虑。以下，是对前面模型的修正，以期使模型对性别比例变化产生足够的灵敏度。出生人口中女性所占比例，对人口的增长起着至关重要的作用6，有必要对其数量进行单独预测。为此，我们在（2）模型的基础上分别设乡村第年龄女性人口数量为，男性人口数量为，镇第年龄女性人口数量为，男性人口数量为，市第年龄女性人口数量为，男性人口数量为，在迁移人口数量中女性所占比例为。根据所给乡村，镇，市男女出生比例数据，进行数据的预处理，得到乡村
12、，镇，市出生人口中女性所占比例，记为。那么乡村，镇，市女性第年0岁人口数量为：相应的男性第年0岁人口数量为：第年1到龄女性以及男性人口的预测只须分别用，替代（2）模型中的 ，并在迁入或迁出的人数前分别乘以迁移率或，以乡村女性第年龄人口数预测为例，有：其他预测以此类推。5.2模型的分析与求解5.2.1模型参数的设定我们根据所建立模型，需要以下几个输入量：总和生育率，婴儿死亡率，生育模型，每年乡村到城镇的人口迁移数为上年人口总数的倍数，乡村到城镇的人口迁移中女性所占比例，各年龄人口的死亡率，出生人口中女性所占比例，总人口数；根据中国统计年鉴2006，我们得到1987年到2005年中每年的总人口数。
13、可以作为模型预测的初值。在题目附录二中，给出了市、镇、乡男女出生比例。通过对其时间序列分析，得出在很短时间内，都稳定到某一值：0.4552732020.4588594280.472078522本论文只关注从2005起，未来45年我国人口的发展趋势。由于社会环境稳定，可以简单地认为未来45年中，各年龄人口的死亡率保持一定。国家政策可以影响5乡村到城镇的人口迁移中女性所占比例及迁移总量的倍率。在以下求解中，我们可以分析不同值下人口增长曲线的变化。同时，国家计划生育政策也将影响总和生育率，从题目附录一中可知，为了保持社会经济的稳定发展，应该把总和生育率控制在1.8左右。由于乡村与城镇总和生育率的不同
14、,我们将各生育率定为如下：2.01.81.6根据2006中国卫生统计年鉴上1991年至2004年婴儿死亡率的变化7，即中国儿童发展纲要（20012010年）提出的每10年降低20个百分点8。经过数据处理，得出,2050年婴儿死亡率10.48.87.5则可以得出，在未来45年里，乡、镇、市总和生育率为17.4,12.9,8.8。我们首先对模型（3）运用Matlab编程对其参数进行分析，由上所定参数，并且取 0.5。此时仅有未知，为分析其对总人口数变化规律的影响，分别令其取0.01,0.001,0.0001。以2005年数据为预测起始数据得到图1：图1 未来45年总人口预测从图中可以看出对总人口数
15、的影响很小，在前期1到15年几乎没有影响。那么可以认为。对于，在确定的情况下分别令其取0.2，0.5，0.8得到图2：图2 未来45年总人口预测从图中可以看出对人口高峰的到来时间有一定影响，但是影响不够显著。所以可以简化认为。以上述给定参数，对我国总人口数从2005年预测到2050年的变化图3（其中，）：图3 未来45年总人口预测从预测数据可以得出总人口在2023年达到峰值14.045亿。市，镇，乡村以及总人口的男女比例从2005年到2050年的变化图4：图4市，镇，乡村以及总人口的男女比例从图中可以看出总人口，市，镇男女比例一直保持在1附近，只有乡村男女比例在2035年后呈指数形增长，但在2
16、050年仍小于1.8。总人口的抚养比从2005年到2050年的变化图5：图5 总人口的抚养比从图中可以看出2005至2010总人口的抚养比一直减小，2010年达最小值：0.38.2010年以后抚养比一直增加，2043年左右达最大值：0.64。人口城镇化比例从2005到2050年的变化图如6：图6 人口城镇化比例总人口考虑60岁以上人口的老龄化进程从2005年到2050如图7：图7 60岁以上人口的老龄化从图中可以看出我国在2035年到2045年形成老龄人口高峰平台，60岁以上老年人口比重在0.3以上。总人口中育龄妇女人数从2005年到2050年变化如图8：图8 总人口中育龄妇女人数从图中可以看
17、出，在2010年后育龄妇女人数一直下降，2045年降到2.4亿以下。5.3对差分方程的修正影响人口增长因素很多，如农村人口向城镇的迁移、婴儿死亡率、总和生育率等。在对长时间的预测中，我们只能假定一些因素是不变的。但是，我们应该清楚地看到，随着时间的推移、大环境的变化必然带来因素的改变。为了解决这个问题，我们需要不断地对模型中一些参数进行修正。在这里，我们采用时间段参数设置法9，将时间区间设置等间距的15年。在区间内，我们认为参数是不变的。在一个区间数据预测完成后，我们在上一区间预测数据的基础上修正因素，进行下一区间预测，以达到最大限度补充信息量，保证预测效果的显著。下面仅在2020年及2035
18、年修正婴儿死亡率，得到以下预测见表1（人口单位/10亿）表1 差分方程修正年份2006200720082009201020112012201320142015一般预测1.31371.32091.32791.33481.34191.34911.35661.36411.37131.378修正预测1.31371.32091.32791.33481.34191.34911.35661.36411.37131.378年份2016201720182019202020212022202320242025一般预测1.3841.38941.39381.39741.40031.4021.40311.40311.4
19、0211.4004修正预测1.3841.38941.39381.39741.40031.40221.40331.40341.40251.4008年份2026202720282029203020312032203320342035一般预测1.39781.39471.39111.38711.38281.37811.37331.36841.36341.3582修正预测1.39831.39531.39171.38771.38351.37881.37411.36921.36421.3591年份2036203720382039204020412042204320442045一般预测1.35281.3472
20、1.34131.33521.32851.32151.31421.30611.29771.2888修正预测1.35381.34841.34261.33661.33011.32321.3161.30811.29991.2911年份20462047204820492050一般预测1.27951.271.25991.24981.2399修正预测1.28191.27261.26261.25261.2428从表中可以看出，参数设置法对人口增长产生了一定影响。如果把各个影响人口增长的因素加到对分段参数设置法修正中，会得到更普配实际的预测。但由于本论文成文仓促，本修正方法没有做完，是此篇论文的一个不足之处。六
21、、 模型的比较与检验为了对差分方程模型进行检验与比较，我们将模型的时间节点朝前移五年，用2001年数据为输入数据，得出2002年到2005年的预测数据，经过同实际统计数据的对比，我们可以看到差分方程模型的优点与不足。另外，对模型进行检验，我们可以比照国家人口发展战略研究报告(2007)的相关数据指标。但是为了更加独立地校验模型，我们另建立了灰色预测模型，以期看到差分方程模型与灰色预测结果的匹配程度。6.1 以普查统计数据检验差分方程模型带入2001年数据，由差分方程模型（3）得出预测结果见表2。表22002年2003年 2004年2005年实际值（万人）12845312922712998813
22、0756预测值（万人）128000128980129890130760绝对误差（万人）45324798-4相对误差0.0035270.0019110.0007540.000031通过实际人口量与预测人口数量的比较，我们可以看到，差分方程模型在总人口的预测上与实际情况符合得很好。6.2灰色预测模型的建立与求解6.2.1灰色数列预测原理灰色数列预测最大的特点是单数据序列的预测。在形式上，只运用预测对象自身的时间序列，而与预测对象相关联的其他因素没有参与运算和建模，灰色系统理论把影响客观系统的诸多因素及它们之间的关系定为灰色量。对这样的灰色量进行预测，就可从自身的时间序列中寻找有用信息，发现内在规律
23、，建立模型进行预测。6.2.2建立模型(1)数据的检验与处理对于给定的原始时间19962005年数据列： 数列的级比 如果所有的级比都落在可容覆盖内，则数列可以作为模型和进行数据灰色预测。经过计算，所以全部落入中。数列可以作为模型和进行数据灰色预测。(2)确定数据序列对做AGO生成，有=AGO，其中 ，则 =122389,246015,370776,496562,623305,750932,87938，1008612,1138600,1269356求均值数列 184202，308395.5，433669，559933.5，687118.5，815158.5，943998.5，1073606，1
24、203978于是建立模型，灰微分方程为10： 相应的白化微分方程为： 记 ，则由最小二乘法，求使得达到最小值的。于是求解方程得 由上可以得到白化方程的时间相应式 而且 得还原方程： (3)对灰色预测模型的改进一般情况下，)模型通过对数据序列长度的不同取舍，得到不同的预测结果，而组成一个预测区间，可供决策选用。但是未来时刻越远，模型预测的灰区间越大。可见预测精度的提高在于不断补充新的信息，使灰度逐步白化。所以，用模型预测时不是建立一个模型一直预测下去，而只用已知数列建立的模型预测的第一个预测值，补充在已知数列之后，同时为不增加数据序列的长度（维度），去掉其第一个已知数据，保持数据序列的等长度，再
25、建立模型，预测下一个值。将预测值再补充到数据序列之后，再去掉该数据序列的第一个数据。这样新陈代谢，逐个预测，依次递补，直到完成预测目的或达到一定的精度要求为止。这种预测方法称“等维灰数递补动态预测”法，即新陈代谢的灰色预测方法 。它可以达到两个目的：及时补充和利用新的信息，提高了灰色区间的白化度。当然随着递补次数的增加，灰度也在增加，信息量减少，因此也不应是无止境的。 每预测一步灰参数做一次修正，模型得到改进。这样灰参数不断修正，模型逐步改进，因而预测值都产生在动态之中。但是，需要注意的是，每次进行新陈代谢后，都需要对模型进行必要的检验，以保证建模的可行性。6.3灰色预测的检验6.3.1灰色预
26、测数列长度的选择在应用“等维灰数递补动态预测”方法，我们知道选择不同的数列长度（维度）对预测模型的影响较大。我们在中国统计年鉴2006中查询到1978年至2005年我国人口数及构成，运用matlab软件编写了“等维灰数递补动态预测”程序（程序代码见附录二）对未来45年的总人口预测。在对总人口数量的预测中，我们采用了三种数列长度。图形表示如图9，（具体数据见附录一）。图9 不同维度的灰色预测运用维度为十的“等维灰数递补动态预测”方法对未来十五年的总人口预测，得到预测结果如下表3：表3年份20062007200820092010人口数（十亿）1.31531.32311.33091.33881.34
27、67年份20112012201320142015人口数（十亿）1.35471.36271.37081.37891.387年份20162017201820192020人口数（十亿）1.39521.40351.41181.42011.4285从上图及具体的预测数据中，我们看到，在未来15年内数据预测与国家人口发展战略研究报告（见题目附录二）的吻合较好。说明等维灰色预测在本论文的人口预测中，比较适合短期预测。6.3.2 灰色预测与差分方程模型结果的比较我们用灰色预测数据与前面建立的差分方程模型得出的数据进行比较。走势图如图10。图10 灰色预测结果与差分模型结果比较从上图中，我们看到人口差分模型结论
28、在短期内是比较合理的。维度为十的灰色预测结果与差分模型结果的比较如表4（人数单位为十亿）：表4年份20062007200820092010灰色预测1.31531.32311.33091.33881.3467差分模型1.31351.32061.32741.33411.341年份20112012201320142015灰色预测1.35471.36271.37081.37891.387差分模型1.3481.35521.36241.36941.3758年份20162017201820192020灰色预测1.39521.40351.41181.42011.4285差分模型1.38151.38661.39
29、081.39411.3967七、 模型的评价与改进本论文通过对差分方程模型的三次改进，着重考虑城镇与乡村间的人口迁移及女性人口比例变化对我国人口增长的影响，并在中长期预测中参考了时间段参数设置法，为减小长期预测的误差带来了一种新的方法。所以本模型较多地从我国实际情况出发，深刻地描述了人口在新时期的发展特点，能较好地反映我国的实际情况。通过以上对模型的比较与检验可以看出，模型在中短期预测中有较高的准确率。在参数设置中，对某些参数进行了灵敏度分析，得出了对人口增长的不同影响程度。但是我们在整个建模过程中没有考虑到文化背景和人口素质的差异会对人口增长以及人口的迁移产生影响，如在乡村由于重男轻女、多子
30、多福等落后思想的存在会导致人口性别比的失调和人口增长快等。没有考虑到可能会由于欧洲等国家老龄化现象普遍，我国生产力水平提高致使劳动力富余而大量向国外输出劳动力而引起的人口增长变化等因素，所以预测结果精度还有待进一步提高。对于这些不足，我们可以采用设置参数来反映这些因素的作用来加以改进。参考文献：1 王彦 马伯强，20世纪80年代以来我国人口发展的数学模型和展望，北京大学学报(自然科学版) ,第39 卷 增刊；2930页，2003。2 宋健等，人口预测和人口控制，北京：人民出版社，1981。3 卢向虎、王永刚，中国“乡城”人口迁移规模的测算与分析( 1979- 2003)，西北人口，第1 期 N
31、o.1(107)2006：14-16页，2006。4 中华人民共和国国家统计局，中国统计年鉴2006，北京：中国统计出版社，2007。5 中国城市发展报告编辑委员会，2006中国城市发展报告,北京：年鉴社，2007。6 Fred Brauer,Carlos Castillo-Chavez,mathematical models in population biology and epidemiology，Springer，New York ：2001。7 国家统计局编辑部, 2006中国卫生统计年鉴,北京：国家统计局出版社，2007。8 国务院新闻办公室，中国儿童发展纲要（20012010年）
32、,北京:中国法制出版社，2001。9 李永胜，人口预测中的模型选择与参数认定，财经科学，总203期：68-72页，2004。10 韩中庚，数学建模方法及应用，北京：高等教育出版社，2006。11王海应 张昌维 ，利用灰色理论预测油田开发指标，新疆石油学院学报，第16卷第3期：2728页，2004。 附录附录一：灰色预测45年的数据维度15维度10维度5维度15维度10维度520061.32691.31891.315320291.57951.53261.506120071.33651.3271.323120301.59161.54271.51520081.3461.33531.330920311
33、.60381.55281.523920091.35561.34391.338820321.61621.56311.532920101.36531.35261.346720331.62881.57331.541920111.37531.36171.354720341.64141.58371.55120121.38561.37081.362720351.6541.59411.560220131.39621.381.370820361.66671.60461.569320141.40721.38911.378920371.67951.61521.578620151.41841.39791.38720
34、381.69241.62581.587920161.42981.40721.395220391.70541.63651.597220171.44121.41651.403520401.71861.64721.606620181.45241.42591.411820411.73181.65811.61620191.46331.43531.420120421.74521.6691.625520201.47361.44471.428520431.75861.67991.635120211.4851.45421.436920441.77221.6911.644720221.49661.46381.44
35、5420451.78581.70211.654420231.50841.47341.453920461.79951.71331.664120241.52011.48321.462520471.81341.72461.673920251.53191.49291.471120481.82731.73591.683720261.54381.50271.479820491.84141.74731.693620271.55561.51261.488520501.85561.75881.703520281.56751.52261.4973附录二：灰色预测的MATLAB程序 total=34743.9554
36、 21468.0445 74544.0000;load dmat M;load bmat b;z=1.8,2.1,2.5;r=0.4;s=0.47 0.45 0.45;yer=10.1,17.3,24.5;%婴儿死亡率yer=10.1,17.3,24.5;bili=0.0105198020.001443894;x1=zeros(3,45,90);women=zeros(3,45,91);men=zeros(3,45,91);%j=1;k=3;for j=1:45 q=sum(total)/100;for k=3:-1:1 if j=1 for i=1:90 if k=3 women(k,j,i
37、+1)=(M(i,k+3*(k-1)+2)/100)*total(k)*(1-(M(i,k+3*(k-1)+3)/1000)-r*q*bili(i); men(k,j,i+1)=(M(i,k+3*(k-1)/100)*total(k)*(1-(M(i,k+3*(k-1)+1)/1000)-(1-r)*q*bili(i);%q*f1*bili(1); elseif k=1women(k,j,i+1)=(M(i,k+3*(k-1)+2)/100)*total(k)*(1-(M(i,k+3*(k-1)+3)/1000)+r*(total(1)/(total(1)+total(2)*q*bili(i)
38、;%q*f1*bili(1); men(k,j,i+1)=(M(i,k+3*(k-1)/100)*total(k)*(1-(M(i,k+3*(k-1)+1)/1000)+(1-r)*(total(1)/(total(1)+total(2)*q*bili(i);%q*f1*bili(1) else women(k,j,i+1)=(M(i,k+3*(k-1)+2)/100)*total(k)*(1-(M(i,k+3*(k-1)+3)/1000)+r*(total(2)/(total(1)+total(2)*q*bili(i);%q*f1*bili(1);men(k,j,i+1)=(M(i,k+3*
39、(k-1)/100)*total(k)*(1-(M(i,k+3*(k-1)+1)/1000)+(1-r)*(total(2)/(total(1)+total(2)*q*bili(i);end x1(k,j,i)=women(k,j,i+1)+men(k,j,i+1);endelsefor i=1:90if k=3 women(k,j,i+1)=(1-M(i,3+k+3*(k-1)/1000)*women(k,j-1,i)-r*q*bili(i);men(k,j,i+1)=(1-M(i,1+k+3*(k-1)/1000)*men(k,j-1,i)-(1-r)*q*bili(i);elseif k
40、=1 women(k,j,i+1)=(1-M(i,3+k+3*(k-1)/1000)*women(k,j-1,i)+(total(1)/(total(1)+total(2)*r*q*bili(i); men(k,j,i+1)=(1-M(i,1+k+3*(k-1)/1000)*men(k,j-1,i)+(total(1)/(total(1)+total(2)*(1-r)*q*bili(i);else women(k,j,i+1)=(1-M(i,3+k+3*(k-1)/1000)*women(k,j-1,i)+(total(2)/(total(1)+total(2)*r*q*bili(i); me
41、n(k,j,i+1)=(1-M(i,1+k+3*(k-1)/1000)*men(k,j-1,i)+(total(2)/(total(1)+total(2)*(1-r)*q*bili(i);end x1(k,j,i)=women(k,j,i+1)+men(k,j,i+1);endendhi=b(:,k)/sum(b(:,k);f=0;for i=1:35f=f+(1-yer(k)/1000)*z(k)*hi(i)*women(k,j,i+15);%(women(k,j,i+15)/(men(k,j,i+15)+women(k,j,i+15)*x1(k,j,i+14);end women(k,j,1)=s(k)*f;men(k,j,1)=(1-s(k)*f;total(k)=sum(x1(k,j,:)+f;endjtotalsum(total)end 19
