【R言語】分散が未知のときの母平均の信頼区間の計算

R言語で分散が未知のときの母平均の信頼区間を算出する関数やその実行例について解説します。

実行例では、母平均の信頼区間だけでなく母平均の差の信頼区間の求め方について見ていきます。2つの母集団の分散が同じときや異なるとき、さらに対応のあるデータの母平均の差の信頼区間の計算の仕方も紹介します。

この記事で扱うプログラミングコードは以下からダウンロードできます。

R言語　分散が未知のときの母平均の信頼区間

1 ファイル 1.71 KB

ダウンロード

母平均および母平均の差の信頼区間差については、次の記事を参照してください。

: 【統計学】母平均の信頼区間
標本が正規分布に従う場合の母平均の信頼区間の導出についてみていきます。母分散が既知または未知である場合についてそれぞれの信頼区間を解説します。母平均の差の信頼区間については以下を参照。信頼区間 ...
続きを見る

: 【統計学】母平均の差の信頼区間
標本が正規分布に従う場合の母平均の差の信頼区間の導出についてみていきます。母分散が既知または未知である場合について二群の母平均の差の信頼区間を解説します。母平均の信頼区間については以下の記事を参照 ...
続きを見る

分散が未知のときの母平均の信頼区間

分散が未知のときの母平均の信頼区間とR言語の母平均んの信頼区間の関数についてまとめました。

母平均と母平均の差の信頼区間の概要について列挙しています。

母平均の信頼区間

$x_1, x_2, \ldots, x_n$を$N(\mu, \sigma^2)$からの$n$個の無作為標本とする。また、標本平均と不偏標本分散を次とする。

\begin{align}\bar{x} &= \cfrac{1}{n}\sum_{i=1}^{n}x_i, \\ u^2&=\cfrac{1}{n-1}\sum_{i=1}^n(x_i-\bar{x})^2. \end{align}

母平均の信頼区間

母平均$\mu$の$100(1 - \alpha)$%信頼区間は次で与えられる。

\begin{align}\left[\bar{x} - \cfrac{u}{\sqrt{n}} t_{n - 1, \alpha / 2} , \bar{x} + \cfrac{u}{\sqrt{n}} t_{n-1, \alpha / 2} \right] .\end{align}

母平均の差の信頼区間

$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$個の無作為標本とし、それぞれの標本平均と不偏標本分散、プールした不偏標本分散を次とする。

\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},\\ u_1^2=\cfrac{1}{n_1-1}\sum_{i=1}^{n_1}(x_{1i}-\bar{x}_1)^2,\ \ u_2^2=\cfrac{1}{n_2-1}\sum_{i=1}^{n_2}(x_{2i}-\bar{x}_2)^2,\\ u = \cfrac{(n_1 - 1)u_1^2 + (n_2- 1)u_2^2}{n_1 + n_2 -2}. \end{gather}

母平均の差の信頼区間

分散が等しいとき

$\sigma_1^2 = \sigma_2^2$のとき、母平均の差$\mu_1 - \mu_2$の$100(1 - \alpha)$%信頼区間は次で与えられる。

\begin{align} \left[\bar{x}_1 - \bar{x}_2 - \sqrt{\cfrac{(n_1 + n_2) u^2}{n_1n_2}} t_{n_1 + n_2 - 2, \alpha/2}, \bar{x}_1 - \bar{x}_2 + \sqrt{\cfrac{(n_1 + n_2) u^2}{n_1n_2}} t_{n_1 + n_2 - 2, \alpha/2}\right] .\end{align}

分散が異なるとき

$\sigma_1^2 \neq \sigma_2^2$のとき、母平均の差$\mu_1 - \mu_2$の$100(1 -\alpha)$%信頼区間は次で与えられる。

\begin{align}\left[ \bar{x}_1 - \bar{x}_2 - \sqrt{\cfrac{ u_1^2}{n_1} + \cfrac{u_2^2}{n_2}} t_{\nu, \alpha/2}, \bar{x}_1 - \bar{x}_2 +\sqrt{\cfrac{ u_1^2}{n_1} + \cfrac{u_2^2}{n_2}} t_{\nu, \alpha/2}\right], \end{align}

ここに、$\nu$は次で与えられる。

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

対応のあるデータのとき

2変量正規母集団$N(\boldsymbol{\mu}, \boldsymbol{\Sigma})$からの大きさ$n$の無作為標本$(x_{11}, x_{21}), \ldots, (x_{1n}, x_{2n})$が得られたとする。ここに

\begin{gather} \boldsymbol{\mu} =\begin{pmatrix}\mu_1 \\\mu_2 \end{pmatrix}, \quad \boldsymbol{\Sigma} = \begin{pmatrix} \sigma_1^2 & \sigma_{12}\\ \sigma_{21} & \sigma_2^2 \end{pmatrix}.\end{gather}

共分散行列$\boldsymbol{\Sigma}$が未知であり対応のあるとき、母平均の差$\mu_1 - \mu_2$の$100(1 - \alpha)$%信頼区間は次で与えられる。

\begin{align} \label{eq6} \left[ \bar{x}_d - \cfrac{u_d}{\sqrt{n}} t_{n- 1, \alpha/2}, \bar{x}_d + \cfrac{u_d}{\sqrt{n}} t_{n - 1, \alpha/2}\right], \tag{6}\end{align}

ここに$n_1 = n_2 = n$、$x_{di} = x_{1i} - x_{2i},\ i = 1, \ldots, n$、

\begin{align}\bar{x}_d &= \cfrac{1}{n}\sum_{i=1}^n x_{di},\\ u_d^2 &= \cfrac{1}{n - 1}\sum_{i=1}^n (x_{di} - \bar{x}_d)^2.\end{align}

関数t.test

関数t.testで母平均および母平均の差の信頼区間を計算することができます。以下、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の引数
x	numeric型のベクトル。1つ目のグループの標本（観測値）。
y	numeric型のベクトル。2つ目のグループの標本（観測値）。
alternative	character型。"two.sided"で両側検定、"greater"で右片側検定、"less"で左片側検定を指定できる。
mu	numeric型。帰無仮説の下での母平均、または母平均の差の値。
paired	logical型。TRUEのとき対応のある2標本検定を行う。
var.equal	logical型。TRUEのとき等分散性を仮定する。
conf.level	numeric型。信頼区間を計算する際の信頼度数。
formula	formula型。lhs ~ rhsと書き、lhsはnumeric型の変数でrhsは2つのlevelsから成るfactor型の変数。
data	データフレーム。検定に用いるデータセット。
subset	どの観測値を検定に用いるか指定するためのベクトル。
na.action	dataがNAを含む場合どうするか。

実行例

母平均の信頼区間の実行例について紹介します。

グラフの作成の際に次のパッケージggplot2を用います。あらかじめインストールしておいてください。

install.packages("ggplot2")
library(ggplot2)

1 2	install.packages("ggplot2") library(ggplot2)

母平均の信頼区間

次の木に関するデータセットtreesを例に母平均の信頼区間の計算を行っていきます。

> dataset <- trees
> head(dataset)
  Girth Height Volume
1   8.3     70   10.3
2   8.6     65   10.3
3   8.8     63   10.2
4  10.5     72   16.4
5  10.7     81   18.8
6  10.8     83   19.7

> dataset <- trees

> head(dataset)

Girth Height Volume

1 8.3 70 10.3

2 8.6 65 10.3

3 8.8 63 10.2

4 10.5 72 16.4

5 10.7 81 18.8

6 10.8 83 19.7

Heightの95%信頼区間を算出するには、次のように関数t.testのオブジェクトの後ろに$conf.intを付けます。

conf_level <- 0.95
ci <- t.test(dataset$Height, conf.level = conf_level)$conf.int

1 2	conf_level <- 0.95 ci <- t.test(dataset$Height, conf.level = conf_level)$conf.int

ciを参照すると信頼区間の下限と上限、また信頼水準が格納されているのが分かります。

> ci
[1] 73.6628 78.3372
attr(,"conf.level")
[1] 0.95

> ci

[1] 73.6628 78.3372

attr(,"conf.level")

[1] 0.95

信頼区間の上限と下限を取得するには、次のようにciの1番目と2番目をインデックスで参照すれば大丈夫です。

ci_lower <- ci[1] #信頼区間の下限
ci_upper <- ci[2] #信頼区間の上限

1 2	ci_lower <- ci[1] #信頼区間の下限 ci_upper <- ci[2] #信頼区間の上限

計算した信頼区間を他のファイルに保存したいときは、次のように関数data.frameでデータフレームに計算した値を格納しwrite.csvなどで他のファイルへ保存します。

resultTable <- data.frame(sample.mean = mean(dataset$Height),
                          c.i.lower = ci[1], c.i.upper = ci[2], row.names = NULL)
write.csv(resultTable, "母平均の信頼区間_分散未知.csv", row.names = FALSE)

resultTable <- data.frame(sample.mean = mean(dataset$Height),

c.i.lower = ci[1], c.i.upper = ci[2], row.names = NULL)

write.csv(resultTable, "母平均の信頼区間_分散未知.csv", row.names = FALSE)

母平均の信頼区間

上で算出した信頼区間をグラフにすると下の画像のようになります。

ggplot(dataset, aes(x = 1, y = Height)) + 
  geom_dotplot(binaxis = "y", stackdir = "center") +
  stat_summary(aes(x = 1, y = Height), fun = mean,
               fun.min = function(x) mean(x), fun.max = function(x) mean(x),
               geom = "crossbar", width = 0.3, color = "#F8766D") +
  stat_summary(aes(x = 1, y = Height), fun.min = function(x) ci[1], fun.max = function(x) ci[2], 
               geom = 'errorbar', width = 0.25, color = "#00BFC4", size = 1) + 
  geom_text(data = data.frame(x = c(mean(dataset$Height), ci), 
                              name = factor(c("sample mean", paste0(100 * conf_level, rep("% c.i.", 2))),
                                            levels = c("sample mean", paste0(100 * conf_level, "% c.i.")))),
            aes(x = 0.7, y = x, color = name, label = paste0(c("sample mean", "c.i.lower", "c.i.upper"), ": ", round(x, 3)))) +
  labs(color = "")

ggplot(dataset, aes(x = 1, y = Height)) +

geom_dotplot(binaxis = "y", stackdir = "center") +

stat_summary(aes(x = 1, y = Height), fun = mean,

fun.min = function(x) mean(x), fun.max = function(x) mean(x),

geom = "crossbar", width = 0.3, color = "#F8766D") +

stat_summary(aes(x = 1, y = Height), fun.min = function(x) ci[1], fun.max = function(x) ci[2],

geom = 'errorbar', width = 0.25, color = "#00BFC4", size = 1) +

geom_text(data = data.frame(x = c(mean(dataset$Height), ci),

name = factor(c("sample mean", paste0(100 * conf_level, rep("% c.i.", 2))),

levels = c("sample mean", paste0(100 * conf_level, "% c.i.")))),

aes(x = 0.7, y = x, color = name, label = paste0(c("sample mean", "c.i.lower", "c.i.upper"), ": ", round(x, 3)))) +

labs(color = "")

母平均の信頼区間（分散未知）

赤線はHeightの標本平均
青線はHeightの母平均の上限と下限

を表しています。

母平均はだいたい73から78の範囲にあることを推定できました。

母平均は母集団のパラメータであり、上の図の赤線の標本平均ではないことに注意してください。

母平均の差の信頼区間（分散が等しいとき、異なるとき）

続いて母平均の差の信頼区間の算出例について見ていきます。

次のデータセットToothGrowthのsuppがVCであるグループとOJであるグループのLenの母平均の差の信頼区間を求めたいと思います。

> dataset <- ToothGrowth
> head(dataset)
   len supp dose
1  4.2   VC  0.5
2 11.5   VC  0.5
3  7.3   VC  0.5
4  5.8   VC  0.5
5  6.4   VC  0.5
6 10.0   VC  0.5

> dataset <- ToothGrowth

> head(dataset)

len supp dose

1 4.2 VC 0.5

2 11.5 VC 0.5

3 7.3 VC 0.5

4 5.8 VC 0.5

5 6.4 VC 0.5

6 10.0 VC 0.5

母平均の差の信頼区間を計算するには、関数t.testの引数を次のように与えます。

ポイント

引数var.equalがTRUEのとき母分散が等しいとき、FALSEのとき母分散が異なるときの母平均の差の信頼区間が計算されます。

conf_level <- 0.95

ci_equalVar <- t.test(len ~ supp, data = dataset, conf.level = conf_level, var.equal = TRUE)$conf.int
ci_notEqualVar <- t.test(len ~ supp, data = dataset, conf.level = conf_level, var.equal = FALSE)$conf.int

conf_level <- 0.95

ci_equalVar <- t.test(len ~ supp, data = dataset, conf.level = conf_level, var.equal = TRUE)$conf.int

ci_notEqualVar <- t.test(len ~ supp, data = dataset, conf.level = conf_level, var.equal = FALSE)$conf.int

また、引数formulaとdataを用いた書き方ではなく、次のように2つのグループのベクトルを与えても同様の結果が得られます。

len_VC <- dataset$len[dataset$supp == "VC"]
len_OJ <- dataset$len[dataset$supp == "OJ"]
t.test(len_VC, len_OJ, conf.level = conf_level, var.equal = TRUE) #formulaを用いない場合

len_VC <- dataset$len[dataset$supp == "VC"]

len_OJ <- dataset$len[dataset$supp == "OJ"]

t.test(len_VC, len_OJ, conf.level = conf_level, var.equal = TRUE) #formulaを用いない場合

母分散が等しいときと異なるときのそれぞれのsuppがVCであるグループとOJであるグループの母平均の差の信頼区間は次の画像の通りです。

resultTable <- data.frame(c.i.lower = c(ci_equalVar[1], ci_notEqualVar[1]), 
                          c.i.upper = c(ci_equalVar[2], ci_notEqualVar[2]),
                          row.names = c("equal variances", "not equal variances"))

write.csv(resultTable, "母平均の差の信頼区間_分散未知.csv")

resultTable <- data.frame(c.i.lower = c(ci_equalVar[1], ci_notEqualVar[1]),

c.i.upper = c(ci_equalVar[2], ci_notEqualVar[2]),

row.names = c("equal variances", "not equal variances"))

write.csv(resultTable, "母平均の差の信頼区間_分散未知.csv")

母平均の差の信頼区間

母平均の差の信頼区間の概要で解説したように、母分散が等しいときと異なるときではt分布の自由度が異なるため、信頼区間の上限と下限が若干異なります。

2つのグループの信頼区間をプロットする際は、次のようにggplotを実行します。

ci_OJ <- t.test(len_OJ, conf.level = conf_level)$conf.int 
ci_VC <- t.test(len_VC, conf.level = conf_level)$conf.int

dataset_ci <- data.frame(mean_len = tapply(dataset$len, dataset$supp, mean),
                         supp = as.factor(c("OJ", "VC")),
                         ci = c(dist(ci_OJ), dist(ci_VC)) / 2)

ggplot() + 
  geom_dotplot(data = dataset, aes(x = supp, y = len, color = supp, fill = supp), 
               binaxis = "y", stackdir = "center") +
  geom_line(data = dataset_ci, aes(x = supp, y = mean_len, group = 1), size = 1, alpha = 0.5) +
  geom_point(data = dataset_ci, aes(x = supp, y = mean_len, group = 1), size = 2, alpha = 0.5) +
  geom_errorbar(data = dataset_ci, aes(x = supp, y = mean_len, ymin = mean_len - ci, ymax = mean_len + ci), width = 0.25, size = 1, alpha = 0.5) +
  geom_text(data = dataset_ci, aes(x = supp, y = mean_len - ci + 1, label = paste0(conf_level, "% c.i.lower: ", round(mean_len - ci, 3)))) +
  geom_text(data = dataset_ci, aes(x = supp, y = mean_len + ci + 1, label = paste0(conf_level, "% c.i.upper: ", round(mean_len + ci, 3)))) +
  geom_text(data = dataset_ci, aes(x = supp, y = mean_len + 1, label = paste0("Sample mean: ", round(mean_len, 3))))

ci_OJ <- t.test(len_OJ, conf.level = conf_level)$conf.int

ci_VC <- t.test(len_VC, conf.level = conf_level)$conf.int

dataset_ci <- data.frame(mean_len = tapply(dataset$len, dataset$supp, mean),

supp = as.factor(c("OJ", "VC")),

ci = c(dist(ci_OJ), dist(ci_VC)) / 2)

ggplot() +

geom_dotplot(data = dataset, aes(x = supp, y = len, color = supp, fill = supp),

binaxis = "y", stackdir = "center") +

geom_line(data = dataset_ci, aes(x = supp, y = mean_len, group = 1), size = 1, alpha = 0.5) +

geom_point(data = dataset_ci, aes(x = supp, y = mean_len, group = 1), size = 2, alpha = 0.5) +

geom_errorbar(data = dataset_ci, aes(x = supp, y = mean_len, ymin = mean_len - ci, ymax = mean_len + ci), width = 0.25, size = 1, alpha = 0.5) +

geom_text(data = dataset_ci, aes(x = supp, y = mean_len - ci + 1, label = paste0(conf_level, "% c.i.lower: ", round(mean_len - ci, 3)))) +

geom_text(data = dataset_ci, aes(x = supp, y = mean_len + ci + 1, label = paste0(conf_level, "% c.i.upper: ", round(mean_len + ci, 3)))) +

geom_text(data = dataset_ci, aes(x = supp, y = mean_len + 1, label = paste0("Sample mean: ", round(mean_len, 3))))

母平均の差の信頼区間

上の画像から分かるように、赤のOJのグループの上限と青のVCのグループの下限の差が母平均の差の信頼区間の上限、OJのグループの下限と青のVCのグループの上限の差が母平均の差の信頼区間の下限に近い値をとっているのが分かります。

各グループの母平均の信頼区間を描画したものであり母平均の差の信頼区間のプロットでないことに注意してください。

母平均の差の信頼区間（対応のあるデータ）

最後に対応のあるデータの母平均の差の信頼区間の計算例を紹介します。

次のデータのsleepについて母平均の差の信頼区間を計算していきます。

> dataset <- sleep
> head(dataset)
  extra group ID
1   0.7     1  1
2  -1.6     1  2
3  -0.2     1  3
4  -1.2     1  4
5  -0.1     1  5
6   3.4     1  6

> dataset <- sleep

> head(dataset)

extra group ID

1 0.7 1 1

2 -1.6 1 2

3 -0.2 1 3

4 -1.2 1 4

5 -0.1 1 5

6 3.4 1 6

対応のあるデータの母平均の差の信頼区間を計算するには、関数t.testの引数を次のように与えます。

ポイント

対応のあるデータの母平均の差の信頼区間を計算したい場合は、次のように引数pairedをTRUEにします。

conf_level <- 0.95
ci_paired <- t.test(extra ~ group, data = dataset, conf.level = conf_level, paried = TRUE)$conf.int

1 2	conf_level <- 0.95 ci_paired <- t.test(extra ~ group, data = dataset, conf.level = conf_level, paried = TRUE)$conf.int

次のように、母平均の信頼区間と同様に母平均の差の信頼区間の上限と下限の取得が可能です。

ci_lower <- ci_paired[1] #信頼区間の下限
ci_upper <- ci_paired[2] #信頼区間の上限

1 2	ci_lower <- ci_paired[1] #信頼区間の下限 ci_upper <- ci_paired[2] #信頼区間の上限

まとめ

R言語で分散が未知のときの母平均の信頼区間を計算する関数やその実行例をみていきました。

母平均の検定のための関数t.testを用いることで分散が未知のときの母平均の信頼区間を計算することができます。

引数var.equalによって分散が等しい、または異なるときの母平均の差の信頼区間、pairedで対応のあるデータの母平均の差の信頼区間が算出できます。

usagi-san

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

2024/04/28

【R言語】F検定テンプレートスクリプト

2024/04/28

【R言語】経験密度関数・経験分布関数のプロット　関数densityとecdfの使い方

2024/04/10

【R言語】ベータ関数とガンマ関数　関数beta, gammaの使い方

usagi-sanの記事をもっと見る

-R言語
-R言語, プログラミング, 統計学, 統計解析

comment コメントをキャンセル

: R言語
【R言語】openxlsxを使ったエクセルファイルの入力・出力
R言語でエクセルファイルを操作する方法を解説します。ここではパッケージopenxlsxを使ってエクセルファイルの入力や出力を行っていきます。エクセルファイル中のデータセットの読み込みから、フォント ...

: R言語
【R言語】Rのインストール方法　Rの始め方
統計解析に特化したプログラミング言語であるR言語のインストールの手順を解説します。画像付きで分かりやすくインストールの手順をまとめました。この記事ではWindowsとmacOS向けにインストール方 ...

: R言語
【R言語】主成分分析　変数の寄与率の計算や次元の削減
R言語で主成分分析を実行する関数やその実行例を解説していきます。多変量解析において各変数の寄与率を調べるためや多変量データの要約に用いられる主成分分析をR言語で実行する流れを詳しくみていきます。実 ...

: R言語
【R言語】アンサリ・ブラッドリー検定　ansari.testの使い方
R言語でアンサリ・ブラッドリー検定を実行する方法を解説していきます。この記事では、非正規母集団からの標本に対する分散の比の検定であるアンサリ・ブラッドリー検定のための関数とその実行例を紹介します。 ...

: R言語
【R言語】csvファイルの読み込み・書き込み
R言語でcsvファイルを操作する方法を解説します。 csvファイルの読み込みや書き込みに用いる関数とその実行例について見ていきます。単純な読み込みと書き込みだけでなく、csvファイルに追記する例も紹 ...

【計算フォーム】ウィルコクソンの順位和検定　マン・ホイットニーのU検定　p値の算出

【R言語】2標本t検定テンプレートスクリプト

【R言語】分散が未知のときの母平均の信頼区間の計算

【統計学】母平均の信頼区間

【統計学】母平均の差の信頼区間

分散が未知のときの母平均の信頼区間

母平均の信頼区間

母平均の差の信頼区間

関数t.test

実行例

母平均の信頼区間

母平均の差の信頼区間（分散が等しいとき、異なるとき）

母平均の差の信頼区間（対応のあるデータ）

まとめ