統計の最も基本的な知識でありながら、理解し難いものに「正規分布」がある。
「釣り鐘状の分布」などと言われているが、実際のところ、
何故、この分布が重要なのか？その理解に苦しむ人は多い。
例えば、庭の植木に生えている花の「花弁」を無作為に1,000枚採ってきたとする。
季節によっては、花など生えていないかもしれない。その時は同じ種類の「葉」でも良い。
とにかく、サイズを測りやすそうな植物の花弁か葉を1,000枚用意する。
そして、その長さを測るとしよう。まぁ、実際には幅でも良いのだが。
当たり前の話ではあるが、無作為に集めた1,000枚の花弁（葉）のサイズが
全て同じということはあり得ない。必ず、微妙な差が生じるのである。
では、この差がどのように現れるのか？
平均的なサイズが最も多く、それより大きくなるか、小さくなるほど数は少なくなる。
このような状況をヒストグラムで表現すると釣り鐘状に値が並ぶ。
このような分布をすることを「正規分布」と呼ぶ。
「自然界の至るところで見ることができる。」と、一般的には言われているが、
実際にそのように分布しているデータを探すのは意外に大変。
試験の結果も、身長のデータも、植物の花弁や葉も、実際には理論通りにはいかない。
では、不必要あるかというと、そういうわけでもない。
正規分布を前提とすることで、多くの統計的な処理が説明しやすくなり、
また、実用上の問題としては、正規分布に従っていなくても良いこともある。
もちろん、正規分布の前提が崩れると危険なものも存在するが…。
それでも、平均と標準偏差が分かれば、ある値の現れる確率を得れるという性質は重要。
それゆえに、現在においても、正規分布を前提とする手法は多い。
正規分布を盲目的に信用しそれを前提に置いてはならないという主張は確かにある。
しかし、だからと言って、知らないで良いわけでは断じてない！
さて、話は逸れたが、実際のデータを使いながら整理してみる。
Rには、様々なサンプルデータがあるので、今回はアヤメ（iris）のデータを使ってみる。
とりあえず、データを準備する。
# アヤメのデータを確認する。以下のコマンドを実行するだけ。
head(iris)
データの中身を見てみると、以下のようになっている。
> head(iris)
Sepal.Length Sepal.Width Petal.Length Petal.Width Species
1 5.1 3.5 1.4 0.2 setosa
2 4.9 3.0 1.4 0.2 setosa
3 4.7 3.2 1.3 0.2 setosa
4 4.6 3.1 1.5 0.2 setosa
5 5.0 3.6 1.4 0.2 setosa
6 5.4 3.9 1.7 0.4 setosa
Sepal はガクのことで、Petal は花弁のことである。
Speciesには、アヤメの種類が入っている。
まずは、全体の状況を可視化して確認したい。
beeswarm()関数とboxplot()関数を使って、アヤメの花弁の長さを可視化してみる。
library(beeswarm)
boxplot((iris$Petal.Length ~ iris$Species))
beeswarm((iris$Petal.Length ~ iris$Species), pch=16, col="red", add=TRUE)
三種類のアヤメ、setosa、versicolor、virginicaの花弁の長さの分布が示されている。
今回は、少々、変わった方法でデータを指定しているが、
これはモデル式と呼ばれる方法による指定。
Petal.Lengthという変数の値をSpeciesという変数で分けるという意味。
非常に有名なデータで、よく、チュートリアルデータとして用いられる。
このデータは、ロナルド・フィッシャーという有名な統計学者が、
判別分析という手法を紹介するためのサンプルデータとして紹介したものである。
元々は、エドガー・アンダーソンという人物が、
アヤメの形状を定量的に分析するために、
北アメリカのガスペ島という所で収集したものらしい。
とのかく、このデータは、3種類のアヤメについて、50株ずつ採集し、
ガクの長さと幅、花弁の長さと幅の4変数について測定されたものである。
さて、今回は、解りやすいようにsetosaの花弁のデータを用いることにする。
なんとなく、正規分布に従っていそうでもある。
X <- split(iris$Petal.Length, iris[,5])[1][[1]]
hist(X.sample)
ここでは、split()関数が新しく登場している。
この関数は、元のデータから、特定の列の値で切り分けるための関数である。
今回の場合は、花弁の長さのデータを、setosaで切り分けている。
最後に、[[1]]を付けているのは、list型からベクトルを取り出すためのオマジナイ。
さて、準備が出来たところで、早速、正規分布の特性について整理していく。
まず、正規分布は、ある値の発生する可能性を「確率」として表すことができ、
その確率を推定するための「関数」を「平均」と「標準偏差」によって決めることができる。
# setosaの花弁の長さの最小と最大をxの取る範囲とする
X.range = seq(min(X), max(X), 0.01)
# setosaの花弁の長さの平均と標準偏差で正規分布を描く
X.dnorm <- dnorm(seq(min(X), max(X), by=0.01),mean(X),sd(X))
plot(x=X.range, X.dnorm, type="l", col="red", xaxt="n", yaxt="n", ylim=c(0,2.6), xlab="", ylab="")
# setosaの元のデータの分布を重ね合わせる。
par(new=TRUE)
plot(density(X, cut=0), xlim=c(min(X), max(X)), ylim=c(0, 2.6))
Rでは、正規分布を描くためのdnorm()関数が用意されていて、
このパラメータに、x軸の値のベクトル、平均、標準偏差を与えてやると、
x軸の値に対応するy軸の値、すなわち、正規分布を描く関数によって導かれた値が出る。
今回は、x軸の値は、setosaの最小値〜最大値の範囲を等間隔に0.01刻みで与えている。
ある分位数を与えていると考えることもできる。Rのマニュアルではそのようになっている。
要するに、y軸の値を求めるときのx軸の値。刻み幅が小さいと鋭角な折れ線となる。
この図では、元のデータとの重ね合わせを行なっているが、今の段階では、
とりあえず、分布密度を取ったヒストグラムに対応するものと理解しておけば良い。
ヒストグラムとは異なり、ある関数を用いて推定を掛けている。
これを見てみると、赤い曲線で描いた正規分布の曲線と、
黒い曲線で描いた元のデータの推定は、なんとなく似た状況を示している。
さて、話を元に戻す。
そうそう、正規分布を描くと、ある値の取る可能性を知ることができる。
例えば、「花弁の長さが1.3センチよりも小さい」という状況は、
どの程度の確率で発生するのか？
そのようなことが、簡単に知ることができるのである。
イマイチ、ピンと来ないかもしれないが、
ある事象の起こり得る確率を求めるのは非常に難しい。
サイコロの目のように、離散的な値であれば良いが、
連続的な値の場合はもっと難しい。
正規分布というのは、「平均」と「標準偏差」だけ与えてやれば、
簡単に確率を教えてくれる。そういう、魔法がかかっているのである。
では、どうすれば、その「確率」とやらを教えてくれるのだろうか？
それを知るためには、まずは、正規分布の「関数」を知らなくてはならない。
ここからは、少し辛抱。厄介な式が出てくるが、要するに、「平均」と「標準偏差」、
そして、「xの値」が分かれば、釣り鐘状の山を描けるような関数であることが解れば良い。
唯一、意味不明なものが、「e」 の記号。これは、指数関数を表している。
説明をし始めると先に進まないので、とりあえずは、この記号は関数であって、
この関数のべき乗の部分が解れば、何かの値が出てくると理解しておく。
いずれ、指数関数の説明はしなければならない。
何度も繰り返すが、この関数は、平均と標準偏差だけが解れば良い。
この二つさえ知ることができれば、確率がでてくる。
xの値は分位数なので気にする必要は無い。
では、なぜ、確率が求められるかというと、
実は、この関数によって描かれる釣り鐘状の山の面積が「1」となるため。
つまり、面積の広さを確率として求めているのである。
したがって、確率を求めたいある範囲があって、それに対応するxの値が解れば、
この正規分布を描く関数を使って、面積を求めることができる。
そして、その面積こそが、真に知りたい確率となる。
ところで、正規分布において平均と標準偏差は何を意味するのか？
まず、この関数において平均は釣り鐘状の分布の最も高い地点に当たる。
そして、厳密に正規分布に従う場合この地点は中央値にも一致する。
そう言えば、対象数が多ければ、平均と中央値が近づくという話をした。
つまり、そのような状況は正規分布を前提としているのである。
一方、標準偏差というのは、要するに、バラツキの指標であって、
バラツキが大きいということは、それだけ、「山の裾」は広がる。
面積が「1」に固定された状態で、山裾が広がるということは、山の形は扁平になり、
逆に、標準偏差が低いということは、裾が狭まるので、山の形は切り立つ。
かなり、遠回りしたが、先程の例の場合の場合はどうなるか？つまり、
「花弁の長さが1.3センチよりも小さい」場合は、どのような確率であろうか？
実際には、かなり複雑な計算となるので、まずはRでの実行から。
pnorm(1.3, mean(X), sd(X))
Rでは、pnorm()関数を用いるとこの計算ができる。実行結果は以下の通り。
> pnorm(max(X), mean(X), sd(X))
[1] 0.1754524
この結果から、17.5%の確率で出現することが解った。
当然、この結果は、あくまで正規分布に従ったという前提の下ではあるが、
正規分布に従っているという前提を与えることで、容易に確率を求めることができる。
では、これはどのような状況なのか、可視化して確認してみる。
# 正規分布を描き直す
X.dnorm <- dnorm(seq(min(X), max(X), by=0.01),mean(X),sd(X))
plot(x=X.range, X.dnorm, type="l", col="red", xlab="", ylab="")
# タイトルと軸ラベルの設定
title(main="Probability of occurence for petal length within 1.3cm")
title(xlab="Length of Petal", ylab="density")
# 面積の部分を塗りつぶす
poly.x <- seq(1,1.3,0.01)
poly.y <- c(0,dnorm(poly.x, mean(X),sd(X)),0)
polygon(c(1,poly.x,1.3), poly.y, col="red")
# 理解しやすいようにオマケ
caption <- paste(round(pnorm(1.3, mean(X), sd(X))*100, 2), "%")
text(1.2, 0.25, caption, cex=1.3, col="white")
結局は、この赤色の面積のことを言っているのである。
ところで、この面積は、正規分布の積分で表され、
今回の例で言うと、1.3に当たる部分を一般化して「a」と置き、
「a」よりも小さい値が発生する確率を、
と置くと次のようになる。
積分は、丁度、f(x)とx軸で囲まれた部分の面積に相当するため、
この式によって、ある値「a」よりも下の値を取る場合の確率が求まる。
ちなみに、面積「1」の場合は、つまり全ての範囲は、以下のようになる。
実は、この積分を計算するのは、少々厄介で、しかも、どうにも解りにくい。
そこで、一般的には、標準正規分布と呼ばれる、
ある特殊な正規分布の形に変換してこの部分の面積を求めることができる。
この話については、次回に。
「釣り鐘状の分布」などと言われているが、実際のところ、
何故、この分布が重要なのか？その理解に苦しむ人は多い。
例えば、庭の植木に生えている花の「花弁」を無作為に1,000枚採ってきたとする。
季節によっては、花など生えていないかもしれない。その時は同じ種類の「葉」でも良い。
とにかく、サイズを測りやすそうな植物の花弁か葉を1,000枚用意する。
そして、その長さを測るとしよう。まぁ、実際には幅でも良いのだが。
当たり前の話ではあるが、無作為に集めた1,000枚の花弁（葉）のサイズが
全て同じということはあり得ない。必ず、微妙な差が生じるのである。
では、この差がどのように現れるのか？
平均的なサイズが最も多く、それより大きくなるか、小さくなるほど数は少なくなる。
このような状況をヒストグラムで表現すると釣り鐘状に値が並ぶ。
このような分布をすることを「正規分布」と呼ぶ。
「自然界の至るところで見ることができる。」と、一般的には言われているが、
実際にそのように分布しているデータを探すのは意外に大変。
試験の結果も、身長のデータも、植物の花弁や葉も、実際には理論通りにはいかない。
では、不必要あるかというと、そういうわけでもない。
正規分布を前提とすることで、多くの統計的な処理が説明しやすくなり、
また、実用上の問題としては、正規分布に従っていなくても良いこともある。
もちろん、正規分布の前提が崩れると危険なものも存在するが…。
それでも、平均と標準偏差が分かれば、ある値の現れる確率を得れるという性質は重要。
それゆえに、現在においても、正規分布を前提とする手法は多い。
正規分布を盲目的に信用しそれを前提に置いてはならないという主張は確かにある。
しかし、だからと言って、知らないで良いわけでは断じてない！
さて、話は逸れたが、実際のデータを使いながら整理してみる。
Rには、様々なサンプルデータがあるので、今回はアヤメ（iris）のデータを使ってみる。
とりあえず、データを準備する。
# アヤメのデータを確認する。以下のコマンドを実行するだけ。
head(iris)
データの中身を見てみると、以下のようになっている。
> head(iris)
Sepal.Length Sepal.Width Petal.Length Petal.Width Species
1 5.1 3.5 1.4 0.2 setosa
2 4.9 3.0 1.4 0.2 setosa
3 4.7 3.2 1.3 0.2 setosa
4 4.6 3.1 1.5 0.2 setosa
5 5.0 3.6 1.4 0.2 setosa
6 5.4 3.9 1.7 0.4 setosa
Sepal はガクのことで、Petal は花弁のことである。
Speciesには、アヤメの種類が入っている。
まずは、全体の状況を可視化して確認したい。
beeswarm()関数とboxplot()関数を使って、アヤメの花弁の長さを可視化してみる。
boxplot((iris$Petal.Length ~ iris$Species))
beeswarm((iris$Petal.Length ~ iris$Species), pch=16, col="red", add=TRUE)
三種類のアヤメ、setosa、versicolor、virginicaの花弁の長さの分布が示されている。
今回は、少々、変わった方法でデータを指定しているが、
これはモデル式と呼ばれる方法による指定。
Petal.Lengthという変数の値をSpeciesという変数で分けるという意味。
非常に有名なデータで、よく、チュートリアルデータとして用いられる。
このデータは、ロナルド・フィッシャーという有名な統計学者が、
判別分析という手法を紹介するためのサンプルデータとして紹介したものである。
元々は、エドガー・アンダーソンという人物が、
アヤメの形状を定量的に分析するために、
北アメリカのガスペ島という所で収集したものらしい。
とのかく、このデータは、3種類のアヤメについて、50株ずつ採集し、
ガクの長さと幅、花弁の長さと幅の4変数について測定されたものである。
さて、今回は、解りやすいようにsetosaの花弁のデータを用いることにする。
なんとなく、正規分布に従っていそうでもある。
hist(X.sample)
ここでは、split()関数が新しく登場している。
この関数は、元のデータから、特定の列の値で切り分けるための関数である。
今回の場合は、花弁の長さのデータを、setosaで切り分けている。
最後に、[[1]]を付けているのは、list型からベクトルを取り出すためのオマジナイ。
さて、準備が出来たところで、早速、正規分布の特性について整理していく。
まず、正規分布は、ある値の発生する可能性を「確率」として表すことができ、
その確率を推定するための「関数」を「平均」と「標準偏差」によって決めることができる。
X.range = seq(min(X), max(X), 0.01)
# setosaの花弁の長さの平均と標準偏差で正規分布を描く
X.dnorm <- dnorm(seq(min(X), max(X), by=0.01),mean(X),sd(X))
plot(x=X.range, X.dnorm, type="l", col="red", xaxt="n", yaxt="n", ylim=c(0,2.6), xlab="", ylab="")
# setosaの元のデータの分布を重ね合わせる。
par(new=TRUE)
plot(density(X, cut=0), xlim=c(min(X), max(X)), ylim=c(0, 2.6))
Rでは、正規分布を描くためのdnorm()関数が用意されていて、
このパラメータに、x軸の値のベクトル、平均、標準偏差を与えてやると、
x軸の値に対応するy軸の値、すなわち、正規分布を描く関数によって導かれた値が出る。
今回は、x軸の値は、setosaの最小値〜最大値の範囲を等間隔に0.01刻みで与えている。
ある分位数を与えていると考えることもできる。Rのマニュアルではそのようになっている。
要するに、y軸の値を求めるときのx軸の値。刻み幅が小さいと鋭角な折れ線となる。
この図では、元のデータとの重ね合わせを行なっているが、今の段階では、
とりあえず、分布密度を取ったヒストグラムに対応するものと理解しておけば良い。
ヒストグラムとは異なり、ある関数を用いて推定を掛けている。
これを見てみると、赤い曲線で描いた正規分布の曲線と、
黒い曲線で描いた元のデータの推定は、なんとなく似た状況を示している。
さて、話を元に戻す。
そうそう、正規分布を描くと、ある値の取る可能性を知ることができる。
例えば、「花弁の長さが1.3センチよりも小さい」という状況は、
どの程度の確率で発生するのか？
そのようなことが、簡単に知ることができるのである。
イマイチ、ピンと来ないかもしれないが、
ある事象の起こり得る確率を求めるのは非常に難しい。
サイコロの目のように、離散的な値であれば良いが、
連続的な値の場合はもっと難しい。
正規分布というのは、「平均」と「標準偏差」だけ与えてやれば、
簡単に確率を教えてくれる。そういう、魔法がかかっているのである。
では、どうすれば、その「確率」とやらを教えてくれるのだろうか？
それを知るためには、まずは、正規分布の「関数」を知らなくてはならない。
ここからは、少し辛抱。厄介な式が出てくるが、要するに、「平均」と「標準偏差」、
そして、「xの値」が分かれば、釣り鐘状の山を描けるような関数であることが解れば良い。
唯一、意味不明なものが、「e」 の記号。これは、指数関数を表している。
説明をし始めると先に進まないので、とりあえずは、この記号は関数であって、
この関数のべき乗の部分が解れば、何かの値が出てくると理解しておく。
いずれ、指数関数の説明はしなければならない。
何度も繰り返すが、この関数は、平均と標準偏差だけが解れば良い。
この二つさえ知ることができれば、確率がでてくる。
xの値は分位数なので気にする必要は無い。
では、なぜ、確率が求められるかというと、
実は、この関数によって描かれる釣り鐘状の山の面積が「1」となるため。
つまり、面積の広さを確率として求めているのである。
したがって、確率を求めたいある範囲があって、それに対応するxの値が解れば、
この正規分布を描く関数を使って、面積を求めることができる。
そして、その面積こそが、真に知りたい確率となる。
ところで、正規分布において平均と標準偏差は何を意味するのか？
まず、この関数において平均は釣り鐘状の分布の最も高い地点に当たる。
そして、厳密に正規分布に従う場合この地点は中央値にも一致する。
そう言えば、対象数が多ければ、平均と中央値が近づくという話をした。
つまり、そのような状況は正規分布を前提としているのである。
一方、標準偏差というのは、要するに、バラツキの指標であって、
バラツキが大きいということは、それだけ、「山の裾」は広がる。
面積が「1」に固定された状態で、山裾が広がるということは、山の形は扁平になり、
逆に、標準偏差が低いということは、裾が狭まるので、山の形は切り立つ。
かなり、遠回りしたが、先程の例の場合の場合はどうなるか？つまり、
「花弁の長さが1.3センチよりも小さい」場合は、どのような確率であろうか？
実際には、かなり複雑な計算となるので、まずはRでの実行から。
pnorm(1.3, mean(X), sd(X))
Rでは、pnorm()関数を用いるとこの計算ができる。実行結果は以下の通り。
> pnorm(max(X), mean(X), sd(X))
[1] 0.1754524
当然、この結果は、あくまで正規分布に従ったという前提の下ではあるが、
正規分布に従っているという前提を与えることで、容易に確率を求めることができる。
では、これはどのような状況なのか、可視化して確認してみる。
X.dnorm <- dnorm(seq(min(X), max(X), by=0.01),mean(X),sd(X))
plot(x=X.range, X.dnorm, type="l", col="red", xlab="", ylab="")
# タイトルと軸ラベルの設定
title(main="Probability of occurence for petal length within 1.3cm")
title(xlab="Length of Petal", ylab="density")
# 面積の部分を塗りつぶす
poly.x <- seq(1,1.3,0.01)
poly.y <- c(0,dnorm(poly.x, mean(X),sd(X)),0)
polygon(c(1,poly.x,1.3), poly.y, col="red")
# 理解しやすいようにオマケ
caption <- paste(round(pnorm(1.3, mean(X), sd(X))*100, 2), "%")
text(1.2, 0.25, caption, cex=1.3, col="white")
結局は、この赤色の面積のことを言っているのである。
ところで、この面積は、正規分布の積分で表され、
今回の例で言うと、1.3に当たる部分を一般化して「a」と置き、
「a」よりも小さい値が発生する確率を、
と置くと次のようになる。
積分は、丁度、f(x)とx軸で囲まれた部分の面積に相当するため、
この式によって、ある値「a」よりも下の値を取る場合の確率が求まる。
ちなみに、面積「1」の場合は、つまり全ての範囲は、以下のようになる。
実は、この積分を計算するのは、少々厄介で、しかも、どうにも解りにくい。
そこで、一般的には、標準正規分布と呼ばれる、
ある特殊な正規分布の形に変換してこの部分の面積を求めることができる。
この話については、次回に。
