リサンプリングとかブートストラップで統計量の分布を作って検定とか推定をするのに使えます。
/* ダミーデータの作成 */
data mult;
call streaminit(070605);
do y = 0, 1;
do i = 1 to 100;
x = rand('Normal') + y;
output;
end;
end;
run;
/* multtest */
/* outsampによりリサンプリングデータを出力 */
proc multtest data = mult outsamp = out
seed = 070605 nsample = 1000 permutation;
class y;
test mean(x);
contrast '0 vs 1' -1 1;
run;
/* 以下適当に整形 */
proc transpose data = out
out = out0(drop=_NAME_ _LABEL_) prefix = c0;
where _class_ = '0';
by _sample_; var x;
run;
proc transpose data = out
out = out1(drop=_NAME_ _LABEL_) prefix = c1;
where _class_ = '0';
by _sample_; var x;
run;
/* 平均値の和を統計量とした場合 */
/* 実際はオリジナルの統計量を計算するときに使うと便利 */
data perm;
merge out0 out1; by _sample_;
stat = mean(c01 -- c0100) + mean(c11 -- c1100);
run;
goptions reset = all;
goptions ftext = 'Times New Roman' ftitle = 'Times New Roman'
htitle = 1.6 htext = 1.6 hsize = 15cm vsize = 15cm;
options linesize=120 pagesize = 300;
/* 平均値の和の分布 */
proc univariate data = perm;
var stat; histogram;
run;
/* データから計算される統計量 */
proc means data = mult;
class y;
var x;
run;
multtest procedureのoutsampオプションでリサンプリングとかブートストラップしたデータセットを出してくれるところがポイント。
outsampオプションはデータセット(HDD)に書き出してくれちゃうので、数万回程度のリサンプリングにもすごい時間がかかる、しょうがないっちゃしょうがないんだけどね。
まあRを使った方が早くて幸せでしょう。
