廿TT

譬如水怙牛過窓櫺 頭角四蹄都過了 因甚麼尾巴過不得

二変量正規分布の密度関数、楕円

2変量正規分布の密度関数

f:id:abrahamcow:20140302202938p:plain

のグラフをかいてみる

 ##2変量正規密度
 normal2dens <- function(x,y,r=0.8){
   det <- 1-r^2
   1/(2*pi*sqrt(det))*exp(-(x^2-2*r*x*y+y^2)/(2*det))
 }
  x<- seq(-3,3,length = 100)
  y<-x
  z <- outer(x,y,r=-0.8,normal2dens)
 persp(x,y,z)


地図みたいな等高線図でかくと

 contour(x,y,z)
 abline(0,1,lty=2)
 abline(v=0)
 abline(h=0)


楕円になっている
なんで楕円になるんだろう? これをかんがえる
等高線図っていうのは,関数を同じ高さのところで輪切りにして見ているので
密度関数 = c (cは任意の定数)
と置いてみる

f:id:abrahamcow:20140303024005p:plain

と整理できる
cは任意なので右辺
 -2 \log ( 2\pi \sigma _1 \sigma _2 \sqrt{1-\rho ^2}\cdot c )
をあらためてcと置く
 \frac{1}{(1-\rho ^2)}\left[\frac{(x_1 -\mu _1)^2}{\sigma ^2_1}-\rho \frac{(x_1-\mu _1)(x_2 - \mu _2)}{\sigma _1 \sigma _2}+\frac{(x_2-\mu _2)^2}{\sigma ^2_2}\right] =c
ここで,
 Y_1 = \frac{X_1 - \mu _1}{\sigma _1},Y_2 = \frac{X_2 - \mu _2}{\sigma _2}
と変数変換すると
 \frac{1}{(1-\rho ^2)}( y^2_1-2\rho y_1 y_2+y^2_2 ) =c
これを行列,ベクトルを使って
 \bf{y}'A\bf{y}=(1-\rho ^2)c
と書ける
 \bf{y}=\left( \begin{array}{c} y_1 \\ y_2 \\ \end{array} \right)
 A=\left( \begin{array}{c} 1 & -\rho \\ -\rho & 1 \\ \end{array} \right)
Aの固有値を求める
 \left| \begin{array}{c} 1 - \lambda & -\rho \\ -\rho & 1- \lambda \\ \end{array} \right| = 0
を解いて
λ=1±ρ
対応する(単位)固有ベクトル
\frac{1}{\sqrt{2}}\left( \begin{array}{c} 1 \\ 1 \\ \end{array} \right)
\frac{1}{\sqrt{2}}\left( \begin{array}{c} -1 \\ 1 \\ \end{array} \right)
で,
 P=\frac{1}{\sqrt{2}}\left( \begin{array}{c} 1 & -1\\ 1&1 \\ \end{array} \right)
により
 P'AP= \left( \begin{array}{c} 1-\rho & 0\\ 0&1-\rho \\ \end{array} \right)
と対角化できる

 \tilde{ \bf{y}} = \left( \begin{array}{c} \tilde y_1 \\ \tilde y_2 \\ \end{array} \right) =P'\bf{y} \Leftrightarrow \bf{y}=P\tilde{ \bf{y}}
と変換して,新しい座標系
 (\tilde y_1 , \tilde y_2)
で方程式をかんがえると,
 \bf{y}'A\bf{y}=\tilde{\bf{y}}'P' A P \tilde \bf{y}=(1-\rho ^2)c
 (1-\rho)\tilde y_1+(1+\rho)\tilde y_1=(1-\rho ^2)c
 \left( \frac{\tilde y_1}{c\sqrt{1+\rho}} \right)^2 +\left( \frac{\tilde y_1}{c\sqrt{1-\rho}} \right)^2=1

これで楕円の方程式の見慣れた形になった
楕円のグラフは

 x<- seq(-3,3,length = 100)
 y<-x
ellipse <- function(x,y,r=0.8){
  x^2/( ( 1+r ) )+y^2/( ( 1-r ) )
}
 z<-outer(x,y,ellipse)
 contour(x,y,z)
 abline(v=0)
 abline(h=0)


だけれども
 P=\frac{1}{\sqrt{2}}\left( \begin{array}{c} 1 & -1\\ 1&1 \\ \end{array} \right)
が原点まわり45°の回転であった
−45°回転させて,元の座標系に戻すと,上の方の密度関数のグラフと一致するのがわかる


多変量正規分布の密度関数は楕円上にある


以上は
An Introduction to Multivariate Statistical Analysis
のp.23あたりに書いてあることを補ってみたものです

An Introduction to Multivariate Statistical Analysis

An Introduction to Multivariate Statistical Analysis

追記(2014年3月3日)

はてなブログに移行したところ、一部はてなTeX 記法が数式を表示してくれなくなってたので、
TeXclipで画像にして貼り直した。

どうも長い式は苦手みたいだ。

2変量正規分布の密度関数

\begin{align*}
\frac{1}{2\pi \sigma _1 \sigma _2 \sqrt{1-\rho ^2}}
\exp \left(
\left{-\frac{1}{2(1-\rho ^2)}\left[\frac{(x_1 -\mu _1)^2}{\sigma ^2_1}-\rho \frac{(x_1-\mu _1)(x_2 - \mu _2)}{\sigma _1 \sigma _2}+\frac{(x_2-\mu _2)^2}{\sigma ^2_2}\right]\right}
\right)
\end{align*}

「密度関数=c と置いてみる」のところ

\begin{align*}
\exp\left(
-\frac{1}{2(1-\rho ^2)}\left[\frac{(x_1 -\mu _1)^2}{\sigma ^2_1}
-\rho \frac{(x_1-\mu _1)(x_2 - \mu _2)}{\sigma _1 \sigma _2}+\frac{(x_2-\mu _2)^2}{\sigma ^2_2}\right]\right)
&=2\pi \sigma _1 \sigma _2 \sqrt{1-\rho ^2}\cdot c \\
-\frac{1}{2(1-\rho ^2)}
\left[
\frac{(x_1 -\mu _1)^2}{\sigma ^2_1}-\rho \frac{(x_1-\mu _1)(x_2 - \mu _2)}{\sigma _1 \sigma _2}+\frac{(x_2-\mu _2)^2}{\sigma ^2_2}
\right]
&= \log ( 2\pi \sigma _1 \sigma _2 \sqrt{1-\rho ^2}\cdot c )\\
\frac{1}{(1-\rho ^2)}
\left[\frac{(x_1 -\mu _1)^2}{\sigma ^2_1}-\rho \frac{(x_1-\mu _1)(x_2 - \mu _2)}{\sigma _1 \sigma _2}
+\frac{(x_2-\mu _2)^2}{\sigma ^2_2}\right] 
&=-2 \log ( 2\pi \sigma _1 \sigma _2 \sqrt{1-\rho ^2}\cdot c )
\end{align*}


R/S‐PLUSによる統計解析入門

R/S‐PLUSによる統計解析入門