STATA入门10随机模拟

《STATA入门10随机模拟》由会员分享，可在线阅读，更多相关《STATA入门10随机模拟（29页珍藏版）》请在装配图网上搜索。

1、STATA 入门 10 随机模拟只要你自己试试模拟随机现象几次，就会加强对概率的了解，比读很 多页的数理统计和概率论的文章还有用。学习模拟，不仅是为了解模拟本 身，也是为更了解概率而了解模拟。10.1 伪随机数生成(0，1)之间均匀分布的伪随机数的函数为 uniform()diuniform()diuniform()diuniform()每次都得到一个大于零小于1 的随机数。如果要生成一位数的随机数(即 0，1，2，3，4，5，6，7，8，9) 可以取小数点后第一位数，通常用下面的命令diint(10 某 uniform()两位随机数(0-99)则取小数点后两位小数，即也可以同时生成多个随机数

2、，然后将该随机数赋给某个变量。要注意 的是，电脑中给出的随机数不是真正的随机数，而是伪随机数，因为它是 按照一定的规律生成的。如果给定基于生成伪随机数的初始数值(即 eteed#)，则对相同的初始数值，生成的伪随机数序列完全一样。某 =begin= =clearetob10gen 某 1=uniform()gen某2二uniform。/注意到某1与某2不一样eteed1234geny1=uniform()eteed1234geny2=uniform()geny3二uniform。/注意到y1与y2 一样，但均与y3不同eteed5634genz1=uniform()eteed1234genz2

3、二uniform。/注意到 z2 与 y1,y2 样，但 z1 与 z2 不同lit某=end=10.2 简单模拟一旦有了可靠的概率模型，模拟是找出复杂事件发生概率的有效工具 一个事件在重复结果中发生的比例，迟早会接近它的概率，所以模拟可以 对概率做适当的估计。例 1：如何执行模拟掷一枚硬币10次，结果中会出现至少3个连续正面或者至少3个连 续反面的概率是多少？思考：（1）猜想这个概率大约是多少？（2）如何从理论上计算出这个概率？（3）如何模拟计算这个概率？第一步：提出概率模型。每一次掷，正面和反面的概率各为 0.5投掷之间，彼此是独立的。也就是说，知道某一次掷出的结果，不会改变任何其他次所掷

4、结果的概率。第二步：分配随机数字以代表不同的结果。随机数字表中的0-9每个数字出现的概率都是 0.1每个数字模拟掷一次硬币的结果。奇数代表正面，偶数代表反面。第三步：模拟多次重复。生成 10 个随机数字记录开心的事件（至少连续三个正面或反面）是否发生，如果发生， 记为1，否则为0 重复 10 次（或者100，1000，1000000次），计算概率 =开心事件发生/总重复次数。真正的概率是 0.826。大部分的人认为连续正面或反面不太容易发生。但模拟结果足以修正我们 直觉错误。某begincaptprogramdropeq3programeq3,rcla/rcla 选项表示计算结果将由 retu

5、rn 返回到 r()verion9drop_all/清空所有数据，不能用clearetob10/将生成10个观察值tempvar某yz/设定某,y,z为临时变量gen某二int（ 10某uniform。）/产生10个随机变量，可能为0, l,，9gen、y=（mod（、某,2）=0）/如果生成的随机变量为奇数，则 y=0;为偶数，y=1genz=0/生成 Z=0forvaluei=3/10replacez=1ify二二 y_n-1&y二二y _n_2ini，/连续 三个变量相等时 z=1umzreturncalarma某二r（ma某）/z取1表示高兴的事发生（三连续）， 否则失败 endimu

6、latema 某=r（ma 某）,rep（10000）nodot: eq3/重复 1 万次，取平均 结果 um某=end=由于上述命令要不停生成变量，进行变量代换，所以计算速度较慢如果使用矩阵，则计算速度会大大加快，命令也更简捷。某 MATA某=begin= =mataA二uniform(10000,10):0.5/每行为一次试验,每列为某次试验中硬 币的正反面结果B二J(10000,l,0)/记录每次试验的结果，先记为0.若高兴 结果出现再改为lfor(j=l;j二row(A);j+)/第j次试验for(i=3;i二col(A);i+)/第j次试验中第i次硬币出现结果if(Aj,(i-2,i

7、-1,i)=(0,0,0)|Aj,(i-2,i-1,i)=(1,1,1)Bj,1=1/若连续为正面或者连续为反面,则记为1breakmean(B)/求开心事件的次数end某=end=10.3 复杂模拟例 2：我们要女儿任务：一对夫妇计划生孩子生到有女儿才停，或者生了三个就停，他 们拥有女儿的概率是多大？思考：理论上，概率是多少？第一步：概率模型每一个孩子是女孩的概率是 0.49，且各个孩子的性别是互相独立的。 第二步：分配数字00，01，02，。，49=女孩49，50，51，。，99=男孩第三步：模拟 从随机数表中生成一对一对的数字，直到有了女儿，或已有 3 个孩子 重复 100 次。6905

8、1648178717648987男女女女女男女男男男 用数学可以计算出来，有女孩的真正概率是 0.867某=begin=captprogramdropgirlprogramgirl,rcladrop_alletob3gen 某二int（100 某 uniform。）geny=（某49）/生女孩y=l,生男孩，y=0umyreturncalarma某二r（ma某）/若有一个女孩，则ma某取1,若三个 全为男孩，取 0endimulatema 某二r（ma 某），rep（10000）nodot: girlum某=end=更快的模拟方式.某=begin= =captprogramdropgirlpr

9、ogramgirlmatA=matuniform（1,3）calargirl=1ifA1,10.49&A1,20.49&A1,30.49calargirl=0/若三个全为男孩，生女孩数为0endimulategirl,rep(10000)nodot:girltab_im某=end= 某=begin=mataA=uniform(10000,3):0.49B=J(10000,1,1)for(i=1;i=row(A);i+)if(Ai,.=(0,0,0)Bi,l=0/若三个全为男孩，生女孩数为0mean(B)end某=end=10.4 多阶段模拟例 3：肾脏移植思考：计算理论概率第一步：采用概率树将

10、信息组织起来第二步：分配数字阶段 1：0=死亡1-9=存活阶段 2：0-5=移植成功6-9=仍需洗肾阶段 3，成功0-6=存活五年7-9=死亡阶段 3：洗肾0-4=存活五年5-9=死亡第三步：模拟数学计算结果为 0.558。begin/某换肾后的五年存活率：撑过手术牛 0.9，术后存活的人中有 0.6 移植成功，0.4 还得回去洗肾；五年后，有新肾的人70%仍然活着，而洗 肾的只有一半仍活着。计算换肾的人活过五年的概率。某/captprogramdropurvprogramurv,rcladrop_alletob1genz=1gen 某 1=int(10 某 uniform()if 某1=0r

11、eplacez=0elegen 某 2=int(10 某 uniform()if 某26replacez=0elegen 某 4=int(10 某 uniform()if 某44replacez=0umzreturncalarma 某二r(ma 某)endimulatema 某二r(ma 某)，rep(10000)nodot: urvum某=end=更简捷的方式begin=captprogramdropurvprogramurvcalarz=1ifuniform()0.1calarz=0eleifuniform()=0.7calarz=0eleifuniform()=0.5calarz=0en

12、dimulatez,rep(10000)nodot:urvum某=end=10.5 商店案例任务：设某种商品每周的需求量某是 10,30上均匀分布的随机变量 而某店进货数量为10,30中的某一整数。商店每销售一单位商品可获利500 元；若供 大于求则削价处理，每处理一单位商店亏损100 元；若供不应求则可从外 部调剂供应，此时每一单位商店仅获利300 元。为使商店所获利润期望值 不少于 9280元，试确定最少进货量。解：由题设，随机变量某的概率密度为110某30f某20其它0设商店进货量为 a 则商店利润 Z 为500 某 100a 某 Z500a300 某 a600 某 100a300 某

13、200a10某aa某30EZZ某f某d某0d 某 10301Z 某 d 某 0d 某1a600 某 100ad 某 102022300a200 某 d 某 20aa12300 某 100a 某 102012150 某 200a 某207.5a2350a525092807.5a2350a403003a622.5a6502026265a26332.5故知最小进货量为21 个单位可使期望利润大于9280元。某=end=captuprogramdropgoodprogramgoodcalar 某=10+int（20 某 uniform。） ifl=某calarz=500 某某-100 某（，-某）el

14、ecalarz=500 某1，+300 某（某-，）endetmoreoffquietlyforvaluei=10/30imula tez,rep(1000)nod ot: goodiquietlyumcalarzi二r(mean)calarlit 某 =end=10.6 练习亚洲随机甲虫的命运(Aiantochaticbeetie)，这种昆虫的繁殖模式 是：(1)20%的虫还没有生雌幼虫之前就死掉了，30%生 1 只雌虫，50%生 2 只雌虫。(2)个别雌虫的繁殖情况互相独立。问：亚洲随机甲虫的前 途会：繁殖很快？勉强保持数目？逐渐灭绝？生日问题：只要一间屋里有 23 个人，则至少有两人同一

15、天生日的概 率会超过 1/2。试模拟两个班 60 名学生同一天过生日的概率，并用你们 班，或者前几届后几届班组验证之。古罗马时代的赌博：掷 4 块绵羊距骨是古罗马最受欢迎的赌博，掷多 次后骨头的四个面挖概率如下：窄而平的一面 0.1 宽而凹的一面 0.4宽而凸的一面 0.4 窄而凹的一面 0.1掷 4 块距骨最好的结果叫“维纳斯”，这是朝上的四个面都不一样的 情况，估计掷出“维纳斯“的概率。中国会有多少人成为可怜的单身汉？不少中国人（尤其是农村）有重 男轻女思想，他们总是想尽办法生男孩，性别失调可能会导致越来越多的 可怜的男性单身汉。假设所有夫妇都直到生出男孩，或者生够3 个孩子才 停止生育，

16、则一个家庭平均会有几个孩子？多少个男孩？多少个女孩？商店的平均获利：一商店经销某种商品每周进货的数量某与顾客对该 种商品的需求量Y是相互独立的随机变量，且都服从区间10,20上的均 匀分布。商店每售出一单位商品可获得利润1000元，若需求量超过了进 货量，商店可从其它商店调剂供应，这时每单位商品可获利润500元，试 计算商店经销该种商品所获利润的期望值。f某,yf某某fYy110 某 2010010y20 其它 0设随机变量Z表示商店每周所获利润，则1000YZ1000 某 500Y 某Y 某 1000 某 500 某 YY 某 EZZ 某 f 某,yd 某 dy11000y 某 d 某 dy

17、500 某 yd 某 dy100D1D2202220y10dyyd 某dy5 某 yd 某10y10102010y20ydy102035y210Y50dy1022011010y2y331015y35y250y210200002010.7 附录亚洲随机甲虫的命运（Aiantochaticbeetle），这种昆虫的繁殖模式 是：（1）20%的虫还没有生雌幼虫之前就死掉了，30%生 1 只雌虫，50%生2 只雌虫。（2）个别雌虫的繁殖情况互相独立。问：亚洲随机甲虫的前途会：繁殖很快？勉强保持数目？逐渐灭绝？某=begin= =captureprogdropbbprogbb,rclalocalbeet

18、le=1whil ebee tle500&bee tl e0drop_alletobbeetletempvar 某 ygen 某二uniform。gen、y=1replace、y=2if、某0.8quietlyumylocalbeetle=r(um)noiilydiplayaerrorbeetleendetmoreoffquietlybb10=end=生日问题：只要一间屋里有23 个人，则至少有两人同一天生日的概率会超过1/2。试模拟两个班60 名学生同一天过生日的概率，并用你们 班，或者前几届后几届班组验证之。某 =begin=captuprogdropbirth progbirth dro

19、p_all etob60 tempvarygeny二int（365 某 uniform。） ortycalarz=0 forvaluei=1/59calarz=1continue,breakendimulatebirthz,rep(100)um某=end=古罗马时代的赌博：掷4 块绵羊距骨是古罗马最受欢迎的赌博，掷多 次后骨头的四个面挖概率如下：窄而平的一面0.1宽而凹的一面 0.4宽而凸的一面0.4窄而凹的一面 0.1掷4 块距骨最好的结果叫“维纳斯”，这是朝上的四个面都不一样的 情况，估计掷出“维纳斯“的概率。某=begin=captuprogdropwnprogwndrop_alleto

20、b4tempvar 某 ygeny二uniform。recodey(0/0.1=l)(0.1/0.2=2)(0.2/0.6=3)(0.6/l=4),gen( 某) ort某，calarz=1forvaluei=1/4if 某，、i，=、某，、i，+lcalarz=0continue,breakendimulatewnz,rep(1000)um某=end=中国会有多少人成为可怜的单身汉？不少中国人（尤其是农村）有重男轻女思想，他们总是想尽办法生男孩，性别失调可能会导致越来越多的可怜的男性单身汉。假设所有夫妇都直到生出男孩，或者生够3 个孩子才停止生育，则一个家庭平均会有几个孩子？多少个男孩？多少

21、个女孩？某=begin=captuprogdropchildprogchilddrop_alletob3tempvarygeny二int（100 某 uniform。）calarboy=0calargirl=0ifyl49calargirl=1ify249calargirl=2/第二个为女孩ify349calargirl=3/三女elecalarboy=1/2 女 1 男elecalarboy=1/1 女 1 男elecalarboy=1/0 女 1 男endimulatechildboy=boygirl=girl,rep(1000)um某=end=如果不是重男轻女，而是生到三个为止，不论是男

22、是女，则begincaptuprogdropchildprogchild,rcladrop_alletob3tempvarybggeny二int(100 某 uniform。)genb=0geng=0replace、b=1if、y=49umbcalarboy=r(um)replace、g=1if、y49um gcalargirl=r(um)endimulatechildboygirl,rep(1000) umend商店的平均获利：一商店经销某种商品每周进货的数量某与顾客对该 种商品的需求量Y是相互独立的随机变量，且都服从区间10,20上的均 匀分布。商店每售出一单位商品可获得利润1000 元，

23、若需求量超过了进 货量，商店可从其它商店调剂供应，这时每单位商品可获利润500 元，试 计算商店经销该种商品所获利润的期望值。当某和Y取连续变量时，某二10, 20，Y=10，20某=begin=captureprogramdroppfprogpfcalar 某=10+10 某 uniform()calary=10+10 某 uniform()if 某=ycalarz=1000 某 yendimulatepfz,rep（10000）um 某=end=结果约为：14153.12 当进货量和销售量均取整数时，即某=10，11，。20；而 y=10，11，。，20 某=begin=capturepr

24、ogramdroppfprogpfcalar 某=10+int（ll 某 uniform。） calary=10+int（ll 某 uniform。） if 某=ycalarz=1000 某 yendimulatepfz,rep(10000)um 某=end= 结果约为：14102.65clearmatA=J(11,11,.)forvaluei=1/11forvaluej=1/11if、i二、jmatAi,j=1000 某(9+i)eleifijmat A、i,、j=500 某(18+i+j)matlitAmata=J(1,11,1)matC=a 某 A 某 a/121matlitC某 15909.091某当进货量和销货量取1 位小数时clearmatA=J(101,101,.)forvaluei=1/101forvaluej=1/101if、i二、jmatA、i,、j=1000 某(10+(iT)/10)eleifijmatA、i,、j=500 某(20+(i-l)/10+(j-l)/10)mata=J(1,101,1)matC=a 某 A 某 a/10201matlitC15841.584某=end=请思考，上述程序有什么问题？为什么答案有如此大的差距？