1、随机性模型与模拟方法,随机变量 蒙特卡罗方法 随机数的生成 模拟,一、随机变量,何谓随机变量？随机变量是一个其值不可预测的变量。虽然一个随机变量在个别试验中其结果不确定，但在大量重复试验中其结果是具有统计规律的。正是随机变量的这种规律性使我们可以利用它来建模。例如我们可以利用下述的数据： 得出一个模型。,时间t（秒） 0 1 2 3 4 5 6 7 8 9变量X 1 0 2 2 1 2 0 1 0 2,是一个离散的随机变量并取值于 0，1和2。我们不可能给出 与 的确定的关系式，但是可以通过数 的不同值出现次数来描述这随机型 的规律列表如下： 这个表给出了随机变量 的变化规律，频率告诉某个特定
2、的事件发生的频繁程度。如果我们需要构造一个含有随机变量的模型，可以假设这个规律总是成立的，模型的假设可以基于这几个数据之上。实际操作时可以把频率分布当作概率函数来处理，但应注意概率是频率的极限值，这两者是有差异的。在处理一个简单的理论模型时，对概率函数,0 1 2频数 3 3 4频率 0.3 0.3 0.4,必须作出合适的选择。例如，假设在上述问题中的随机变量取三个值时等于可能的，这样其概率函数为这个例子说明在处理随机变量的模型时有以下两种选择： （1）使用一个理论模型。这在任何一本概率统计的书上都可以找到一些标准的理论模型如二项分布等。每一个都基于一定的假设之下成立的，所以在选用时要特别注意
3、其假设条件。 （2）使用基于实际数据的频率表，并不去套用不准理论模型。,0 1 2,使用前者的好处在于能精确地叙述变量的概率，在处理问题时可以充分发挥数理统计的作用。但这一好处把所求模式制约在了处理简单情形。随着复杂性的增加，数学就变的太难。使用后者的好处在于模型时基于观测到的数据而不是基于假设之上。增加复杂性并不成为一大障碍，但我们不再能利用数理统计而得求助于模拟以及模型的统计结果。在建立随机性模型时，首先要注意，将要处理的是离散还是连续的随机变量。 1、离散随机变量 离散随机变量的理论模型是由概率函数 来刻画的。这个式子说明随机变量 取值 时的概率。对于离散型的随机变量有下面三种重要的分布
4、,（01）分布 设随机变量 只可能取0、1两个值，它的分布规律是 则称 服从（01）分布。对于一个随机实验，如果它的样本空间只包含两个元素，即 ，我们总能在 上定义一个服从（01）分布的随机变量来描述这个随机实验的结果。例如，对新生儿的性别进行登记，检查产品的质量是否合格等都可以用（01）分布的随机变量来描述。,（2）二项分布 设实验 只有两个可能的结果，将 独立地重复地进行 次，则称这一串重复的独立实验为 重贝努利实验。它是一重和重要的数学模型，有着广泛的应用。若用 表示 重贝努利实验中事件 发生的次数， 是一个随机变量，它服从如下的二项分布 特别，当 时二项分布就是（01）分布。,（3）泊
5、松分布 设随机变量 所有可能的取值为 而取各个值的概率为 其中， 是常数，则称 服从参数为 的泊松分布。可以证明当 很小时，以 为参数的二项分布，当 时趋于以 为参数的泊松分布，其中,2、连续的随机变量,理论模型的连续型随机变量可以由概率密度函数 来描述，对所有的 存在 ，且 ，随机变量落在区间 的概率可由 来给出，在连续型随机变 量中下述两种是重要的 。,(1)均匀分布 设连续型随机变量 具有概率密度则称 在区间（a，b）上服从均匀分布。 在区间（a，b）上服从均匀分布的随机变量 ，具有下述意义的等可能性，即它落在区间（a，b）中任意等长度的子区间内的可能性是相同的，或者说它落在子区间内的概
6、率只依赖于子区间的长度而与子区间的位置无关。 （2）正态分布 设连续型随机变量 的概率密度为其中 为常数，则称 服从参数为 的,正态分布。 连续型随机变量的值如同离散的一样可以用频率表给出，但不同的是离散的随机变量每个频率对应于随机变量的一个值，而对于随机变量每一个频率对应于随机变量的一个取值范围。,二、蒙特卡罗方法,蒙特卡罗方法是计算模拟的基础，其名字来源于世界著名的赌城摩纳哥的蒙特卡罗。其思想来源于著名的蒲丰投针问题。1777年法国科学家蒲丰提出了下述著名问题：平面上画有等距离 的一些平行线，取一根长度为 的针，随机地向有平行线的平面上掷去，求针与平行线相交的概率。我们用几何概型来解决这一
7、问题。设M为针落下后的中点， 表示中点M到最近一条平行线的距离， 表示针于平行线的交角，如图2.18所示。那么基本时间区域,图2.18,它为平面上的一个矩形，其面积为 。 为使针与平行线（与 最后的一条平行线）相交，其充要条件是 的面积为 ，这样针与平行线相交的概率为 设一共投掷 次（ 是一个事先选好的相当大的自然数），观察到针和直线相交的次数为 。,从上式我们看到，当比值 不变时， 值始终不变。取 为 的近似值，我们可以算出 的近似值。可以想象当投掷次数越来越多时计算的结果就越来越准确。下表时这些实验的有关资料 (此处把 折算为1）：,由此可以看出蒙特卡罗方法的基本步骤：首先，建立一个概率模
8、型，使它的某个参数等于问题的解。然后按照假设的分布，对随机变量选出具体的值（这一过程又叫着抽样），从而构造出一个确定性的模型，计算出结果。再通过几次抽样实验的结果，的到参数的统计特性，最终算出解的近似值。 蒙特卡罗方法主要用再难以定量分析的概率模型，这种模型一般的不到解析的结果，或虽然又解析结果，但计算代价太大以至不可用。也可以用在算不出解析结果的定性模型中。 用蒙特卡罗方法解题，需要根据随机变量遵循的分布规律选出具体的至，即抽样。随机变量的抽样方法很多，不同的分布采用的方法不尽相同。在计算机上的各种分布的随机数事实上都是按照一定的确定性方法产生的伪随机数。,三、随机数的生成,我们知道对于丢硬
9、币的随机结果可以用以下的离散随机变量的改里函数来描述 如果我们需要模拟随机变量的以个值或一个集合，可以用丢硬币然后记录其其结果的方法来得到，然而这具又相当的局限性，这里我们用数学程序来产生拟随机变量。即看上去是随机出现的，但并非真正的随家便朗，它们产生于一个梯推公式。不过这些拟随机数并没有明显的规律，当给于适当的伸缩之后，它们非常接近于在 区间的均匀分布。,X 0 1P(x) 0.5 0.5,这种方法的思想是，设计一个把 和 之间的整数映射到它们自身上的函数 ，然后从 开始，依次计算 例如通过下面的公式可以产生这样的一组随机变量 给定任意一个初值，如 代入公式得 ，然后用 去除得 ；同样 代入
10、公式，可以得 ，重复这一过程可以得到我们所需要的一组随机变量。在程序设计和软件包中通常用 来表示由这样，我们用它来表示从 上的均匀分布所产生的随机变量。,我们可以从它构造出另外的随机变量。例如，可以从 给出区间 上的连续均匀分布的随机变量。如果我们要生成带参数 的指数分布，可以用 。 如果我们要生成平均值未零，标准差为 1 的正态分布，可以用下列公式 和 来给出 的两个值，令 或 可以生成 型的正态分布。,为了得到离散的随机变量，我们把 分成若干部分。例如设计一个离散的随机变量有下列的概率函数。取一个RND值：如果 ，则 ；如果 ，则 ； 如果 ，则 。对于连续的随机变量除了取生成的随机变量是
11、每类的中点外，我们可以用同样的思想进行列表分类。如,0 1 2 0.3 0.3 0.4,0-10 10-15 15-20频率 0.2 0.5 0.3,的一个 值将平移到 。一个更细致的方法是用线性插值而不是取中点，即 给出 。 从已知的 模拟一个连续随机变量的理论分布，可以用以下方法：,（1）逆累积分布函数法 如果随机变量的 是 ， 则累积分布函数是 。如果把它作为一个随机变量， 是 上的均匀分布。从 上的均匀分布取一个 值，解方程 得对应得 的值，例如，设 累积分布函数为 解 得 。这就是我们所要的由这个分布所生成的 的值,(2) 排除法 对于这种方法我们需要用两个 值来生成一个 值。设 的
12、值在区间 外为 ，而 的最大值是 。 我们可以通过如下的步骤生成 的值。 从 上的均匀分布生成 和 ； 用 计算 ； 计算 ； 用 算出 ； 如果 ，则接受 ，否则排除 回到 。 对于上面的例子，我们取,四、模拟,模拟是现象的模型所产生的再现。所谓数学模拟就是用模型使现象再现。因此，表示现象的部分或总体的基本方程和表示自然规律的数学模型全是数学模拟。然而，狭义地讲主要指的是数字模拟。它是将复杂现象作出可以用数字计算机表达的数学模型，从数值上进行各种实验。各种方法随着计算机的进步已广泛地应用起来。因此我们所说的模拟主要是指数学模拟。,例 2.18 一列火车大约在下午1点离开 站其规律如下； 火车
13、从 到 途中所需要的平均时间为 分，由 分钟的标准差。如果你要赶的是这趟火车的下一站 ,而你到达 的站的时间分布为 问你能赶上这列火车的概率是多少？,离站时间 13.00 13.05 13.10 概率 0.7 0.2 0.1,时间 13.28 13.30 13.32 13.34概率 0.3 0.4 0.2 0.1,为了回答这个问题，我们需要一些随机数。这里我们将采用上面给出的那些随机数，即 等。而我们所要模拟的是 火车离站的时间 ； 火车途中的时间 ； 你到达车站 的时间 。 这样你赶上火车的条件是 。为模拟这个问题只需要生成 ， 和 的值，然后检验这条件。但如何得到 的值是不明显的，因并不知
14、道这个分布。这样，假设一个模型，取平均值为30，标准差为2的正态分布，由所给的条件知 ， 为离散的，而 为连续的随机变量。,以分为时间单位，从 的下午以点起算，构造的模型如下 其中 。 计算结果为 ， 和 ，这样 。在这种场合你比火车提前到达4分钟。但需要指出，这并不是说我们已经回答了这个问题，要回答这个问题我们要作多次这样的模拟，记下这些结果，算出能赶上火车的频率。通过足够多次的模拟之后我们就可以看出能赶上火车的概率。,一般用在模拟建模时，一次模拟的成功并不能说明什么问题，更不能说我们的主要工作已经完成。你必须多次的进行模拟，然后分析其结果。分析的种类要看模型的对象，而这在模拟的一开始就应该
15、清楚的。在实验的模拟模型的对象是在变化的，但常常包括一下几种： 对系统的长期性态作出统计； 比较系统的可选择对象的安排； 研究参数变化的影响； 研究模型假设的影响； 找出系统最优方案； 上面的例子是相当平凡的，根本不能作为用模拟解决问题的例子。下面我们仅举两个简单的例子以理顺模拟模型的思路。,例2.19 某个理发店中有两名理发员 和 ，顾客随机地来理发，据统计 的顾客仅需要剪发，化时 分钟；而有 的顾客即需要剪又需要又需要吹风，许花时 分。 对任意一个模拟，首先要作的是 找出能完全描述任意时刻的系统的状态变量的集合。 给出能从时刻 的状态变量算出时刻 的新的状态变量的程序。 这个例子中有三个状
16、态变量：在等待的顾客的人输（离散的非负整数）；是否 正在工作（是或否）；是否 正在工作（是或否）。,一次模拟式由始于 ，结束于 的状态变量的值的一系列演算组成的。一个事件是时间中的一点，在这个时刻一个随机变量改变了它的值。在这个例子中的事件有: 一个顾客到达； 开始服务； 结束服务； 开始服务； 结束服务。 一个元素是一个离散，或者是系统的长期部分，或者是进入和离开，这里的元素是顾客和两个理发员。 对研究一个模拟模型来说，有两种程序类型： （1）时间切片 考察状态变量和在时间切片中（通常是等时间的切片）元素的位置。在每一个时间切片中状态变量可变可不变。,（2）事件序列 考察在每一事件的系统，并
17、不考虑时间之间的时间。 这两种途径我们有时称为“时间传动”和“事件传动”模型。一般我们用时间传动模型于连续的的确定型系统，事件传动模型于离散的概率模型，但这不时绝对。在这个例子中我们将用“事件传动”。 对于时间切片模型，我们必须决定时间切片的大小，为简单计我们将取 分钟。问题的描述并不包括任何有关顾客到达率的信息。假设在任何一分钟顾客到达的概率是 。实际上有两种不同类型的顾客，取决于是否要吹风。我们通过取服务时间的平均值，即 分，构造一个粗糙的模型。,为了描述一个顾客是否到来这个随机变量，我们用一个硬币将作为一个随机数的生成器，用 表示反面， 表示正面。设扔出的序列是 。用 表示一个顾客到达，
18、且取初始状态为顾客，运行前 分钟，就有下表的结果：,时间（分） 到达？ A在工作 B在工作 排队 0 否 否 否 0 1 是 是 否 0 2 否 是 否 0 3 是 是 是 0 4 是 是 是 0 5 否 是 是 0,到 这里，人们将要为我们希望知道什么。通常我们感兴趣的是平均队伍的长度，最长的队伍，顾客等待的平均时间以及两个理发员的忙 闲 程度等，注意到这里有两种不同的平均，即一个是关于时间，而另一个是关于顾客的平均，为回答上述为她我们设 是任意时刻的排队的顾客数。顾客和时间的关系通常可由图 给出。,6 是 是 是 0 是 是 是 0 是 是 是 0 否 是 是 010 否 是 是 0,它是
19、一个右连续的阶梯函数是合理的，这是由于只有新顾客到来或有顾客完成服务后离去，函数值才发生变化， 关于时间的平均是 ，其中 图下额面积，设 表示一个时间区间，在其上 保持常数（这里 本身是变量）。当我们进行模拟时我们累积其和 。用 记在进行模拟期间到达的顾客数。这样我们所要的两个平均分别为,队长平均 等待时间的平均 下面是用来说明累积排队时间 的记录（注意这里仅给出 变化的时间）：,时间 Q 0 0 0 0 011.584 1 0 0 0 12.935 0 1.351 1.351 1.35117.290 1 4.355 0 1.35117.935 0 0.645 0.645 1.99618.67
20、6 1 0.741 0 1.99623.156 0 4.480 4.480 6.47625.217 1 2.061 0 6.47625.327 2 0.110 0.110 6.58625.935 1 0.608 1.216 7.80227.431 2 1.406 1.406 9.208,这里所执行的总时间 。在这期间 。队伍的最长长度 。累计排队时间 队伍的平均长度 。平均等待时间 。 在我们结束模拟时还有两个顾客，一个是排队的，而另一个是新来的。 服务的总时间是 分。因此 忙碌的概率是 。 服务的总时间是 分或忙碌时间为 。,评注1 模拟一个系统的目的不是为了模仿一个现实系统，而是通过解决问
21、题达到优化系统的目的。例如在这个例子中可以分析诸如增雇一个理发员或改变服务时间等对系统的影响。 评注2 在我们的模型中，为使问题简单我们已经作了一些假设： （1）假设了在任何一分钟有一个顾客到达的概率是 。 （2）默认在同一分钟内的顾客数 。 （3）如果两个理发师均空闲，顾客可以任意选。 （4）排队的原则是安先后的秩序。如果有预约可以先服务。,（5）我们的模型中允许一下情形出现，一个顾客的来到，发现有很多人在等就走啦。也可能是一个顾客在等了一段时间之后等不及了就离开了。意味这允许其中的一个理发员有短暂的休息。 例2.20 倒媒台的操作方案 某煤矿公司有一个大型煤台，用于向运媒列车装煤。该倒煤台
22、的容量是 列标准列车。装满一个空的倒煤台需要一个小组 个小时的时间，费用是 。为提高装煤速度可以以 的代价动用第二个小组。铁道部门每天向这个倒煤台发三列空的标准车。这些列车可在上午 点倒下午 点之间的任何时刻到达。给一列标准车装满煤需要 小时，向倒煤台装煤和从倒煤台向列车装煤不能同时进行。如果列车到达后因等待装煤二停滞，铁道部门将征收每车 的滞费。,此外，每星期四上午 点到下午 点之间还有一列大容量列车到达，其容量为标准的列车的 倍，滞期费为 。请问 （1）安标准规则操作可使装煤费最低？费用使多少？ （2）如果标准列车能在指定的时间到达，什么样的调度安排最经济 ？ 这个问题中列车的到达时间使随
23、机因素，适合于建立概率模型，同计算机模拟加以解决。 首先，模型中需要考虑的费用由两部分组成。一部分使装煤小组向倒煤台装煤的费用，记为 ，另一部分使列车等待装煤的滞期费 。因每天要装的煤数量使固定的， 的大小只受是否使用大二小组影响。,通过使用第二小组，有可能减少 。模型的主要任务是将总费用 降到最低。故 是模型的目标函数。 其次，由于理论上的困难，很难得到最优方案。考虑到这是一个每天重复发生的为她，重要的的是提供一组简单明确的规划，使煤矿公式可以根据规则方便地获得接近最优的解。因此，我们将在方案的优化程度和简明性之间做一个折中。 设： 为装满列车 所需的煤量； 为倒煤台中剩下煤量； 表示当前时
24、间。其中 和 均以 小时向列车装的煤量为单位。,根据题意写出下面一些应该遵循的规则： 有列车等待时，用两个小组装煤节省的滞期费大于增加的装煤费用，此时应使用第二个小组。 当同时有两列或三列标准列车等待装煤时。应将已装煤量最多的车排在前面先装，已装煤量最少的车排到最后面。可以证明，这样安排滞期费最少。 当同时有大容量车 和标准车 等待时，先装 后装 的滞期费 先装 时的滞期费为 当 时，先装 ，否则先装 。,设当前待装的车为 ，则用两个小组装倒煤台直到 或 为止，然后装列车。 周四时，装标准车和大容量列车共需要 小时。即便是倒煤台在周四上午 点以前就已提前装满，当天用两个小组装倒煤台仍需 小时，
25、合计 小时故最快也要到周五早上 点才能完成周四的任务，且此时倒煤台为空。为保证周五正常工作，应马上开始装倒煤台。由以上分析知，周时间最紧张，就始终用两个小组。 非周四，在此刻 无列车等待，设已知下一列 车到达时间为 。若 ，则时间充 足，可以用一个小组装倒煤台至满或下一列车来。否则用两个小组。,非周四，不知道列车的到达时间。设在时刻 倒煤台中尚有煤量 ，没有列车等待，当天还有 列标准车未到达。假设列车到达时间服从独立的均匀分布，则存在 ，当 时用一个组装煤即可，否则要用两个组。 的选择应满足最小的原则，因其解析解难以求出，故采用计算机模拟的方法。首先任意取一个 值 ，注意到 ，在上述约束条件下
26、以一定步长 取 各种组合，分别用计算机模拟求出平均费用，找出使平均费用最少的一组 ， 和 值， 作为在该组合给定 下的函数值。选取一系列不同的 的值重复以上过程，就可以得到函数 在各点上的值。,在以上规则指导下，我们用时间切片法进行模拟，流程图如图 所示。模拟的结果如下：年滞期费用 元，年度总费用 元。,开始,模拟时钟 0,初始化系统状态和时间队列,找出最近的下次事件。推出模拟时钟到该事件的时刻,计算新的系统状态,产生未来事件（可能没有，也可能有一个或多个）并加入事件列队,结束条件满足吗？,计算并输出统计结果,结束,当标准列车的到达时间可以确定时，分三种情况考虑（记 ， 和 为三列车的到达时间
27、，且 ） 非周四，周五 可以推出，滞期费用为 ，且不使用第二小组，当且仅当 成立。不等式组 的解不唯一，任何一小组即可，如取 。 周四 标准列车的到达时间应尽量和大容量列车错开，故取 。再此前提下用模拟的方法确定 ，得 时费用最少。,周五 因周四工作量大，将积压到周五（周四的最后一列车最快能再周五早上四点装完，最慢要拖到六点）。 为减少等待， 发 车 的 时 间 尽 量 靠 后，故取 。 在计算机模拟中事件序列比较常用，我们不在举例说明，这个方法的流程图可以用图 来说明。,开始,，初始化统计数据,用独立的均匀分布随机数 产生三列标准车的到达时刻周四时产生大容量车的到达时刻,在时刻 倒煤台处有车等待 吗?,当天列车到齐了吗？,装满倒煤台，时间 推进到第二天,模拟天数足够了吗？,输出模拟结果,结束,按规则 选出待装列车,按规则 装车 或装倒煤台,倒煤台是满的吗?,按规则 装倒煤台,
