真是够乌龙的 ... ※ 引述《cjlu (Joe)》之铭言： > ※ 引述《[email protected] (羊咩咩)》之铭言： > > 演算法的步骤为 > > 1. 产生 U~uniform(-1,1) > > V~uniform(-1,1) > > 2. let W=U^2 + V^2 > > 3. If W>1 , 回到步骤一重新来一次 > > 这样的话就可以产生随机变数 服从常态分配了 这个程序产生的根本就不是常态分配的随机乱数! 竟然一时不察, 还写了下述的程式, 那只能纯视为是一种程式设计的练习! > > 可是在写程式的时候 要如何写 > > 才能够在一个function内产生n个服从常态分配的样本呢?? > > 对於软体方面还是新手 不熟悉该怎麽写回圈...>< > > 不写回圈的话 一次只能产生一个样本了...Orz > ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ > 若是为了程式设计的练习 [但愿不是], > 下面这个被包装 (wrap-up) 为 function 的可供参考: > rNorm <- function(n=10) { > z <- numeric(n) > k <- 1 > while(k <= n) { > w <- sum(runif(2, -1, 1)^2) > if(w <= 1) { > z[k] <- w > k <- k + 1 > } > } > z > } > 例如: > > rNorm(5) > [1] 0.4193 0.9392 0.2459 0.7765 0.2454 > > rNorm(10) > [1] 0.5804 0.2529 0.7186 0.6608 0.1869 0.6637 0.3839 0.6172 0.5858 0.4746 > 只是这麽作非常没有效率! > 用下面这个小 function 去比较 nNorm() 与 R 内设提供的 rnorm() > 在产生较大数量的样本数时其耗费的时间差 (单位是秒) > compN2n <- function(n) { > system.time(rNorm(n))[3] - system.time(rnorm(n))[3] > } > > compN2n(1000) > [1] 0.03 > > compN2n(5000) > [1] 0.13 > > compN2n(10000) > [1] 0.25 > > compN2n(50000) > [1] 1.2 > 若只需要产生一次大小为 n 的常态样本, 那也就罢了; > 但若这是为一个重复执行的模拟所需, 那就很没效率了. > Accept-reject algorithm 是个蛮重要的方法, > 但是去练习这方法的题目或许应可更具有实用价值会比较好. -- 夫兵者不祥之器物或恶之故有道者不处君子居则贵左用兵则贵右兵者不祥之器非君子 之器不得已而用之恬淡为上胜而不美而美之者是乐杀人夫乐杀人者则不可得志於天下 矣吉事尚左凶事尚右偏将军居左上将军居右言以丧礼处之杀人之众以哀悲泣之战胜以 丧礼处之道常无名朴虽小天下莫能臣侯王若能守之万物将自宾天地相合以降甘露民莫 之令而自均始制有名名亦既有夫亦将知止知止可以不殆譬道之在 140.116.152.108海