外れ値を含むデータに対する効果量：正規化四分位区間と中央値を使ったほんのちょっとだけ頑健な方法

背景

こんなデータを作ります。対応なしの2群のデータ（x, y）です。

xは，平均50，標準偏差10の正規分布に従う乱数39個です。

yは，平均55，標準偏差10の正規分布に従う乱数39個です。

なので，理論的には，こんな感じになり，

効果量（標準化平均差）は0.50くらいになるはずです。

ただし，xには，値が0というケースを一個混ぜます（n = 40になります）。

ただし，yには，値が100というケースを一個混ぜます（n = 40になります）。

つまり，外れ値を紛れこませるわけです。

これで標準化平均差をもとめてみます。

Rだと，

x<-c(rnorm(39,50,10),0)

y<-c(rnorm(39,55,10),100)

こんな感じでしょうか。

標準化平均差は，今回は⊿でいいでしょう。xを統制群とみて，

delta<-(mean(y)-mean(x))/sd(x) 

 としましょう。ここでは不偏標準偏差推定値を使っています。

このように，各変数ひとつずつでも外れ値が紛れ込んでいると，特に小標本の場合，標準化平均差は大きく違った値を取ってしまいます。

ノンパラメトリック・ブートストラップ法（B = 1,000，パーセンタイル法）で95%信頼区間を構成すると，

delta<-numeric(0)
for(i in 1:1000){
subx<-sample(20,x,replace=T)
suby<-sample(20,y,replace=T)
delta[i]<-(mean(suby)-mean(subx))/sd(subx)
}
quantile(delta,c(.025,.975))

手元のシミュレーションだと，95% CI [0.05, 1.20]ってかんじでした。点推定値は.59。うわ～。

正規化四分位区間と中央値を使った方法

こういうデータに使うノンパラメトリックな代用品はほんとうにいろいろあるのだけど，一番簡単，ストレートなのは，順序統計量とかで代用する方法ね。

平均値は中央値で代用しよう。

標準偏差は正規化四分位区間で代用しよう。

正規化四分位区間とは，四分位区間（IQR）に係数0.7413をかけたもの。

標準正規分布に従う変数の四分位区間は1.349なので，この係数は1/1.349ってことね。

これで標準偏差とおなじスケールになるね。

なので，⊿のノンパラメトリックな代用としては，

(yの中央値-xの中央値)/xの正規化四分位区間

っていうことね。きれいな正規分布に従うデータなら，標準化平均差と基本的に同じ値がでるはず。

計算するには，

delta2<-numeric(0)
for(i in 1:1000){
subx<-sample(x,20,replace=T)
suby<-sample(y,20,replace=T)
NIQRx<-(quantile(x,.75)-quantile(x,.25))*.7413
delta2[i]<-(median(suby)-median(subx))/NIQRx
}
quantile(delta2,c(.025,.975))

こんなかんじね。

手元のシミュレーションだと，95% CI [-0.18, 1.08]という感じ。点推定値は，.49くらいだった。区間の幅は広いけど，数値例が外れ値によって効果量の値を大きくしようと設定されているのに対して，あきらかにパラメトリックな方法より下方に修正されている点で，まあいいかな。もちろん，歪んだ分布の効果量には，ブートストラップ法もBCa法がいいよ，とかそういう論文もあるし，パーセンタイル法じゃないほうがたいていはいいよね。

さらに，ノンパラメトリックな方法っていっぱいあるからね，マンホイットニーのU統計量から変換するやつとか。でもこの方法は標準化平均差と同じようなスケールで数字がでるからいいよね。わかりやすいよね。