myscratchpad

Van Selst, M., & Jolicoeur, P. (1994). A solution to the effect of sample size on outlier elimination. The quarterly journal of experimental psychology, 47A, 631-650.
より。
結論としては基準変動型のいずれかがお勧めらしい


## 仮想反応時間データ生成
rt <- round(sample(rnorm(10000, mean=300, sd=10)+(300*rexp(10000)), 20))
hist(rt)
win.graph(); plot(rt)
mean(rt); sd(rt)
range(rt)

# 基準変動型手続きに使うデータ
ssz <- c(4,5,6,7,8,9,10,11,12,13,14,15,20,25,30,35,50,100) # 抽出するサンプルの数
nrv <- c(1.458,1.68,1.841,1.961,2.050,2.12,2.173,2.22,2.246,2.274,2.31,2.326,2.391,2.410,2.4305,2.45,2.48,2.50) # 基準変動型非再帰手続き用
mrv <- c(8, 6.20, 5.30, 4.80, 4.475,4.25,4.11,4.00,3.92,3.85,3.80,3.75,3.64,3.595,3.55,3.54,3.51,3.50) # 基準変動型修正再帰手続き用
## 表にしてみる
tb <- rbind(nrv, mrv)
colnames(tb) <- ssz
tb

## 線形補間で予測値を得る関数を生成
nrvfun <- approxfun(ssz, nrv) # 基準変動型非再帰手続き用
mrvfun <- approxfun(ssz, mrv) # 基準変動型修正再帰手続き用


## 基準変動型非再帰手続き (Non-recursive with moving criterion)
mrt <- mean(rt)
sdrt <- sd(rt)
sszrt <- length(rt)
crco <- nrvfun(sszrt) # 線形補間関数による基準値
cv <- c(mrt-(crco*sdrt), mrt+(crco*sdrt))
rt2 <- rt[cv[1]<=rt & rt<=cv[2]]
mean(rt); sd(rt); length(rt)
mean(rt2); sd(rt2); length(rt2) # 最終出力
range(rt); range(rt2); cv; crco


## 基準変動型修正再帰手続き (Modified recursive with moving criterion)
rttmp <- rt[rt<max(rt)]
mrttmp <- mean(rttmp)
sdrttmp <- sd(rttmp)
sszrt <-  length(rt)
crco <- mrvfun(sszrt) # 線形補間関数による基準値
cv <- c(mrttmp-(crco*sdrttmp), mrttmp+(crco*sdrttmp))
minv <- min(rt); maxv <- max(rt)
ucv.l <- cv[1]; ucv.u <- cv[2]
svct <- rt
    while((minv < ucv.l | maxv > ucv.u)&(length(svct)>3)) {
        if (minv < ucv.l) {svct<- svct[svct>minv]}
        if (maxv > ucv.u) {svct<- svct[svct<maxv]}
        cat(ucv.l, ucv.u, crco, "\n") # 最後に使われた基準値とカットオフポイントを返す
        if (length(svct) < 2) break
        maxv <- max(svct)
        minv <- min(svct)
        svct2 <-svct[svct<maxv]
        if (length(svct2) < 2) break
        ucv.l <- mean(svct2)-(crco*sd(svct2))
        ucv.u <- mean(svct2)+(crco*sd(svct2))
    }
rt2 <- svct
mean(rt); sd(rt)
mean(rt2); sd(rt2) # 最終出力
cv; crco
range(rt); length(rt)
range(rt2); length(rt2)

## ハイブリッド手続き (Hybrid procedure)
crco <- 2.5; crco.nr <- crco # 基準値は勘で
## まず非再帰手続きを行う
mrt <- mean(rt)
sdrt <- sd(rt)
sszrt <- length(rt)
cv <- c(mrt-(crco*sdrt), mrt+(crco*sdrt)); cv.nr <- cv
rt.nrc <- rt[cv[1]<=rt & rt<=cv[2]]
## 修正再帰手続きを行う。この際、非再帰手続きの基準値に1プラスする
crco <- crco.nr+1; crco.mr <- crco
rttmp <- rt[rt<max(rt)]
mrttmp <- mean(rttmp)
sdrttmp <- sd(rttmp)
cv <- c(mrttmp-(crco*sdrttmp), mrttmp+(crco*sdrttmp)); cv.mr <- cv
minv <- min(rt); maxv <- max(rt)
ucv.l <- cv[1]; ucv.u <- cv[2]
svct <- rt
    while((minv < ucv.l | maxv > ucv.u)&(length(svct)>3)) {
        if (minv < ucv.l) {svct<- svct[svct>minv]}; cv.mr[1] <- ucv.l
        if (maxv > ucv.u) {svct<- svct[svct<maxv]}; cv.mr[2] <- ucv.u
        cat(ucv.l, ucv.u, crco, "\n") # 最後に使われた基準値とカットオフポイントを返す
        if (length(svct) < 2) break
        maxv <- max(svct)
        minv <- min(svct)
        svct2 <-svct[svct<maxv]
        if (length(svct2) < 2) break
        ucv.l <- mean(svct2)-(crco*sd(svct2))
        ucv.u <- mean(svct2)+(crco*sd(svct2))
    }
rt.mrc <- svct
## 非再帰手続の平均と修正再帰手続の平均を平均する
mean(c(mean(rt.nrc),mean(rt.mrc))) # 最終の代表値
mean(rt); sd(rt); length(rt)
mean(rt.nrc); sd(rt.nrc); length(rt.nrc); cv.nr; crco.nr
mean(rt.mrc); sd(rt.mrc); length(rt.mrc); cv.mr; crco.mr
range(rt); range(rt.nrc); range(rt.mrc)