项目
03/14/2016

2016 年 3 月

第 31 卷，第 3 期

C# - 离散事件模拟：人口增长示例

从古至今，多个科学的发展得益于模拟能力。模拟人体的医学模型加强人体解剖学研究。计算机模拟游戏（如“魔兽世界”）重新创建整个幻想世界，“飞行模拟器”帮助在地面上培训飞行员。不同的模拟程序探索恐怖袭击的反应、流行性疾病和其他可能的危机。甚至电影“侏罗纪公园”中的模拟恐龙都暗示了模拟及其潜在性的广泛应用。

模拟是一项技术，可在其中由设计的模型模仿现实生活系统或过程。模型封装了该系统的所有功能和行为；模拟是随着时间的推移执行本系统。设计一个模拟有以下几个阶段：

定义要建模的系统，其中包括研究手头的问题、识别环境的属性和指定要达到的目标。
制定模型，其中包括定义其所有变量及其逻辑关系以及创建必要的流程图。
定义模型将需要的数据，以产生预期结果。
生产计算机化的模型实施。
验证实施的模型是否满足设计。
通过比较验证模拟器确实呈现了模拟的真实系统。
尝试使用此模拟器生成所需的数据。
分析和解释模拟器的结果，并基于这些结果作出决定。
将创建的模型和模拟器记录为工具。

模拟一般包括连续处理或离散事件。要模拟天气系统，例如，随着所有元素不断变化跟踪持续发生。因此，根据时间变量放置的温度变量会以连续的曲线来表示。相反，飞机起飞或降落随着时间点发生，因此，模拟可以只考虑那些精确的时刻或事件并抛弃其余一切。这种类型的模拟称为离散事件模拟 (DES)，这正是我将在本文中讨论的内容。

离散事件模拟

DES 将系统或进程建模成一段时间（即，从一事件时间到另一事件时间）内有序排列的各个事件。因此，DES 模拟中的时间通常比实际时间要短得多。

开发 DES 时，需要考虑以下六个主要元素：

对象代表实际系统的元素。它们具有属性，它们与事件相关，它们消耗资源，它们随着时间的推移进入和离开队列。在前面提到过的飞机起飞和降落方案中，对象是飞机。在医疗保健系统中，对象可能是病人或器官。在仓库系统中，对象是库存产品。对象应该是相互或与系统交互，而且在模拟期间可以随时创建对象。

属性是特定于每个对象存储的特征（大小、登录时间、性别、价格等），以便确定对在模拟中可能发生的各种方案的反应；此类值可以修改。

事件是可以在系统中发生的事情，尤其是对象，如飞机的着陆、产品到达仓库、特定疾病的出现等等。事件可以按任何顺序发生或再发生。

资源是为对象提供服务的元素（例如，机场中飞机的起落跑道、仓库中的存储单元和诊所中的医生）。当某个对象需要一个已被占用的资源时，该对象必须排队等待直到该资源可用。

队列是在其中组织对象以等待释放当前占用的资源的管道。队列可以有最大的容量，而且可以有不同的调用方法：先进先出 (FIFO)、后进先出 (LIFO) 或基于某些标准或优先级（疾病级数、燃料消耗等等）。

时间（因为它发生在现实生活中）在模拟中是必不可少的。要测量时间，时钟在模拟开始时就要启动，然后可以用来跟踪特定的时间段（出发或到达时间、运输时间、某些症状花费的时间等等）。这种跟踪十分重要，因为通过它您可以知道下一个事件应何时发生。

因为模拟编程很复杂，为了便于开发，已经有许多创造包含模拟范例所有要求的语言方面的尝试。SIMULA 就是一种这样的语言，在 20 世纪 60 年代由 Ole Johan 和 Kristen Nygaard 发明，首次引入面向对象编程 (OOP) 的概念，是当今领先的编程范例。如今，更注重于创建可将创建模拟时程序员需要的内容包含在内的程序包、框架或库。这些库应从 C#、C++、Java 或 Python 等普通语言调用。

在“离散事件模拟编程语言的历史记录”中，Nance 提出了任何 DES 编程语言应满足的六个最低要求 (bit.ly/1ZnvkVn)：

随机数字生成。
进程转换器，允许统一的随机变量以外的变量。
列表处理，方便对象的创建、操作和删除。
统计分析，提供模型行为的描述性摘要。
报告生成，协助大型数据集的演示和促进决策制定。
时间流机制。

这些要求都可以在 C# 中得以满足，因此在本文中我将演示以 C# 开发的人口增长的离散事件模拟。

人口增长的 DES

在种群（动物和植物）如何随着时间和空间变化以及与其所在环境交互的研究中，人口增长是要考虑的许多方面之一。人口是生活在同一时期的相同物种的生物群体，有时位于同一空间或依据具体特征进一步划分。

为什么研究人口增长很重要？更好的了解人口如何增长或缩减使科学家能够更好地预测生物多样性保护、资源利用、气候变化、污染、医疗、运输需求等诸多方面中的的变化。还可以深入了解生物之间以及与环境的交互，这是在考虑人口是繁荣还是衰退时需要考虑的重要方面。

在本文中，我将演示人口增长的 DES。目标是观察人口随着时间推移如何演变，以及通过分析最终结果（人口规模、老年人、年轻人等等）得到一些统计信息。人口将以 m 位男性和 n 位女性个体开始，每位都有相应的年龄。显然，m 和 n 必须大于零，否则模拟毫无意义。此模型中的对象是个体。

达到参数 λ = 18 岁的泊松函数分布的年龄后，一个个体能与另一个个体开始交往。（您不需要马上完全理解泊松分布、正态分布或指数分布；下一节将会进行解释。）

单身且有行为能力的异性个体彼此接触的概率小于 50%，甚至仅在年龄差不超过 5 年时发生。图 1 显示了两个个体之间的关系结束概率的一个示例。

图 1 关系结束的概率

平均年龄	可能性
14-20	0.7
21-28	0.5
29 以上	0.2

在参数 λ= 8 年的指数函数分布的某个时间后，关系中涉及的个体可以有孩子。

一个女人如果其年龄符合参数 µ = 28 和 σ²= 8 年的正态（钟形）分布函数便可以怀孕。每个女人都有许多孩子，按参数 µ = 2 和 σ²= 6 年的正态函数分布。（参数 µ 代表平均年龄，而 σ²是年龄差异的度量值。）

每个个体都有一个预期寿命，按参数 λ = 70 岁的泊松函数分布，其中 λ 代表平均预期寿命。

在前面的描述中可以确定以下几个事件：

开始一段关系。
结束一段关系。
怀孕。
生孩子。
死亡。

每个事件都伴随着一个离散随机变量，确定事件将在哪个时刻发生。

概率分布

离散随机变量的值的集合是有限的或可数无穷的。即，值可列为有限或无限序列的值 1, 2, 3 …。对于一个离散随机变量，其概率分布是可将概率分配给每个可能值的任意的图、表或公式。所有概率之和必须为 1，且每个单独的概率必须介于 0 和 1 之间。例如，当投掷公平合理的骰子（各面几率均等）时，代表可能结果的离散随机变量 X 将获得以下概率分布：X(1) = 1/6, X(2) = 1/6, …, X(6) = 1/6。各面几率均等，因此分配给各面的概率均为 1/6。

参数 l 和 µ 表示其相应分布中的平均值（期望值）。平均值表示随机变量平均采用的值。换句话说，它是总和 E=[(每个可能的结果) × (对应结果的概率)]，其中 E 表示平均值。如果是骰子，平均值将为 E = 1/6 + 2/6 + 3/6 + 4/6 + 5/6 + 6/6 = 3.5。请注意，结果 3.5 实际上是骰子可以采用的所有可能值之间的中间位置；它是多次投掷骰子时的期望值。

参数 σ2 表示分布的方差。方差表示随机变量的可能值的分散性，并且它总是非负的。小方差（接近 0）表明值相互之间很接近而且接近平均值；高方差（接近 1）表明值之间的最大距离以及与平均值之间的距离。

泊松是一个离散分布，表示每时间单位事件数量的概率（参见图 2）。当事件的概率小且发生几率大时，通常应用此分布。一本书中印刷错误的数量、到达某个商业中心的客户数量、到达交通灯的汽车数量和给定年龄组中每年死亡人数都是泊松分布的应用示例。

图 2 泊松分布，参数 λ = 18

指数分布表示泊松过程中事件之间的时间（参见图 3）。例如，如果您正在处理一个描述在某段时间内到达某个商业中心的客户数量的泊松过程，您可能会对指明在第一个客户到达前已过去了多少时间的随机变量感兴趣。指数分布可以实现这个目的。它还可以应用于物理过程，例如用于表示粒子的生命周期，其中 λ 表示粒子寿命增长的速率。

图 3 指数分布，参数 λ = 8

正态分布描述倾向于围绕一个中心值的概率，没有左偏置或右偏置，如图 4 所示。正态分布是对称的，而且具有钟形密度曲线，在平均值处有一个峰值。平均值左右两侧各有 50% 的分布。标准偏差表示钟形曲线的传播或周长；标准偏差越小数据越集中。平均值和标准偏差必须都表示为正态分布的参数。许多自然现象都密切遵循正态分布：血压、人的身高、测量误差等等。

图 4 正态分布，参数 µ = 28 和 σ2 = 8 年

现在，我将向您展示如何以一种流行、成熟且简洁的语言（如 C#）实施拟议的 DES。

实施

为了开发此模拟，我将利用 OOP 范例的所有优点。我的理念是尽可能获得可读性最高的代码。科学应用程序往往是复杂且难以理解的，有必要使它们尽可能清晰，以便其他人可以理解您的代码。我将利用图 5 中所示的结构在 C# 中开发控制台应用程序。

图 5 模拟应用程序的结构

逻辑路径是必不可少的；注意命名空间结构如何维护其自己的常识： Simulation.DES.PopulationGrowth。

Events 文件夹包含一个 Event.cs 文件，用于定义一个表示模拟中的所有可能事件的枚举类型：

namespace DES.PopulationGrowth.Events
{
public enum Event
  {
Capable Engaging,
Birth EngageDisengage,
GetPregnant,
ChildrenCount,
TimeChildren,
    Die
  }
}

Objects 文件夹包含与个体相关的每个类；这些是将从 OOP 受益最多的代码段。有三个与个体相关的类： Individual、Male 和 Female。第一个是抽象类且其他类都从它继承；即，男性和女性都是个体。大多数编码发生在 Individual 类中。图 6 显示了它的属性和方法。

图 6 Individual 类的属性和方法

以下是每个属性的描述：

Age：个体的年龄。
Couple：如果个体是单身，值为 null；否则，值为他/她的配偶（另一个体）。
Engaged：如果个体涉及到一段关系，值为 ture；否则，值为 false。
LifeTime：个体死亡的年龄。
RelationAge：个体可以开始一段关系的年龄。
TimeChildren：模拟中个体可以有孩子的时间（不是年龄）。

以下是其对应方法的描述：

Disengage：结束两个个体之间的关系。
EndRelation：根据图 1 中的概率表确定关系是否应结束。
FindPartner：找到某个个体的适宜合作伙伴。
SuitablePartner：确定一个个体是否适合开始一段关系（年龄差异和异性）。
SuitableRelation：确定某人是否可以开始一段关系；即，是单身且年龄适合开始交往。
ToString：重写以将个体信息表示为字符串。

类首先定义所有属性，之后是只设置年龄的构造函数：

public abstract class Individual
  {
public int Age { get; set; }
public int RelationAge{ get; set; }
public int LifeTime{ get; set; }
public double TimeChildren{ get; set; }
public Individual Couple { get; set; }
protected Individual(int age)
{
  Age = age;
}

SuitableRelation 和 SuitablePartner 方法及 Engaged 属性只是逻辑测试，相当简单：

public bool SuitableRelation()
{
return Age >= RelationAge&& Couple == null;
}
public bool SuitablePartner(Individual individual)
{
return ((individual is Male && this is Female) ||
(individual is Female && this is Male)) &&Math.Abs(individual.Age - Age) <= 5;
}  
public bool Engaged
{
get { return Couple != null; }
}

Disengage 方法通过先后将配偶和个体上的 Couple 字段设置为 null 结束关系。它也将他们生孩子的时间设置为 0，因为他们不再有关系。

public void Disengage()
{
Couple.Couple = null;
Couple = null;
TimeChildren = 0;
}

EndRelation 方法基本上是概率表的翻译，用于确定一段关系结束的可能性。它使用一个统一的随机变量，该变量产生一个范围 [0, 1] 内的随机值，相当于产生一个接受率。分布字典在模拟类（稍后描述）中创建，并保留项对（事件、概率分布），从而将每个事件与其分布相关联：

public bool EndRelation(Dictionary<Event, IDistribution> distributions)
{
var sample =
  ((ContinuousUniform) distributions[Event.BirthEngageDisengage]).Sample();
if (Age >= 14 && Age <= 20 && sample <= 0.7)
return true;
if (Age >= 21 && Age <= 28 && sample <= 0.5)
return true;
if (Age >= 29 && sample <= 0.2)
return true;
return false;
}

图 7 中显示的 FindPartner 方法为实例个体查找适宜的合作伙伴。它接收人口列表、模拟的当前时间和分布字典作为输入。

图 7 FindPartner 方法

public void FindPartner(IEnumerable<Individual> population, int currentTime,
  Dictionary<Event, IDistribution> distributions)
{
foreach (var candidate in population)
if (SuitablePartner(candidate) &&
candidate.SuitableRelation() &&
((ContinuousUniform) distributions[Event.BirthEngageDisengage]).Sample() <= 0.5)
{
// Relate them
candidate.Couple = this;
Couple = candidate;
// Set time for having child
var childTime = ((Exponential)  distributions[Event.TimeChildren]).Sample()*100;
// They can have children on the simulated year: 'currentTime + childTime'.
candidate.TimeChildren = currentTime + childTime;
TimeChildren = currentTime + childTime;
break;
}
}

最后，ToString 方法定义某个体的字符串表示形式：

public override string ToString()
{
Return string.Format("Age: {0} Lifetime {1}", Age, LifeTime);
}

Male 类非常简单：

public class Male: Individual
{
public Male(int age) : base(age)
{
}
public override string ToString()
{
return base.ToString() + " Male";
}
}

Female 类更为复杂，因为它包括处理怀孕、生育等。类定义首先声明所有属性和构造函数：

public class Female :Individual
{
public bool IsPregnant{ get; set; }
public double PregnancyAge{ get; set; }
public double ChildrenCount{ get; set; }
public Female(int age) : base(age)
{
}
}

以下是 Female 类的属性：

IsPregnant：确定这个女人是否已怀孕。
PregnancyAge：确定一个女人可以怀孕的年龄。
ChildrenCount：指示要出生的孩子的数量。

以下是它包含的方法：

SuitablePregnancy：确定这个女人是否可以怀孕。
GiveBirth：表示一个女人生育，为人口增加新个体。
ToString:override：用于将女性信息表示为字符串。

SuitablePregnancy 是一个测试方法，当实例女人满足怀孕的每个条件时输出 true：

public bool SuitablePregnancy(intcurrentTime)
  {
return Age >= PregnancyAge && currentTime <= TimeChildren && ChildrenCount > 0;
}

图 8 中显示的 GiveBirth 方法是代码，可将新个体添加到人口并进行初始化。

图 8 GiveBirth 方法

public Individual GiveBirth(Dictionary<Event, IDistribution> distributions, int currentTime)
  {
var sample =
  ((ContinuousUniform) distributions[Event.BirthEngageDisengage]).Sample();
var child = sample > 0.5 ? (Individual) newMale(0) : newFemale(0);
// One less child to give birth to
ChildrenCount--;
child.LifeTime = ((Poisson)distributions[Event.Die]).Sample();
child.RelationAge = ((Poisson)distributions[Event.CapableEngaging]).Sample();
if (child is Female)
{
(child as Female).PregnancyAge =
  ((Normal)distributions[Event.GetPregnant]).Sample();
(child as Female).ChildrenCount =
  ((Normal)distributions[Event.ChildrenCount]).Sample();
}
if (Engaged && ChildrenCount> 0)
{
TimeChildren =
  currentTime + ((Exponential)distributions[Event.TimeChildren]).Sample();
Couple.TimeChildren = TimeChildren;
}
else
TimeChildren = 0;
IsPregnant = false;
return child;
  }

首先生成一个统一示例以确定新个体将为女性还是男性，概率各为 50% (0.5)。ChildrenCount 值递减 1，表明这个女人少生一个孩子。代码的其余部分涉及个体初始化和重置一些变量。

ToString 方法更改 Females 表示为字符串的方式：

public override string ToString()
  {
return base.ToString() + " Female";
}

因为与个体相关的所有函数都被放置在单独的类中，所以 Simulation 类现在更加简单。属性和构造函数位于代码的开头：

public class Simulation
  {
public List<Individual> Population { get; set; }
public int Time { get; set; }
private int _currentTime;
private readonly Dictionary<Event, IDistribution> _distributions ;

属性和变量是自我描述的，而且之前描述了一些。Time 表示模拟将持续的时间（以年为单位），而 _currentTime 表示模拟的当前年份。本示例中的构造函数（如图 9 所示）比较复杂，因为它包括每个个体的随机变量的初始化。

图 9 Simulation 类构造函数

public Simulation(IEnumerable<Individual> population, int time)
{
Population = new List<Individual>(population);
Time = time;
_distributions = new Dictionary<Event, IDistribution>
{
{ Event.CapableEngaging, new Poisson(18) },
{ Event.BirthEngageDisengage, new ContinuousUniform() },
{ Event.GetPregnant, new Normal(28, 8) },
{ Event.ChildrenCount, new Normal(2, 6) },
{ Event.TimeChildren, new Exponential(8) },
{ Event.Die, new Poisson(70) },
                                 };
foreach (var individual in Population)
{
// LifeTime
individual.LifeTime = ((Poisson) _distributions[Event.Die]).Sample();
// Ready to start having relations
individual.RelationAge = ((Poisson)_distributions[Event.CapableEngaging]).Sample();
// Pregnancy Age (only women)
if (individual is Female)
                {
(individual as Female).PregnancyAge = ((Normal) _distributions[Event.GetPregnant]).Sample();
(individual as Female).ChildrenCount = ((Normal)_distributions[Event.ChildrenCount]).Sample();
}
}
}

最后，Execute 方法（如图 10 所示）是所有模拟逻辑发生的地方。

图 10 Execute 方法

public void Execute()
{
while (_currentTime< Time)
{
// Check what happens to every individual this year
for (vari = 0; i<Population.Count; i++)
{
var individual = Population[i];
// Event -> Birth
if (individual is Female&& (individual as Female).IsPregnant)
Population.Add((individual as Female).GiveBirth(_distributions,
  _currentTime));
// Event -> Check whether someone starts a relationship this year
if (individual.SuitableRelation())
individual.FindPartner(Population, _currentTime, _distributions);
// Events where having an engaged individual represents a prerequisite
if (individual.Engaged)
{
// Event -> Check whether arelationship ends this year
if (individual.EndRelation(_distributions))
individual.Disengage();
// Event -> Check whether a couple can have a child now
if (individual is Female &&
(individual as Female).SuitablePregnancy(_currentTime))
(individual as Female).IsPregnant = true;
}
// Event -> Check whether someone dies this year
if (individual.Age.Equals(individual.LifeTime))
{
// Case: Individual in relationship (break relation)
if (individual.Engaged)
individual.Disengage();
Population.RemoveAt(i);
}
individual.Age++;
_currentTime++;
}
}

外部循环中的迭代表示模拟中经过的年数。内部循环遍历在特定年份内可能发生在那些个体上的事件。

为了了解人口随着时间推移如何发展，我建立了一个新的控制台应用程序，如图 11 所示。

图 11 查看随着时间推移的人口

static void main()
{
var population = new List<Individual>
{
new Male(2),
new Female(2),
new Male(3),
new Female(4),
new Male(5),
new Female(3)
};
var sim = new Simulation(population, 1000);
sim.Execute();
// Print population after simulation
foreach (var individual in sim.Population)
Console.WriteLine("Individual {0}", individual);
Console.ReadLine();
}

图 12 显示了 1,000 年后的人口。

图 12 1,000 年后的人口

在本文中，我开发了一个离散事件模拟，以了解随着时间推移人口如何发展。面向对象的方法被证明是非常有用的，可用于获得可读的、简明代码，如有必要，读者可以尝试和改进这些代码。

Arnaldo Pérez Castaño* 是古巴的计算机科学家。他是一系列编程书籍的作者—“JavaScript Fácil”、“HTML y CSS Fácil”和“Python Fácil”(Marcombo S.A)。他的专长包括 Visual Basic、C#、.NET Framework 和人工智能，并作为一个自由职业者通过 nubelo.com 提供服务。电影和音乐是他的一些爱好。通过 arnaldo.skywalker@gmail.com* 与他联系。

衷心感谢以下 Microsoft 技术专家对本文的审阅： James McCaffrey

C# - 离散事件模拟： 人口增长示例

离散事件模拟

人口增长的 DES

概率分布

实施

其他资源

C# - 离散事件模拟：人口增长示例