確率変数

表が出る確率が \(p \in [0, 1]\) 、裏が出る確率が \(1-p\) のコインがある。これを投げて出た目を \(x\) とする。 ただし、\(x = 1\) は表、\(x = 0\) は裏を意味するとしよう。

このとき、\(x\) という値は、どのような値を取りうるだろうか？当然、0か1であろう。 しかし、\(x\) の値は、\(p\neq 0\)、\(p\neq 1\) であれば、毎回コインを投げるたびに変わる。 すなわち、\(t\) 回目に出るコインの目を便宜的に \(x_t\) とすると、確率 \(p\) で、\(x_t = 1\)、確率 \(1-p\) で、\(x_t = 0\) となる。 つまり、コインを投げる前には \(x_t\) の値は決定せず、確率的にのみ値を知ることができる。 このような変数を、確率変数と呼び、確率変数 \(X\) は \(X = x_t\) と書くことができる。

確率分布

上記のコイントスにおいて、\(p(X = x_t) = f(x_t)\) と書くことができる。つまり、\(x_t = 1\) のとき、\(f(x_t = 1) = p\) とすると、\(p(X = 1) = p\)。同様に、\(x_t = 0\) のとき、\(f(x_t = 0) = 1 - p\) とすると、\(p(X = 0) = 1 - p\)。このとき、\(f(x_t)\) を確率分布と呼ぶ。

確率分布は、\(f(x_t) \ge 0\) かつ、\(x_t\) の取りうる場合の数すべての足し算 \(\sum_{x_t}\) に対して、\(\sum_{x_t} f(x_t) = 1\) を満たす。つまり、必ず0より大きく、場合の数すべての確率分布を足し算したときに1になるものを確率分布と呼ぶ。先程のコイントスの例では、\(f(x_t = 1) + f(x_t = 0) = p + 1-p = 1\) であり、\(p \in [0, 1]\) にて確率分布の定義を満たす。

確率分布の例

しばしば利用される確率分布の例を挙げる。以下、表記を少し変更して確率分布を \(p(x)\)、コイントスに関わる確率を \(q \in [0, 1]\) とする。表 \((x=1)\) が出る確率を \(q\) とすると、上記のコイントスにて利用した確率分布をベルヌーイ分布と呼び、

\(p(x) = q^x (1-q)^{1-x}\)

として定義される。例えば \(p(x=1) = q\)、\(p(x=0) = 1 - q\) となり、たしかに確率分布の定義を満たす。

ベルヌーイ分布に従う確率変数は、Rにおいては、

remove(list = ls()) # 変数初期化

q = 0.1 # 確率
nobs = 10 # 何回コインをふるか
A = rbinom(nobs, 1, q)
A
##  [1] 0 1 0 0 0 0 1 0 0 0

でベルヌーイ分布に従う確率変数を生成することができる (乱数のシードごとに答えは変わる)。

ベイズの定理の数値計算

\(A\), \(B\) を、\(\{0, 1\}\) の2値のいずれかをとるものとする。加えて、

\(p(A=1) = 0.01\) (i.e., \(p(A=0) = 0.99\))、

\(p(B=1 | A=1) = 0.99\) (i.e., \(p(B=0 | A=1) = 0.01\))、

\(p(B=0 | A=0) = 0.99\) (i.e., \(p(B=1 | A=0) = 0.01\))

が条件として与えられているものとする。このとき、\(p(A = 1 | B = 1)\) を計算してみる。

解析計算の詳細は、前回の講義スライドや、演習の答えを参照のこと。

remove(list = ls()) # 変数初期化

## ディレクトリ適当に合わせてください

numobs = 10000 # 観測数
A = matrix(0, 1, numobs) # 初期化
B = matrix(0, 1, numobs) # 初期化

p = 0.01 # 確率の定義
A = rbinom(numobs, 1, p) # Aの生成 (二項分布に従う確率変数 → 意味は後述)
for(i in 1:numobs){
  if(A[i] == 0){
    B[i] = rbinom(1, 1, p) # A = 0ならば、B = 1となる確率は0.01
  }
  else if(A[i] == 1){
    B[i] = rbinom(1, 1, 1-p) # A = 1ならば、B = 1となる確率は0.99
  }
}
length(which(A==1))/length(A) # A = 1となる確率
## [1] 0.0108
numseq = which(B == 1) # B = 1となっている番号
length(which(A[numseq] == 1))/length(A[numseq]) # B = 1の中で、A = 1となっているものの数
## [1] 0.5118483

確率分布の例その2

ベルヌーイ分布と同様の条件において、\(N \in \mathbb{N}\)、\(k(\le N) \in \mathbb{N}\) のとき、

\(p(k| N, q) = _NC_k q^k (1-q)^{N-k}\)

を二項分布と呼ぶ。これは、\(N\) 回コイントスを行い \(k\) 回表がでる確率を意味する。

remove(list = ls()) # 変数初期化

q = 0.1 # 確率
N = 4 # 何回コイントスを行うか(Nに相当)
tri = 10000 # 何回繰り返すか
A = rbinom(tri, N, q) # 二項分布からのサンプリングをtri回繰り返す。
count = matrix(0, N, 1) # 初期化

# 数値解の計算と理論解
for(i in 1:(N+1)){
  count[i] = length(which(A == (i-1)))
  print(choose(N, (i-1))*(q^(i-1))*((1-q)^(N-(i-1))))
}
## [1] 0.6561
## [1] 0.2916
## [1] 0.0486
## [1] 0.0036
## [1] 1e-04
# 数値解
t(count)/tri
##      [,1]   [,2]   [,3]   [,4] [,5]
## [1,] 0.65 0.2979 0.0486 0.0035    0

連続値の確率分布

確率分布の概念はコイントスなどの離散値(表か裏か、など飛び飛びの整数など)のみではなく、連続値(0.00001, 0.01, …などの実数全体)にも拡張可能である (かなり不正確な議論ではあるが、実用上は下記の理解でおおよそ問題ない)。このとき、確率分布 \(p(x)\) を、下記のように定義する。

確率変数 \(X\) が、\(X \in (x, x + \Delta x)\) に含まれる確率は、 \(p(x)\Delta x\) で与えられる( \(\Delta x \ll 1\) )。このとき、\(p(x)\) を確率変数 \(x\) の確率分布(確率密度関数)と呼ぶ。 ただし、\(p(x) \ge 0\) かつ、\(x\) の定義域 \(D\) 全体での積分 \(\int_D\) に対して \(\int_D p(x)dx = 1\) を満たすものとする。離散での確率分布から、積分になったものである。

例えば、\(x \in (0, 1)\) において一定の値をとり(i.e., この範囲にて満遍なく値を取りうる)、その他の範囲では0になる確率分布は、

curve(dunif, -1, 2) # 確率密度関数

と表現される。この確率分布の積分値が1になることは容易に確認可能である(正方形の面積)。 また、このような、とある領域においてどのような値も満遍なく実現される確率分布を、一様分布と呼ぶ。

この講義ではあまり利用しないが、知っておくべき知識としては、 \(F(x) = P(X\le x) = \int_{-\infty}^x p(x)dx\) となる \(F(x)\) を分布関数と呼ぶ。つまり、\(F(x)\)は実現する確率変数が\(x\)よりも小さくなる確率を表している。分布関数を微分したものは確率分布になる。 \(a \le x \le b\)において\(p(x) = \frac{1}{b-a}\) となる一様分布の場合であれば、 \(F(x) = 0 \ (x < a), \frac{x-a}{b-a} \ (a\le x \le b), 1 \ (x \ge b)\) である。

curve(punif, -1, 2)

乱数の生成

これまたベルヌーイ分布や二項分布と同様に、確率分布に従う乱数を生成可能 (実用的にかなり頻繁に使われる)。一様分布であれば、

par(mfrow=c(1,2)) 
hist(runif(1000,0,1))

正規分布

確率分布の中でも、最も重要な確率分布である。その確率分布は、 \(p(x | \mu, \sigma) = \frac{1}{\sqrt{2\pi\sigma^2}} \exp \left( - \frac{1}{2\sigma^2}(x - \mu)^2 \right)\) であたえられる ( \(x \in (-\infty, \infty)\) )。

正規分布はガウス分布とも呼ばれ、その確率分布は下記のような形をしている。

curve(dnorm, -5, 5) # 確率密度関数

つまり、中心周りに左右対称に分布しており、中心から離れるほどに確率分布の値は小さくなっていく。

分布関数は、

curve(pnorm, -5, 5)

であり、これは誤差関数とも呼ばれる。この講義では扱うことはないが、分野によってはかなり頻繁に顔を出す関数である。

乱数も、

par(mfrow=c(1,2)) 
hist(rnorm(1000,0,1))

として生成可能であり、これまたかなり頻繁に利用される。

期待値

ベルヌーイ分布

\(p(x| q) = q^x(1-q)^{1-x}\)

を考える。つまり、表( \(x=1\) )がでる確率が \(q\) のコイントスを考えるに等しい。このとき、\(x\) の平均を求めてみる。

\(q = 0.5\) のとき、\(x\) の平均は0.5になることは予想するに難くない。2回に1回表が出そうに漠然と思うであろう。つまり、2回コイントスをしたときの表が出る平均回数は1であろう。 そのため、1回で表が出る平均回数は0.5であることも想像できる。

一般的に、ベルヌーイ分布の平均は、\(x\) の定義域を \(D\) とすると、

\(\sum_D xp(x| q) = \sum_D xq^x(1-q)^{1-x} = 0\times q^0 \times (1-q)^{1-0} + 1\times q^1 \times (1-q)^{1-1} = q\)

となる。\(q = 0.5\) のとき、上記の記述と一致する。

さらにベルヌーイ分布の分散を計算してみよう。表が出る回数の分散はあまり直感的にはわからないであろう。 分散を計算するために、ひとまず \(x^2\) の平均を計算する。

\(\sum_D x^2p(x| q) = \sum_D x^2q^x(1-q)^{1-x} = 0\times q^0 \times (1-q)^{1-0} + 1\times q^1 \times (1-q)^{1-1} = q\)

分散は、2乗の平均から平均の2乗を差し引いたものである。したがって、

\(\mathrm{Var}(x) = \sum_D x^2p(x| q) - \left( \sum_D xp(x| q) \right)^2 = q - q^2 = q(1-q)\)

となる。ここまでの内容で重要なことを整理する。確率分布が \(p(x)\) として与えられている確率変数 \(x\) の平均は、\(x\) が離散変数のとき、

\(\sum_D x p(x)\)

で与えられる。分散を計算するときには、\(x^2\) の平均も計算した。一般的に、任意の関数 \(f(x)\) について、

\(\sum_D f(x) p(x)\)

を計算することができる。この \(\sum_D f(x) p(x)\) を、一般的に期待値と呼ぶ。

連続値の期待値

\(x\) が連続変数のとき、

\(\int_D f(x) p(x) dx\)

を期待値と呼び、特に \(\int_D x p(x) dx\) を平均と呼ぶ。

例えば、\(x \in (a, b)\) にて定義される一様分布では、\(p(x) = \frac{1}{b-a} \ (a \le x \le b)\)、それ以外は \(p(x) = 0\) である。このとき平均は、

\(\int_a^b x \frac{1}{b-a} dx = \frac{1}{b-a} \frac{1}{2} (b^2 - a^2) = \frac{1}{2}(b + a)\)

となる。例えば \(a = -1, b = 2\) のとき、平均は0.5である。以下の図より、長方形の真ん中であることがわかる。 つまり、一様分布に従う確率変数を多数生成したときに、平均は0.5になるということである。

curve(dunif(x, -1, 2), -2, 3) # 確率密度関数

mean(runif(10000, -1, 2)) # 一様分布に従う確率変数を多数生成して、平均
## [1] 0.5211511

多次元分布

複数の確率変数に対する確率分布も同様にして定義することが可能である。以下は連続値での議論だが、離散値でも同様。定義域\(D\)内における任意の確率変数\(X_1, X_2\)に対して、

\(f(X_1, X_2) \ge 0\)

であり、なおかつ

\(\int \int_D p(x_1, x_2) dx_1dx_2 = 1\)

を満たす\(p(X_1, X_2)\)を確率分布として定義する。同様にして、\(N\)次元の確率変数へと拡張可能である。

多次元分布特有の性質を挙げていく。１つ目は周辺化(全確率)

\(\int_D p(x_1, x_2) dx_2 = f(x_1)\)

である。注目したい確率変数以外はすべて足し算することで消し去るイメージである。これはベイズの定理で説明した全確率と同様のものである。この周辺化の式の両辺を\(x_1\)で積分すると、1になっていることにも注目。

この周辺化を利用することで、\(x_1\)の平均\(m_1\)

\(m_1 = \int\int_D x_1p(x_1, x_2) dx_1dx_2 = \int_Dx_1p(x_1)dx_1\)

や、期待値

\(E[f(x_1)] = \int\int_D f(x_1)p(x_1, x_2) dx_1dx_2 = \int_Df(x_1)p(x_1)dx_1\)

が計算できる。

２つ目の多次元分布特有の性質は、共分散である。共分散\(\mathrm{Cov}(x_1, x_2)\)は、各々の平均を用いて

\(\mathrm{Cov}(x_1, x_2) = \int \int_D (x_1 - m_1)(x_2 - m_2)dx_1dx_2\)

として計算できる。同様にして、相関係数も計算できる。

独立な確率変数の多次元分布

この講義では、主に完全にランダムかつ独立に観測した複数の確率変数\(X_1, X_2, ..., X_N\)の確率分布\(p(X_1, ..., X_N)\)を考えることが多い。これは、以前に学んだ独立の性質を利用して、

\(p(X_1, ..., X_N) = p(X_1)p(X_2) ... p(X_N)\)

と書くことができる。また、\(X_1\)の平均\(m_1\)は、

\(m_1 = \int\int ... \int_D x_1 p(x_1, x_2, ..., x_N)dx_1dx_2, ..., dx_N = \int_D x_1 p(x_1) dx_1\)

として計算できる。つまり、一見ややこしそうな確率分布を扱うものの、計算自体は一つの確率変数を扱うことと大差ない。

演習問題

  1. 任意の連続値の確率分布において、確率変数が \(x = 0\) となる確率はいくつか。

  2. \(x \in (a, b)\) にて定義される一様分布の分散を計算せよ。

  3. \(p(x, y) = x + y \ (0 \le x \le 1, 0 \le y \le 1)\)、\(p(x, y) = 0 \ (\mathrm{otherwise})\)

となる確率分布を考える。この確率分布において

*確率分布の定義を満たしていること、

*\(x\)と\(y\)の平均、

*\(x\)と\(y\)の標準偏差、

*\(x\)と\(y\)の共分散、

*\(x\)と\(y\)の相関係数

を計算せよ。