例えば,2 群の平均値の差の検定を行うとき,データ数と平均値と標準偏差しかわかっていないときに,検定はできないだろう(あるいは,なんとかして,検定を行いたい)と考える人がいる。
2 群の平均値の差の検定は,t.test 関数を使えば行えるが,そのためには確かに,元データが必要である。
しかし,平均値の差の検定(いわゆる t 検定)の解説を見れば検定のために必要な情報は,各群のデータ数,平均値,不偏分散の,合計 6 個の情報だけでよいことがわかる。どんな分布に従っていようと無関係なのである。また,不偏分散は標準偏差を二乗するだけであるから,結局,「二群について,データ数と平均値と標準偏差がわかっていれば,平均値の差の検定は行える」ということになる。
このために,専用の関数を書くこともできるが,もっと汎用の方法がある。それは,ある群のデータをでっち上げる方法である。ただ単にでっち上げるのではなく,データ数,平均値,標準偏差が指定されたものと同じになるデータを作るのである。
データがどんな分布に従うかは問わないということを上に述べたが,まあ,せめて母分布が正規分布に従うデータを作成しよう。R には,正規乱数を生成する関数 rnorm がある。平均値 mean,標準偏差 sd の正規分布に従う n 個のデータ x を生成するには x <- rnorm(n, mean, sd) とするだけでよい。しかし,この段階では x の平均値と標準偏差は指定されたものとは異なる。なんといっても「乱数」なのであるから。正確に指定された平均値と標準偏差を持つ正規乱数を作るには適当に生成した正規乱数を一度標準化する。標準化した正規乱数は,平均値が 0,標準偏差が 1 になるというのは誰もが知っていることである。目的とする正規乱数はこの標準化された正規乱数に,指定された標準偏差を掛け,指定された平均値を加えればよいのである。そうすれば,できる正規乱数の平均値と標準偏差は,指定された平均値と標準偏差にぴたりと一致する。(ぴたりと一致するというところに注意が必要。つまり,指定された平均値と標準偏差は,示されたとおりの数値ではない,丸められた数値であるということ。そう言う点では,作られるデータは正確に元のデータと同じであるとは言えないということである)
例: データ数が 24 で,平均値が 123.45 標準偏差が 6.78 の正規乱数を作る。 > options(digits=20) > x <- rnorm(24) > cat("mean =", mean(x), " S.D.=", sd(x), "\n") mean = -0.2130160848930967 S.D.= 0.950709790266453 # 乱数なので毎回異なった値になるが,指定された値でないことは明らか > x <- scale(x) # 標準化を行う関数 scale を適用する > cat("mean =", mean(x), " S.D.=", sd(x), "\n") mean = -2.313868136446147e-17 S.D.= 1 # 平均値は -2.313868136446147×10の(-17)乗ということで,殆ど 0 ということ > x <- x*6.78+123.45 # 標準化の逆 > cat("mean =", mean(x), " S.D.=", sd(x), "\n") mean = 123.45 S.D.= 6.78 # 表示されているとおりの平均値と標準偏差を持つということである
さて,ここまで来れば答えはすぐそこ。各群ごとに上に述べた方法でデータをでっち上げ,そのでっち上げたデータを t.test 関数を使って検定すればよい。上の手順を二回使うし,他でも利用できそうな関数であるから,以下のように定義しておくとよい。
gendat1 <- function(n, mean, sd) return(scale(rnorm(n))*sd+mean)これを使って,以下のように検定を行うことができる。
例: データ数 平均値 標準偏差 A 群 24 123.45 6.78 B 群 35 126.28 5.43 > (A <- gendat1(24, 123.45, 6.78)) # A 群のデータ [1] 126.4254447481493 119.2677005077077 109.9809264790594 125.2044488778086 118.5402943426712 126.7757599329591 [7] 133.1987272298903 131.5923794221886 113.4644342955914 119.5887977878966 119.7263229865866 116.3732560196951 [13] 124.9521671141838 129.0740980418013 122.4442664845299 122.3567964723745 133.8439545000280 129.6773031494384 [19] 129.5823721125123 126.7703578980750 110.9939206252493 126.5980379505317 129.2924096537701 117.0758233673015 > (B <- gendat1(35, 126.28, 5.43)) # B 群のデータ [1] 139.9765197560334 132.1464449725871 136.6031115798909 125.9914781879804 123.5525306533870 130.6480113694447 [7] 125.2960720595223 124.8400514520517 126.2640875155946 124.8183936737128 124.5685462179979 125.6605823013066 [13] 128.7106647045467 116.3830616911681 120.6424874454661 130.1845143553381 131.3087061533496 124.2051838899657 [19] 122.4152007594897 118.5603630994196 125.6499879613931 130.4590540768779 120.0483665135104 125.7390187089629 [25] 128.3522876051743 126.9468484651417 126.3746312355432 116.9733911064371 123.4602969373899 122.7968444006815 [31] 124.9700778037051 124.2663388696768 140.6244929523544 124.3357595465487 126.0265919783500小数点以下の桁数が異様に長いが,平均値と標準偏差は指定されたものと等しい
> cat("mean =", mean(A), " S.D.=", sd(A), "\n") mean = 123.45 S.D.= 6.78 > cat("mean =", mean(B), " S.D.=", sd(B), "\n") mean = 126.28 S.D.= 5.43この 2 つの「生データ」を使って,t.test 関数で検定を行うことができる。
> t.test(A, B) Welch Two Sample t-test data: gendat1(24, 123.45, 6.78) and gendat1(35, 126.28, 5.43) t = -1.7041, df = 42.164, p-value = 0.09571 alternative hypothesis: true difference in means is not equal to 0 95 percent confidence interval: -6.180952255076989 0.520952255076993 sample estimates: mean of x mean of y 123.45 126.28専用の関数を使った結果は,以下のようになる。同じ結果が得られることが確認できる。
> my.t.test(24, 123.45, 6.78^2, 35, 126.28, 5.43^2) ウエルチの方法による,二群の平均値の差の検定 data: n1 = 24, mean1 = 123.45, variance1 = 45.9684 n2 = 35, mean2 = 126.28, variance2 = 29.4849 t = 1.7041, df = 42.164, p-value = 0.09571
この方法を見てくると,検定というのが元のデータのほんの限られた情報しか使っていないことに改めて気づく。