R言語

【R言語】関数t.testを用いた2標本t検定・母平均の差の検定

  1. HOME >
  2. R言語 >

【R言語】関数t.testを用いた2標本t検定・母平均の差の検定

スポンサーリンク

母平均の差についての検定を行う方法を解説していきます。

R言語を用いることで、母平均の差の検定である2標本t検定が簡単に実行できます。

1標本の時と同様に、検定結果などは自分で計算する必要がなく非常に楽です。

また統計学に詳しくなくても、検定結果の各項目について詳しく解説していきます。

今回扱うプログラミングコードは以下よりダウンロードできます。

テンプレート用として、ぜひぜひダウンロードしてください。

1標本t検定が見たいという方は次のリンクをクリックしてください。

【R言語】関数t.testを用いた母平均の検定 1標本t検定

今回はR言語で、母平均についての検定であるt検定を行う方法を解説していきます。 t検定を実行する関数やその実行例を紹介します。 実行例では、実際にデータセットを用いて母平均の検定を行うだけでなく、検定 ...

続きを見る

また、2標本t検定のより詳しい内容は次の記事を参照してください。

統計学初学者向けの教材です
¥2,970 (2022/06/01 17:41時点 | Amazon調べ)

2標本t検定

検定の概要

2標本t検定の概要を以下にまとめました。以下、母分散が等しい場合と等しくない場合の2通りの母平均の差の検定の手順です。

分散が等しいとき

\(x_{11}, x_{12}, \ldots, x_{1n_1}\)を\(N(\mu_1, \sigma^2)\)からの\(n_1\)個の無作為標本、\(x_{21}, x_{22}, \ldots, x_{2n_2}\)を\(N(\mu_2, \sigma^2)\)からの\(n_2\)個の無作為標本とする。また、標本\(x_{11}, x_{12}, \ldots, x_{1n_1}\)、\(x_{21}, x_{22}, \ldots, x_{2n_2}\)に対応する確率変数を\(X_{11}, X_{12}, \ldots, X_{1n_1}\)、\(X_{21}, X_{22}, \ldots, X_{2n_2}\)とする。このとき、次の仮説を考える。

\begin{align}\left\{\begin{array}{cc}H_0:\ \mu_1 - \mu_2 = 0\\H_1:\ \mu_1 - \mu_2 \neq 0\end{array}\right. .\end{align}

この仮説検定の検定統計量は次で与えられる。

\begin{align}\cfrac{\bar{X}_1-\bar{X}_2-(\mu_1-\mu_2)}{\sqrt{S^2 (1/n_1+ 1/n_2)}}\sim t_{n_1+n_2-2},\quad \mathrm{under}\ H_0, \end{align}

ここに、\(\bar{X}_1\)、\(\bar{X}_2\)、\(S^2\)はそれぞれ

\begin{gather}\bar{X}_1=\cfrac{1}{n_1}\sum_{i=1}^{n_1}X_{1i},\ \ \bar{X}_2=\cfrac{1}{n_2}\sum_{i=1}^{n_2}X_{2i},\ \ S^2=\cfrac{1}{n_1+n_2-2}\bigl\{(n_1-1)S_1^2+(n_2-1)S_2\bigr\}\end{gather}

であり

\begin{align}S_1^2=\cfrac{1}{n_1-1}\sum_{i=1}^{n_1}(X_{1i}-\bar{X}_1)^2,\ \ S_2^2=\cfrac{1}{n_2-1}\sum_{i=1}^{n_2}(X_{2i}-\bar{X}_2)^2.\end{align}

また、有意水準\(\alpha\)の棄却域は次で与えられる。

\begin{align}( -\infty, - t_{n_1 + n_2 -2, \alpha/ 2}) \cup (t_{n_1 + n_2 -2, \alpha/ 2}, \infty), \end{align}

ここに、\(t_{n, \alpha}\)は自由度\(n\)のt分布の上側\(\alpha\)点である。

分散が異なるとき

\(x_{11}, x_{12}, \ldots, x_{1n_1}\)を\(N(\mu_1, \sigma_1^2)\)からの\(n_1\)個の無作為標本、\(x_{21}, x_{22}, \ldots, x_{2n_2}\)を\(N(\mu_2, \sigma_2^2)\)からの\(n_2\)個の無作為標本とする。また、標本\(x_{11}, x_{12}, \ldots, x_{1n_1}\)、\(x_{21}, x_{22}, \ldots, x_{2n_2}\)に対応する確率変数を\(X_{11}, X_{12}, \ldots, X_{1n_1}\)、\(X_{21}, X_{22}, \ldots, X_{2n_2}\)とする。このとき、次の仮説を考える。

\begin{align}\left\{\begin{array}{cc}H_0:\ \mu_1 - \mu_2 = 0\\H_1:\ \mu_1 - \mu_2 \neq 0\end{array}\right. .\end{align}

この仮説検定の検定統計量は次で与えられる。

\begin{align} \cfrac{\bar{X}_1-\bar{X}_2-(\mu_1-\mu_2)}{\sqrt{S_1^2/n_1+S_2^2/n_2}}\sim t_{\nu},\quad \mathrm{under}\ H_0 \end{align}

ここに、\(\bar{X}_1\)、\(\bar{X}_2\)、\(S_1^2\)、\(S_2^2\)はそれぞれ

\begin{gather}\bar{X}_1 = \cfrac{1}{n_1}\sum_{i=1}^{n_1}X_{1i},\ \ \bar{X}_2 =\cfrac{1}{n_2}\sum_{i=1}^{n_2}X_{2i},\\ S_1^2=\cfrac{1}{n_1-1}\sum_{i=1}^{n_1}(X_{1i}-\bar{X}_1)^2,\ \ S_2^2=\cfrac{1}{n_2-1}\sum_{i=1}^{n_2}(X_{2i}-\bar{X}_2)^2\end{gather}

であり、\(\nu\)は次で与えられる。

\begin{align}\nu = \cfrac{(s_1^2/n_1+s_2^2/n_2)^2}{s_1^4/n_1^2(n_1-1)+s_2^4/n_2^2(n_2-1)}\end{align}

また、有意水準\(\alpha\)の棄却域は次で与えられる。

\begin{align}( -\infty, - t_{\nu, \alpha/ 2}) \cup (t_{\nu, \alpha/ 2}, \infty), \end{align}

関数t.test

関数t.testでt検定を実行することができます。以下、t.testとその引数です。

# S3 method for default
t.test(x, y = NULL, alternative = c("two.sided", "less", "greater"), mu = 0, paired = FALSE, var.equal = FALSE, conf.level = 0.95, …)
# S3 method for formula
t.test(formula, data, subset, na.action, …)
関数t.testの引数
xnumeric型のベクトル。1つ目のグループの標本(観測値)。
ynumeric型のベクトル。2つ目のグループの標本(観測値)。
alternativecharacter型。"two.sided"で両側検定、"greater"で右片側検定、"less"で左片側検定を指定できる。
munumeric型。帰無仮説の下での母平均、または母平均の差の値。
pairedlogical型。TRUEのとき対応のある2標本検定を行う。
var.equallogical型。TRUEのとき等分散性を仮定する。
conf.levelnumeric型。信頼区間を計算する際の信頼度数。
formulaformula型。lhs ~ rhsと書き、lhsはnumeric型の変数でrhsは2つのlevelsから成るfactor型の変数。
dataデータフレーム。検定に用いるデータセット。
subsetどの観測値を検定に用いるか指定するためのベクトル。
na.actiondataがNAを含む場合どうするか。

実行例

次のパッケージggplot2とgridExtraをインストールしましょう。ヒストグラムや箱ひげ図のプロットの際に用います。

2標本のt検定は、関数t.testを用いることで、行うことができます。データセットとして次のToothLengthを使います。

ToothLengthは歯の長さlenや薬の種類suppなどの変数から構成されています。

supp別にlenのヒストグラムや経験密度、箱ひげ図をプロットすると下の画像のようになります。

ヒストグラムと経験密度

また、supp別のlenの散布図は以下の通りです。

散布図

suppの種類によってlenの異なることが分かります。VCよりもOJの方がlenが大きいことが図から読み取れます。

実際に、t検定で統計的に2つのグループに差があるのかを確かめてみようと思います。

2群の分散が等しいとき

まずは、2つのグループの母分散が等しいと仮定した場合の検定について解説します。

「suppがVCのグループとOJのグループでlenが異なるか」という仮説検定を実行したいと思います。

この仮説検定を実行するには、関数t.testの引数を次のように設定します。

引数var.equal=TRUEにすることで、2グループの等分散性を仮定したt検定を実行することができます。

また、1行目のようにt.testの引数formulaを指定する方法と5行目のように引数xとyを指定する方法の2通りがあります。

testResult を参照すると、コンソールにt検定の結果が表示されます。

検定統計量やp値、信頼区間などの値が書いてあるのが分かります。

ポイント

今、「p-value = 0.06039」であることから、有意水準0.05の下で2つのVCとOJによってlenに違いがあるとはいえないという結論が得られました。ただし、結果の解釈に注意が必要です(有意水準次第で解釈が変わってきます)。

これらの検定結果を取得したいときは、次のようにtestResultのうしろに$をつけます。検定統計量、自由度、p値、母平均の信頼区間、標本平均、帰無仮説の下での母平均の差の値などを取得できます。

また、次のようにデータフレームに格納すると、csvファイル等に保存するときに便利です。

2群の分散が等しくないとき

次に、2群の分散が等しくないときのt検定、すなわちウェルチのt検定の実行方法について紹介します。

ウェルチのt検定を行うには、次のようにt.testの引数var.equalをFLASEにします。関数t.testの引数var.equalはデフォルトでFALSEであるため、var.equal = FALSEは書かなくても大丈夫です。

testResultを参照すると、次のように検定結果が格納されていることが分かります。

ポイント

分散が等しいときと同じく、「p-value = 0.06063」であることから、有意水準0.05の下で2つのVCとOJによってlenに違いがあるとはいえないという結論が得られました。有意水準0.1だと、「VCとOJ間でlenに差がないとはいえない」という結果が得られます。

分散が等しい場合と同様に次のようにすることで、このウェルチのt検定の結果を取得することが可能です。これらの検定結果をデータフレームに格納する例を以下に載せます。

最後に、検定結果をcsvファイルへ保存する例を紹介します。

分散が等しい場合と異なる場合の2種類のt検定の結果をcsvファイルへ出力したいと思います。

まず、2種類の検定結果が格納されたデータフレームを次のように縦に結合させて、名前を各検定名に変更します。あとは結合させたデータフレームを関数write.csvでcsvファイルへ出力すればおkです。

上記を実行すると次の画像のcsvファイルを作ることができます。

2標本t検定結果

片側仮説検定

補足として、片側検定の実行方法についても紹介します。

片側検定は次のように、t.testの引数alternativeを指定することで行うことができます。

alternative="greater"で右片側検定、alternative = "less"で左片側検定が実行できます。

上のように引数formulaで実行した場合、suppのlevelsの順番がグループの番号に対応します。つまりsuppのlevelsはOJ、VCの順なので、2行目の左片側検定の対立仮説は「OJのグループのlenの母平均<VCのグループのlenの母平均」となることに注意が必要です。

testResultを参照すると次のようになります。

ポイント

「p-value = 0.03032 < 0.05」より、suppがOJのグループのlenよりVCのグループのlenは等しいとは言えないという結論が出ました。

OJのグループよりもVCのグループの方がlenは大きいということが分かりました。

以下、検定結果をデータフレームに格納する例です。

まとめ

R言語で2標本t検定を実行する関数やその実行例をみていきました。

関数t.testを用いることで2標本の母平均の差を実行することができます。

引数var.equalによってStudentのt検定、Welchのt検定を選択することができます。

スポンサーリンク

  • この記事を書いた人
  • 最新記事

usagi-san

統計学とゲームとかをメインに解説していくよ。 数式とかプログラミングコードにミスがあったり質問があったりする場合はコメントで受け付けます。すぐに対応します。

-R言語
-, , ,