我偶爾會去逛一個叫做「統計生活館」的網站。這個大概是我覺得國內最完整的一個在介紹統計學的網站。裡面有一個討論區,經常有人會上去問些問題。如果我知道答案,也會在那邊幫點小忙。



前幾日有位仁兄問了一個問題。他用亂數生成的方法模擬了一個資料庫包含1068筆資料,然後用 SAS 的 PROC CLUSTER 做 Cluster analysis,最後再把樹狀圖畫出來。這個流程其實是很簡單的,但是他要重複上述的動作 10000 次。他問程式該怎樣寫。



我本來想說用個 Macro 來搞定,但是卡在怎樣在亂數模擬變數那一階段請 SAS 生出 10000 個資料庫,因為 Do loop 不能掛在 Data procedure 的外面。後來想到一個簡單的方法,不用 Macro 就可以搞定。



原始程式如下:



data simulation;

do no=0001 to 1068;

array Q[*] Q01-Q25;

do j=1 to 25;

k=rand('uniform');

Q(j)=ranbin(1,1,k );

end;

output ;

end;

proc distance data=simulation method=djaccard absent=0 out=distjacc;

var anominal(Q01-Q25);

run;

proc cluster data=distjacc method=ward rsq rmsstd pseudo outtree=tree;

var _numeric_;

run;

proc tree data=tree out=out;

run;






修改後的程式如下:



data simulation;

do index=1 to 10000;

do no=0001 to 1068;

array Q[*] Q01-Q25;

do j=1 to 25;

k=rand('uniform');

Q(j)=ranbin(1,1,k );

end;

output ;

end;

end;

run;



proc distance data=simulation method=djaccard absent=0 out=distjacc;

by index;

var anominal(Q01-Q25);

run;



proc cluster data=distjacc method=ward rsq rmsstd pseudo outtree=tree;

by index;

var Dist1-Dist1068;

run;



proc tree data=tree;

by index;

run;


差別主要在生成資料那一步驟,只需要在外面多掛一個迴圈,讓裡面生亂數的迴圈再跑 10000 次,然後資料裡面會多了一個叫做 index 的變數。這個變數從 1 到 10000,就等於是我把 10000 個資料庫全部放在一起。之後的任何 Procedure 就可以用 by statement 來搞定。稍微需要注意的是,在 PROC CLUSTER 的 Var statement 中不能再用通用變數 _numeric_,而要改成Dist1-Dist1068。原因是若用 _numeric_ 通用變數的話,會把 index 也放進去分析,這樣會導致距離矩陣不是 n x n,而是 (n+1) x n。目前唯一還不算解決的是最後的 Proc tree 程序,不能再使用 out option,我想這是可以解決的,只是懶得再想而已。:P














-----
arrow
arrow
    全站熱搜

    cchien 發表在 痞客邦 留言(0) 人氣()