Jiangtang's profile技止于此BlogListsNetwork Tools Help

Blog


    1/12/2008

    SAS和蒙特卡罗模拟(2):随机数基础

    **************************************************************************************************************************

    SAS和蒙特卡罗模拟(1):开篇》——简介,通过例子建立起蒙卡的直观概念;参考软件包及书目

    **************************************************************************************************************************************************

    SAS for Monte Carlo Simulations (2): Random Numbers

    最简单、最基本、最重要的随机变量是在[0,1]上均匀分布的随机变量。一般地,我们把[0,1]上均匀分布随机变量的抽样值称为随机数,其他分布随机变量的抽样都是借助于随机数来实现的。以下谈的都是所谓“伪随机数”(Pseudo Random Numbers)。产生随机数,可以通过物理方法取得(很久很久以前,兰德公司就曾以随机脉冲源做信息源,利用电子旋转轮来产生随机数表),但当今最为普遍的乃是在计算机上利用数学方法产生随机数。这种随机数根据特定的迭代公式计算出来,初值确定后,序列就可以预测出来,所以不能算是真正的随机数(就成为“伪随机数”)。不过,在应用中,只要产生的伪随机数列能通过一系列统计检验,就可以把它们当成“真”随机数来用。

    现在大多软件包内置的随机数产生程序,都是使用同余法(Congruential Random Numbers Generators)。“同余”是数论中的概念。

    0.预备知识:同余

    捡回小学一年级的东西:"4/2=2"读作“4除以2等于2”,或者,“24等于2”。还有求模的符号mod(number,divisor),其中,number是被除数(在上式中,为4),divisor是除数(上式中的2)。这样的约定对SAS和Excel都通用,如mod(4,13)=4,mod(13,4)=1。

    现在我们可以开讲“同余”了。设m是正整数,用m去整数a、b,得到的余数相同,则称“a与b关于模m同余”。上面的定义可以读写成,对整数a、b和正整数m,若mod(a,m)=mod(b,m),则称“a与b关于模m有相同的余数(同余)”,记做a≡b(mod m)(这就是同余式)。举个例子,mod(13,4)=1,mod(1,4)=1,则读成13和1关于模4同余,记做13≡1(mod 4)。当然,同余具有对称性,上式还可以写成1=13(mod 4)。

    a≡b(mod m)的一个充要条件是a=b+mt,t是整数,比如13=1+3*4。a=b+mt可以写成(a-b)/m=t,即m能整除(a-b)。

    这些就够了。更多基础性的介绍,可以参考《同余(数据基础)》。

    1.乘同余法

    同余法是一大类方法的统称,包括加同余法、乘同余法等。因为这些方法中的迭代公式都可以写成上面我们见过的同余式形式,故统称同余法。常用的就是下面的乘同余法(Multiplicative Congruential Generator.)。符号不好敲,做些约定,如R(i)就是R加一个下标i。

    乘同余法随机数生成器的同余形式如下:R(i+1)=a*R(i) (mod m)。这个迭代式可以写成更直观的形式,R(i+1)=mod[a*R(i),m],其中初值R(0)称为随机数种子。因为mod(x,m)总是等于0到m-1的一个整数,所以最后把R(i+1)这个随机序列都除以m,就可以得到在[0,1]上均匀分布的随机数。下面用电子表格演示一遍,假设随机数种子R(0)=1,a=4,m=13:

      A B C D E
    1 i R(i) a*R(i) mod[a*R(i),m] mod[a*R(i),m]/m
    2 0 1 =B2*4 =mod(C2,13) =D2/13
    3 1 =D2 =B3*4 =mod(C3,13) =D3/13
    4 2 =D3 =B4*4 =mod(C4,13) =D4/13

    上面显示的是公式(Excel中,公式与计算结果的转换,用快捷键ctrl+~实现),结果看起来是这样的,E列就是我们想要的在[0,1]上均匀分布的随机数:

      A B C D E
    1 i R(i) a*R(i) mod[a*R(i),m] mod[a*R(i),m]/m
    2 0 1 4 4

    0.307692308

    3 1 4 16 3

    0.230769231

    4 2 3 12 12

    0.923076923

    上表中,a*R(1)=16,R(2)=3,m=13,一个同余式就是R(2)=a*R(1) (mod m),3=16(mod 13),或者说,16-3能够被13整除——同余法就是这意思了。说“乘同余法”,要点在于a*R(i)是乘法形式。在SAS系统中,a=397,204,094,m = 2^31-1=2147483647(一个素数),随机种子R(0)可以取1到m-1之间的任何整数。

    2.伪随机数的检验

    现在内置于各大软件包的随机数生成器都经过了彻底的检验,我们当然可以安心地使用这个黑箱。或则我们也可以多了解些。一个好的“伪”随机数列,应该看起来就跟从[0,1]上均匀分布的随机数列中随机抽取出来的一样。两个比较直观的检验有:

    1. 均匀性检验——所有的数是否均匀地分布在[0,1]区间;
    2. 独立性(或不相关)检验——序列中的数不存在序列相关,每个数都跟其前后出现的数独立无关。一个例子是如0.1、0.2、0.3、0.4、0.5,这里每一个相继的数都正好比它前面的一个数大0.1,这样这个数列就似乎有了某种“格局”。

    具体的检验方法就略之了,更多请参见徐钟济(1985)。

    3.生成其他概率分布的随机数

    前面提到过,我们把[0,1]上均匀分布随机变量的抽样值称为随机数,其他分布随机变量的抽样都是借助于随机数来实现的。对其他连续型分布,(累积)分布函数为F(x),r是在[0,1]上均匀分布的随机变量,另r=F(x),解出的x就是该连续型分布的随机抽样。由于x可以写成r的反函数,该变换就称作“逆变换”(Inverse Transformation method)。对离散型分布,思路类似,不过繁琐些,具体可以见徐钟济(1985)、Ross(2006)和Evans等(2001)。大多数相关软件包都会提供各种分布的随机数生成函数,知道有这回事就可以。最重要的,是我们陈述一类问题时,知道需要用哪种概率分布来描述Evans等(2001):

    常用的连续分布

    1. 均匀分布(Uniform Distribution)描述在某最小值和最大值之间所有结果等可能的随机变量的特性。
    2. 正态分布(Normal Distribution)是对称的,具有中位数等于均值的性质,就是我们熟悉的钟形形状。各种类型的误差常常是正态分布的。最后,作为中心极限定理的推论,大量具有任意分布的随机变量的平均数也呈正态分布。
    3. 三角形分布(Triangular Distribution)由三个参数来定义,最大值、最小值和最可能值。临近最可能值的结果比那些位于端点的结果有更大的出现机会。
    4. 指数分布(Exponential Distribution)常用来构建在时间上随机重现的事件的模型,如顾客到达服务系统的时间间隔、机器元件失效前的工作时间等。它的主要性质是“无记忆性”(Memoryless) ,即当前时间对未来结果没有影响。例如,不论机器原先已经运转了多长时间,它继续运转直至出现故障的时间间隔总有同样的分布。
    5. 对数正态分布(Lognormal Distribution):若随机变量x的自然对数是正态的,则x就服从对数正态分布。对数正态分布的常见例子是股票价格。

    常用的离散分布

    1. 贝努利分布(Bernoulli Distribution)描述的是只有两个以常数概率出现的可能结果的随机变量的特性;
    2. 二项分布(Binomial Distribution)给出每次试验成功概率为p的n次独立重复贝努利实验的模型;
    3. 泊松分布(Poisson Distribution)用于建立某种度量单位内发生次数的模型的一种离散分布,如某时段内某事件发生的次数。

    其他有用的分布如伽玛分布(Gamma Distribution)、威布尔分布(Weibull Distribution)、贝塔分布(Beta Distribution)略之。

    4.下期预告

    上次落了些东西。说到蒙卡与C++,在数量金融领域,不能不提的一本书就是Mark Joshi的C++ Design Patterns and Derivatives Pricing (Cambridge Uni. Press, 2004),面试中一定要说自己读过的,这样人家以为你至少会用C++做一个蒙特卡罗模拟。这个系列只讲SAS,下次先讲如何用SAS生成随机数,然后就是具体的项目。

    ---------参考资料---------

    1. Xitao Fan, etc..Monte Carlo Studies: A Guide for Quantitative Researchers. SAS Institute Inc.,2002
    2. Sheldon M. Ross. Simulation(4th Ed).Elsevier Inc..2006
    3. J.R.Evans等《模拟与风险分析》,洪锡熙译,上海:上海人民出版社,2001
    4. 徐钟济《蒙特卡罗方法》,上海:上海科学技术出版社,1985

    Comments (3)

    Please wait...
    Sorry, the comment you entered is too long. Please shorten it.
    You didn't enter anything. Please try again.
    Sorry, we can't add your comment right now. Please try again later.
    To add a comment, you need permission from your parent. Ask for permission
    Your parent has turned off comments.
    Sorry, we can't delete your comment right now. Please try again later.
    You've exceeded the maximum number of comments that can be left in one day. Please try again in 24 hours.
    Your account has had the ability to leave comments disabled because our systems indicate that you may be spamming other users. If you believe that your account has been disabled in error please contact Windows Live support.
    Complete the security check below to finish leaving your comment.
    The characters you type in the security check must match the characters in the picture or audio.

    To add a comment, sign in with your Windows Live ID (if you use Hotmail, Messenger, or Xbox LIVE, you have a Windows Live ID). Sign in


    Don't have a Windows Live ID? Sign up

    No namewrote:
    您需要二手液晶显示屏废旧液晶屏么?我们是不折不扣的二手液晶屏、旧液晶屏大批发商,长期大量供应可再利用的旧液晶屏。我公司提供的各种尺寸的二手液晶屏, 不同厚薄如笔记本屏,均已经过我们严格的分类,检验,测试流程。请访问协力液晶屏www.sceondhandlcd.com[ghfdabdjiaihhjc]
    Nov. 21
    Nov. 9
    Nov. 3

    Trackbacks

    The trackback URL for this entry is:
    http://johnthu.spaces.live.com/blog/cns!2053CD511E6D5B1E!369.trak
    Weblogs that reference this entry
    • None