R Language
順列テストの実行

かなり一般的な関数

ビルトインの歯の成長データセットを使用します。モルモットにオレンジジュースとビタミンCを与えた場合、歯の成長に統計的に有意な差があるかどうかに興味があります。

完全な例は次のとおりです。

teethVC = ToothGrowth[ToothGrowth$supp == 'VC',]
teethOJ = ToothGrowth[ToothGrowth$supp == 'OJ',]

permutationTest = function(vectorA, vectorB, testStat){
  N = 10^5
  fullSet = c(vectorA, vectorB)
  lengthA = length(vectorA)
  lengthB = length(vectorB)
  trials <- replicate(N, 
                      {index <- sample(lengthB + lengthA, size = lengthA, replace = FALSE)
                      testStat((fullSet[index]), fullSet[-index])  } )
  trials
}
vec1 =teethVC$len;
vec2 =teethOJ$len;
subtractMeans = function(a, b){ return (mean(a) - mean(b))}
result = permutationTest(vec1, vec2, subtractMeans)
observedMeanDifference = subtractMeans(vec1, vec2)
result = c(result, observedMeanDifference)
hist(result)
abline(v=observedMeanDifference, col = "blue")
pValue = 2*mean(result <= (observedMeanDifference))
pValue

CSVを読み込んだ後、関数を定義します

permutationTest = function(vectorA, vectorB, testStat){
  N = 10^5
  fullSet = c(vectorA, vectorB)
  lengthA = length(vectorA)
  lengthB = length(vectorB)
  trials <- replicate(N, 
                      {index <- sample(lengthB + lengthA, size = lengthA, replace = FALSE)
                      testStat((fullSet[index]), fullSet[-index])  } )
  trials
}

この関数は2つのベクトルを取り、それらの内容を一緒にシャッフルし、シャッフルされたベクトルに対して関数testStatを実行します。 teststatの結果が戻り値であるtrialsに追加されます。

これはN = 10^5回行う。値Nは関数のパラメータである可能性が非常に高いことに注意してください。

これにより、新しいデータセット、 trials 、実際には2つの変数間に関係がない場合に生じる可能性のある一連の手段が残されます。

次に、テスト統計を定義します。

subtractMeans = function(a, b){ return (mean(a) - mean(b))}

テストを実行します。

result = permutationTest(vec1, vec2, subtractMeans)

我々の実際の観測された平均差を計算する：

observedMeanDifference = subtractMeans(vec1, vec2)

我々のテスト統計のヒストグラム上での我々の観察の様子を見てみましょう。

hist(result)
abline(v=observedMeanDifference, col = "blue")

私たちの観察された結果がランダムなチャンスによって非常に起こりそうなようには見えません...

2つの変数の間に関係がない場合、元の観測結果の可能性であるp値を計算します。

pValue = 2*mean(result >= (observedMeanDifference))

それをちょっと打ち砕いてみましょう：

result >= (observedMeanDifference)

次のようなブール値ベクトルを作成します。

FALSE TRUE FALSE FALSE TRUE FALSE ...

resultの値がobservedMean以上になるたびにTRUEます。

ファンクションmean 、このベクトルを解釈する1についてTRUEと0のためFALSE 、そして私たちの割合与える1ミックスの年代を、回数たちのシャッフルベクトルが違いを意味する。すなわち突破か、我々が観察するもの等しかったです。

最後に、テスト統計の分布が非常に対称的であるために、2を掛けます。実際に観測された結果よりも "極端な"結果を知りたいのです。

残っているのは、 0.06093939と判明したp値を出力すること0.06093939です。この値の解釈は主観的ですが、私はビタミンCがオレンジジュースよりもかなり多く歯の成長を促進するように見えると言います。