用SAS模擬隨機數據 求PIE值

剛剛看到一本好書《統計模擬》作者叫羅斯[英文：Sheldon M. Ross. Simulation(4th Ed).Elsevier Inc..2006 ]. 顧名思義，這是一本描述怎么利用模擬一些符合統計學理論的數據，用途很廣，也就是說實際中的任何數據的分布都符合某種統計學模型，于是在沒有得到真實數據之前，數據分析師可以通過模擬數據來研究這些現實中的問題。如果通過模擬來研究未知問題，可以說得上是研究境界很高了?？偛荒苣玫揭恍嶋H數據，畫個好看的圖，就覺得自己可畫遍天下了吧。

由于自己不是統計出生，但是受過統計學老師的循循教誨，凡事從簡單開始。于是goolge了一下，當當中有這本書的中文介紹：

本書系統闡述了統計模擬的一些實用方法和技術。在對概率的基本知識進行了簡單的回顧之后，介紹如何利用計算機產生隨機數以及如何利用這些隨機數產生任意分布的隨機變量、隨機過程等。然后討論了一些分析統計數據的方法和技術。如Bootstrap（自助法）、方差縮減技術等。接著講述了如何利用統計模擬來判斷所選的隨機模型是否擬合實際的數據。最后介紹MCMC及一些最新發展的統計模擬技術和論題，如隨機序列函數和隨機子集函數的評估。本書在每章的最后還提供了不同難度的習題。本書可作為高等院校數學、統計學、科學計算、保險學、精算學等專業的教材，也可供工程技術人員和應用工作者參考。

一看有很多不懂的術語，頓時心生敬仰，后面寫著可以供“工程技術人員”參考，很顯然，我可以是這本書的讀者。全書近300頁，在今天剩下不多的時間內，很顯然，我決定不去看這本書，哪怕是一個字。還是老師的教導，從簡單開始。google到一個有趣的問題：用統計模擬計算圓周率pie值。Forcode提供了一種用excel求解pie的詳細過程，然后有人用Mathematica計算出來了。Hujiangtang很仔細的闡述了什么是隨機數？什么是蒙特卡羅模擬？為什么選擇SAS做蒙特卡羅模擬？SAS怎么做出來這些，還可以做哪些哪些分布…… 其中引用這個用統計模擬計算pie的例子和上述用excel和mathematica的方法，可是就是不提供解決這個問題的SAS代碼，在我看來，這是不可想象，我很仔細的找了半天，未果，于是我自己編了下列代碼，用SAS來實現模擬pie值。

借用別人寫的原理和圖來說明一下，

原理很簡單：

1）生成隨機數——生成n個均勻落在正方形內的點；

2）對落在正方形內的n個點，數一數正好落在圓里面的點的個數，假設為k（另外n-k個點就落在圓外面的正方形區域內）。數據分析師培訓

3）k/n就可以大致認為是圓的面積與正方形的面積之比，另其等于pai/4，就可以求出圓周率∏的估計值。n越大，算出來的pai-值越精確

SAS實現代碼：

%let num_s=10000;
data ex;
do i= 1 to &num_s.;
x=ranuni(0); y=ranuni(0);
area=sqrt(x**2+y**2);
output;
end;
run;

proc sort;
by area;
run;

data ex2;
set ex;
y1=.;y2=.;
num=_n_;
if area<=1 then y1=y ; else y2=y;
if area>1 and lag(area)<=1 then
do;
pie=(num*4)/&num_s.;
sas_pie=constant(‘pi’);
put pie= sas_pie=;
end;
run;
FILENAME file “c:\simulation.png”;
goptions reset=all hsize=8cm vsize=8cm  noborder device=png gsfname=file;;
symbol1 v=dot color=’red’ height=0.3;
symbol2 v=dot color=’blue’ height=0.3;
proc gplot;
plot  y1*x y2*x/overlay noaxis;
run;
quit;

模擬的點數為 100時pie=3；1000時，pie=3.068； 10000時， pie=3.1392；100000時pie=3.13543；1000000時，pie=3.141524，10000000時，pie=3.1418008。當然SAS也提供了pie的精確值供使用，使用函數constant就可以：sas_pie=constant(‘pi’); 值為：3.1415926536。上面的示例只是為了演示統計模擬的使用，其實SAS提供了很多種隨機函數，幾乎所有分布的隨機函數供大家使用，這里就不詳述了.CDA數據分析師培訓