t 検定と分散分析

2標本の t 検定と2群しかないDataでの一元配置分散分析のP値が一致することは有名な話で、それは F 統計量が t 統計量の2乗にあたるからである。

計算機統計の基本として、P値、例数設計と検出力、信頼区間の算出が基本ﾃｸﾆｯｸであるが、 t 検定の P 値を R で1000回シミュレーションするコードと、SAS で 1000 回シミュレーションするコードをつくったのだが、SAS での結果がおかしい。SAS では PROC GLM を用いて PROC TTEST を用いなかったのだが、まさかこれが？と一瞬疑ってしまった。

調べてみたら、Data生成のときの単純なコーディングミスで一安心した。output を 2回書いていたので、群内でつねに重複するDataが発生していた。こうなると群内平均平方が小さくなってしまうので有意差がつきやすくなってしまう。1000回のシミュレーションで、150回くらい有意になってしまっていた。

計算機シミュレーションでは小さいミスでも結論が異なってしまうことを、いまさらながら実感する。

R のコードはこれ。

Sim.ttest <- function(rep) {
alpha <- 0
for (i in 1:rep) {
   A <- rnorm(6, mean=1.5, sd = 0.15)
   B <- rnorm(6, mean=1.5, sd = 0.15)
   D <- data.frame(x = A, y = B)
   out <- t.test(x = D$x, y=D$y, var.equal=TRUE)
   if (out$p.value < 0.05) {
   alpha <- alpha + 1
  }
}
return(alpha)
}

SAS のコードはこれ。

Data d1;
array mu(2) mu1-mu2;
mu1 = 1.5; mu2 = 1.5;
do i = 1 to 1000;
  do group = 1 to 2;
  do j = 1 to 6;
   if group = 1 then
    y = RAND('NORMAL', 0, 0.15) + mu1;
  if group = 2 then
   y = RAND('NORMAL', 0, 0.15) + mu2;
   output;
  end;
end;
end;
RUN;

ods NORESULTS;
ods listing close;
ods output OverallANOVA = out;
PROC GLM data=d1;
by i;
class group;
model y = group;
RUN;

DATA result;
set Out;
if (ProbF < 0.05) then flag = 1;
else                 flag  =0;
if (ProbF > 0) then output;
RUN;

ods listing;
PROC FREQ data=result;
tables flag / norow nopercent;
RUN;
QUIT;