【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/10

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

2024/01/20

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

2024/01/10

【統計学】標準偏差の不偏推定量

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

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

comment コメントをキャンセル

: R言語
R言語　解析結果を自動で色塗りできる関数【エクセル操作不要】
こんにちは、usagi-sanです。 R言語で統計解析を行う際、検定結果がsummary(ans)の中にリストとして格納されているため、検定結果の操作が難しいと思う方が多いと思います。また結果の表を ...

: R言語
【R言語】コクラン=マンテル=ヘンツェル検定
R言語でコクラン=マンテル=ヘンツェル検定を実行する関数やその実行例を紹介します。各層のオッズ比が1であるかの検定であるコクラン=マンテル=ヘンツェル検定をRで行っていきます。また、カイ二乗検定で ...

: R言語
【R言語】1標本t検定テンプレートスクリプト
R言語の1標本t検定のテンプレートファイルを紹介します。 t検定のスクリプトファイルは以下のgithubまたは、ダウンロードリンクより入手できます。 t検定の実行したいときのテンプレートファイルとして ...

: R言語
【R言語】一般化線形モデル　関数glmの使いかた
一般化線形モデルを行う方法を解説していきます。 R言語で回帰分析を行う方法については、以下の記事で紹介しています。一般化線形モデルを用いることで、母分布が正規分布でなくても線形モデルを構築することが ...

: R言語
【R言語】加重平均の計算　重み付きの平均　記述統計
R言語で重み付きの平均である加重平均を計算する関数やその使用例について紹介していきます。実行例では、加重平均の具体的な計算例について見ていき、どのような場面で用いられるのか理解を深めていきます。こ ...

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

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

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

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

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

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

母平均の信頼区間

母平均の差の信頼区間

関数t.test

実行例

母平均の信頼区間

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

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

まとめ