見出し画像

6-4 和と差の確率変数の性質 ～ ２変量正規分布を３次元グラフにして動かす

ネイピア DS

今回の統計トピック

２変量正規分布に挑戦します！
グラフは３次元・立体です！
EXCELとPython の立体表現をご堪能ください。

赤青の3Dメガネのイラスト：「いらすとや」さんより

公式問題集の準備

「公式問題集」の問題を利用します。お手元に公式問題集をご用意ください。
公式問題集が無い場合もご安心ください！
「知る」「実践する」の章で、のんびり統計をお楽しみください！


問題を解く

📘公式問題集のカテゴリ

標本分布の分野
問4 和と差の確率変数の性質（データなし）

試験実施年月
統計検定2級 2017年6月 問9（回答番号20）

問題

公式問題集をご参照ください。

解き方

題意
正規分布に従う独立の確率変数$${X,\ Y}$$の和と差$${X+Y, \ X-Y}$$について、同時分布である「2変量正規分布」と周辺分布に関する知識を確認する問題です。

・確率変数$${X}$$と$${Y}$$は互いに独立です。
・確率変数$${X}$$は正規分布$${N(0, \sigma^2_1),\ \sigma^2_1 >0}$$に従います。
・確率変数$${Y}$$は正規分布$${N(0, \sigma^2_2),\ \sigma^2_2 >0}$$に従います。
・$${U=X+Y,\ V=X-Y}$$は確率変数です。
・確率変数$${U,\ V}$$の期待値・分散・相関係数を式展開して解答します。

公式問題集の記述を改変

２変量正規分布
確率変数$${X}$$は正規分布$${N(0, \sigma^2_1)}$$に従い、確率変数$${Y}$$は正規分布$${N(0, \sigma^2_2)}$$に従います。
この正規分布のパラメータより、確率変数$${X}$$の期待値$${E[X]=0}$$、分散$${V[X]=\sigma^2_1}$$、確率変数$${Y}$$の期待値$${E[Y]=0}$$、分散$${V[Y]=\sigma^2_2}$$です。

$${U=X+Y}$$、$${V=X-Y}$$も確率変数であり、正規分布に従う確率変数$${U,\ V}$$の同時分布は２変量正規分布になります。
$${U,\ V}$$の２変量正規分布の相関係数が０のとき、$${U,\ V}$$は独立です。

２変量正規分布のグラフは「知る」の章で紹介いたします！

結婚式のイラスト「ウェディングベル」：「いらすとや」さんより

期待値・分散・共分散の公式を用いて計算
公式に当てはめてさくさく計算しましょう。
計算式が続きます。何卒何卒～。

いろいろなお辞儀をする人と「よろしくお願いします」のイラスト：「いらすとや」さんより

①$${\boldsymbol{U,\ V}}$$の平均（期待値）の計算（問題Ⅰ・Ⅲ関係）

■確率変数$${X, Y}$$の線形結合$${aX+bY+c}$$の公式
期待値$${E[aX+bY+c]=aE[X]+bE[Y]+c}$$

公式に$${U=X+Y}$$、$${V=X-Y}$$、$${E[X]=E[Y]=0}$$を当てはめます。

$$
\begin{align*}
E[U]&=E[X+Y]\\
&=E[1X+1Y] \\
&=1E[X]+1E[Y]\\
&=E[X]+E[Y] \\
&=0+0\\
&=0\\
  \\
E[V]&=E[X-Y]\\
&=E[1X+(-1)Y] \\
&=1E[X]+(-1)E[Y]\\
&=E[X]-E[Y] \\
&=0+0\\
&=0\\
\end{align*}
$$

【問題Ⅰの答え】
確率変数$${U,\ V}$$の平均は両方とも0であり、等しいです。

紙に何かを書く人のイラスト（会社員）：「いらすとや」さんより

②$${\boldsymbol{U,\ V}}$$の分散の計算（問題Ⅱ・Ⅲ関係）

■確率変数$${X, Y}$$が独立のときの公式
$${V[aX+bY+c]=a^2V[X]+b^2V[X]}$$

公式に$${U=X+Y}$$、$${V=X-Y}$$、$${V[X]=\sigma^2_1}$$、$${V[Y]=\sigma^2_2}$$を当てはめます。

$$
\begin{align*}
V[U]&=V[X+Y]\\
&=V[1X+1Y] \\
&=1^2V[X]+1^2V[Y]\\
&=V[X]+V[Y]\\
&=\sigma^2_1+\sigma^2_2\\
  \\
V[V]&=V[X-Y]\\
&=V[1X+(-1)Y] \\
&=1^2V[X]+(-1)^2V[Y]\\
&=V[X]+V[Y]\\
&=\sigma^2_1+\sigma^2_2\\
\end{align*}
$$

確率変数$${U,\ V}$$の分散$${V[U],\ V[V]}$$は両方とも$${\sigma^2_1+\sigma^2_2}$$であり、等しいです。

確率変数$${U,\ V}$$の平均（期待値）は$${0}$$、分散は$${\sigma^2_1+\sigma^2_2}$$です。
また、$${U,\ V}$$の周辺分布も正規分布に従います。

つまり、確率変数$${U}$$は正規分布$${N(0, \sigma^2_1+\sigma^2_2)}$$に従い、確率変数$${U}$$も同様に正規分布$${N(0,\ \sigma^2_1+\sigma^2_2)}$$に従います。

【問題Ⅲの答え】
$${ \sigma^2_1+\sigma^2_2}$$の値にかかわらず、$${U,\ V}$$は同じ分布に従うことが分かりました。

紙に何かを書く人のイラスト（会社員）：「いらすとや」さんより

③$${\boldsymbol{U,\ V}}$$の共分散と相関係数の計算（問題Ⅱ関係）

確率変数$${X, Y, W, Z}$$の変数変換$${U=aX+bW,\ V=cY+dZ}$$
■共分散の分配法則的な公式
$${\mathrm{Cov}[aX+bW,\ cY+dZ]}$$
$${=ac\mathrm{Cov}[X,Y]+ad\mathrm{Cov}[X,Z]+bc\mathrm{Cov}[W,Y]+bd\mathrm{Cov}[W,Z]}$$

公式に$${U=X+Y}$$、$${V=X-Y}$$、$${V[X]=\sigma^2_1}$$、$${V[Y]=\sigma^2_2}$$を当てはめます。
また、１つの確率変数の共分散は分散なので、$${\mathrm{Cov}[X,X]=V[X],\ \mathrm{Cov}[Y,Y]=V[Y]}$$です。

$$
\begin{align*}
\mathrm{Cov}[U,V]&=\mathrm{Cov}[X+Y,X-Y]\\
&=\mathrm{Cov}[1X+1Y, 1X+(-1)Y]\\
&=\mathrm{Cov}[X,X]-\mathrm{Cov}[X,Y]+\mathrm{Cov}[Y,X]-\mathrm{Cov}[Y,Y]\\
&=V[X]-\mathrm{Cov}[X,Y]+\mathrm{Cov}[X,Y]-V[Y]\\
&=V[X]-V[Y]\\
&=\sigma^2_1-\sigma^2_2
\end{align*}
$$

確率変数$${U,\ V}$$の共分散$${\mathrm{Cov}[U,V]}$$は$${\sigma^2_1-\sigma^2_2}$$です。

ちょっと休憩です。

チャイのイラスト：「いらすとや」さんより

続いて相関係数$${r[U,V]}$$を求めましょう。
公式に分散$${V[U]=V[V]=\sigma^2_1+\sigma^2_2}$$、共分散$${\mathrm{Cov}[U,V]=\sigma^2_1-\sigma^2_2}$$を当てはめます

■相関係数の公式
相関係数$${r[X,Y]=\cfrac{\mathrm{Cov}[X,Y]}{\sqrt{V[X]V[Y]}}}$$

$$
\begin{align*}
r[U,V]&=\cfrac{\mathrm{Cov}[U,V]}{\sqrt{V[U]V[V]}}\\
  \\
&=\cfrac{\sigma^2_1-\sigma^2_2}{\sqrt{(\sigma^2_1+\sigma^2_2)\times(\sigma^2_1+\sigma^2_2)}}\\
  \\
&=\cfrac{\sigma^2_1-\sigma^2_2}{\sqrt{(\sigma^2_1+\sigma^2_2)^2}}\\
  \\
&=\cfrac{\sigma^2_1-\sigma^2_2}{\sigma^2_1+\sigma^2_2}\\
\end{align*}
$$

相関係数が０のとき、２変量正規分布の確率変数$${U,\ V}$$は独立です。
そして、相関係数が０になるのは、次式の分子$${\sigma^2_1-\sigma^2_2=0}$$のときです。

$$
\begin{align*}
r[U,V]&=\cfrac{\sigma^2_1-\sigma^2_2}{\sigma^2_1+\sigma^2_2} \\
 \\
0&=\cfrac{\sigma^2_1-\sigma^2_2}{\sigma^2_1+\sigma^2_2} \\
 \\
0&=\cfrac{0}{\sigma^2_1+\sigma^2_2}
\end{align*}
$$

これは$${\sigma^2_1=\sigma^2_2}$$のときでもあります。

【問題Ⅱの答え】
$${\sigma^2_1=\sigma^2_2}$$のときのみ、確率変数$${U,\ V}$$は独立です。

目標を達成した人のイラスト：「いらすとや」さんより

解答

⑤ Ⅰ・Ⅱ・Ⅲは正しい です。

難易度　ややむずかしい

・知識：２変量正規分布、変数変換
・計算力：数式組み立て（高）、数式計算（中）
・時間目安：2分

知る

おしながき

公式問題集の問題に接近してみましょう！
今回は、２変量正規分布のイメージを確認しましょう。
難しい数式は抜きにして。

２変量正規分布

📕公式テキスト：2.9.3 ２変量正規分布（106ページ～）
２つの確率変数のそれぞれ値とその確率の対応関係を同時確率分布と呼びます。
確率変数$${X=2}$$で$${Y=1}$$のときの確率が◯◯%、みたいな感じです。
２変量正規分布は正規分布に従う２つの確率変数の同時確率分布です。

確率変数が１つの正規分布＝平面（２次元）
確率変数が１つ＝変数が１つ＝変量が１つ＝１変量です。
平均$${0}$$、分散$${1}$$の正規分布に従う確率変数$${X}$$と$${Y}$$の確率密度関数のグラフをそれぞれ別個に見てみましょう。

よく見るベル型の曲線がそれぞれプロットされます。
横軸が確率変数$${X,\ Y}$$の値、縦軸が確率密度関数の値です。
確率変数の値が１次元です。
シンプルです。

２変量の正規分布＝立体（３次元）

ようこそ、２変量正規分布の３次元の世界へ！
正規分布に従う確率変数$${X}$$と$${Y}$$の同時確率密度関数（２変量正規分布）のグラフを見てみましょう。

グラフが３次元になりました！
三角コーンみたいです。
お好みの色彩でプロットすると３次元のグラフに親近感を持てるかもです！

光るカラーコーンのイラスト：「いらすとや」さんより

１変量の正規分布では確率変数が横軸の１次元でした。
２変量の正規分布では底面の２次元が確率変数に相当します。
左側が確率変数$${X}$$の値、右側が確率変数$${Y}$$の値です。
そして高さが$${X}$$と$${Y}$$の同時確率密度関数$${f(x, y)}$$なのです。

そしてそして、立ち上がる両壁にベル型の曲線がまるで影のように描画されています。
これらを周辺分布と呼びます。

左側の壁には周辺確率変数$${Y}$$の周辺確率密度関数$${f_y(y)}$$、右側の壁には周辺確率変数$${X}$$の周辺確率密度関数$${f_x(x)}$$をプロットしています。

この２変量正規分布を平面に表すと次のようなグラフになります。
横軸が確率変数$${X}$$、縦軸が確率変数$${Y}$$、「等高線」で同時確率密度関数$${f(x,y)}$$を表しています。
三角コーンを真上から見下ろすような見え方です。

２変量の正規分布のパラメータ

１変量の正規分布のパラメータは平均と分散です。
２変量正規分布のパラメータも平均と分散なのですが、確率変数が２つあるので形状が変わります。

平均

確率変数$${X}$$と$${Y}$$それぞれの平均を設定します。
$${\mu= \left( \begin{array}{c} \mu_x \\ \mu_y \end{array} \right) }$$です。
ベクトルです。

分散

確率変数が複数あるので、確率変数どうしの関係性が生じます。
共分散が出現するのです。

$${\Sigma = \left( \begin{matrix} \sigma^2_x & \sigma_{xy}\\ \sigma_{xy} & \sigma^2_y \end{matrix} \right) }$$です。

$${\sigma^2_x, \sigma^2_y}$$は確率変数$${X, Y}$$のそれぞれの分散、$${\sigma_{xy}}$$は確率変数$${X, Y}$$の間の共分散です。
行列の形状になります。
この行列を分散共分散行列、または、共分散行列と呼びます。

さまざまな平均・共分散行列の２変量正規分布

平均$${\mu= \left( \begin{array}{c} 0 \\ 0 \end{array} \right) }$$、共分散行列$${\Sigma = \left( \begin{matrix} 30 & 0 \\ 0 & 30 \end{matrix} \right) }$$

平均$${\mu= \left( \begin{array}{c} 0 \\ 0 \end{array} \right) }$$、共分散行列$${\Sigma = \left( \begin{matrix} 10 & 0 \\ 10 & 30 \end{matrix} \right) }$$

平均$${\mu= \left( \begin{array}{c} 0 \\ 0 \end{array} \right) }$$、共分散行列$${\Sigma = \left( \begin{matrix} 1 & 0 \\ 0 & 30 \end{matrix} \right) }$$

共分散行列の値が変わることによって、２変量正規分布の形状が平たくなったり、斜めにズレたりしました。

３次元グラフの回転

３次元グラフの視点を変えるにはグラフを回転させます。
PythonとJupyter Notebookがあればマウスでグリグリするだけで、簡単にグラフを回転させてることができます。
Pythonのコードは「Pythonで作成してみよう！」の章で紹介いたします。

アクセルジャンプをしているスケート選手のイラスト：「いらすとや」さんより


２変量正規分布の同時確率密度関数の公式

$$
f(x,y)=\cfrac{1}{2\pi \sqrt{\sigma^2_x \sigma^2_y (1- \rho^2)}} \exp \left[ -\cfrac{1}{2(1-\rho^2)} \left( \cfrac{(x-\mu_x)^2}{\sigma_x}-2\rho \left( \cfrac{x-\mu_x}{\sigma_x} \right) \left(\cfrac{y-\mu_x}{\sigma_y} \right)+\left( \cfrac{y-\mu_y}{\sigma_y} \right)^2  \right) \right]
$$

２変量正規分布の同時確率密度関数のパラメータは、確率変数$${X}$$の①平均$${\mu_x}$$・②分散$${\sigma^2_x}$$、確率変数$${Y}$$の③平均$${\mu_y}$$・④分散$${\sigma^2_y}$$、$${X,Y}$$の⑤相関係数$${\rho}$$の５つです。

公式問題集の問題では、２変量正規分布の相関係数$${\rho}$$が０のときに２つの確率変数$${X,\ Y}$$が独立である性質を利用しています。

また、確率変数$${X,Y}$$の共分散が$${\sigma_{xy}}$$のとき、相関係数$${\rho=\cfrac{\sigma_{xy}}{\sigma_x \sigma_y}}$$です。

実践する

2変量正規分布のグラフを描いてみよう

「知る」の２変量正規分布のグラフを描いてみましょう。
EXCELとPythonを活用して３次元グラフをプロットしましょう。


電卓・手作業で作成してみよう！

「知る」の章の２変量正規分布のグラフを見ていただき、可能であれば手作業で・・・、３次元グラフを描くのは厳しいですよね。。。
EXCEL、Pythonのグラフを見てイメージを膨らましていきましょう！


EXCELで作成してみよう！

データ数が多い場合、やはり手作業では非効率になります。
パソコンを利用して、手早く作表できるようになれば、実務活用がしやすくなるでしょう。

3-D等高線グラフ

３次元用のグラフをプロットしましょう。

EXCELには、２変量正規分布の同時確率密度関数を計算する機能がありません。
そこで、次のサイトの情報を利用させていただき、グラフを作ってみましょう。
ご協力ありがとうございます！

まず、サイトの方法で「２変量正規分布の同時確率密度関数の表」を作り、「3-D等高線グラフをプロット」します。

表の作成

「２変量正規分布の同時確率密度関数の公式」をEXCELの計算式で書きます。
パラメータの値を変更して、さまざまな形状のグラフを確認してみてくださいね。

サンプルファイルの計算式をご活用ください。
式に誤りがないかチェックをしていますが、万が一、ミスがありましたら、すみません。

グラフの作成

グラフに描画する表の範囲を選択して、
メニュー＞「挿入」＞「グラフ」の右隅の矢印＞「グラフの挿入」画面で「3-D等高線」を選択します。

グラフ完成イメージ

３次元グラフの回転

「グラフエリアの書式設定」で「3-D回転」の調整ができます。
この例では、「X方向に回転」を上のグラフの「20°」から「40°」に変えました。
グリグリ回して、さまざまな角度から形状を確認してみてください！

EXCELサンプルファイルのダウンロード
こちらのリンクからEXCELサンプルファイルをダウンロードできます。


Pythonで作成してみよう！

プログラムコードを読んで、データを流したりデータを変えてみたりして、データを追いかけることで、作表ロジックを把握する方法も効果的でしょう。
サンプルコードを揃えておけば、類似する作表作業を自動化して素早く結果を得ることができます。

今回は、「知る」の章で用いた２変量正規分布のグラフをプロットするPythonコードです。

コード作成にあたっては、次のサイトの情報を参考にさせていただきました。
ありがとうございます！

①インポート

  • scipy.statsのmultivariate_normalで２変量正規分布の同時確率密度関数の値を取得します。

  • Axes3Dで３次元グラフを描画します。

import numpy as np
from scipy.stats import multivariate_normal
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

②matplotlibのモードの選択
上のインタラクティブモードの設定を利用すると3Dグラフをマウスで回転できます。
なお、Google Colabなどでインタラクティブモードを利用できない場合には、下の非インタラクティブモードをご利用ください。

# インタラクティブモード:マウスで3Dグラフを回転する場合
%matplotlib notebook

# 非インタラクティブモード：Google Colabなど、上記のインタラクティブを使えない場合
# %matplotlib inline

③グラフ描画関数の定義
プロットのコードを繰り返し利用するので関数化しました。

■３次元グラフ

  • 入力パラメータは「平均ベクトルmean_vec」と「共分散行列cov_mtx」です。

  • multivariate_normal.pdfで２変量正規分布の同時確率密度関数の値を取得します。

  • plot_surfaceで３次元グラフをプロットします。

  • ２つのcontourで周辺分布の周辺確率密度関数の値（壁のベル型の曲線）をプロットします。

# 2変量正規分布の3次元グラフを描画する関数
def bi_norm_plot(mean_vec, cov_mtx):
    # x,y軸の設定
    side = np.arange(-20, 20, 0.5)
    x, y = np.meshgrid(side, side)
    xy = np.dstack((x, y))

    # z軸：2変量正規分布の確率密度関数の算出
    z = multivariate_normal.pdf(xy, mean_vec, cov_mtx) # x,y軸、平均、分散共分散

    # プロット
    surfacecolor = 'OrRd' # coolwarm YlOrRd
    linecolor = 'binary'
    fig = plt.figure()
    
    # ax: 3次元グラフの描画
    ax = fig.add_subplot(111, projection='3d')
    ax.plot_surface(x, y, z, rstride=1, cstride=1, cmap=surfacecolor)  
    ax.contour(x, y, z, zdir='x', offset=min(side), levels=1, cmap=linecolor)
    ax.contour(x, y, z, zdir='y', offset=max(side), levels=1, cmap=linecolor)
    ax.set_xlabel('$X$')
    ax.set_ylabel('$Y$')
    ax.set_zlabel('$f(x,y)$', labelpad=1)

    plt.show()

■等高線グラフ

  • 入力パラメータは「平均ベクトルmean_vec」と「共分散行列cov_mtx」です。

  • multivariate_normal.pdfで２変量正規分布の同時確率密度関数の値を取得します。

  • contourで等高線をプロットします。

  • contoufrで等高線の間を塗りつぶします。

  • colorbarで等高線グラフの右横のカラーバーを表示します。

# 2変量正規分布の等高線グラフを描画する関数
def bi_norm_plot2(mean_vec, cov_mtx):
    # x,y軸の設定
    side = np.arange(-20, 20, 0.5)
    x, y = np.meshgrid(side, side)
    xy = np.dstack((x, y))

    # z軸：2変量正規分布の確率密度関数の算出
    z = multivariate_normal.pdf(xy, mean_vec, cov_mtx) # x,y軸、平均、分散共分散

    # プロット
    surfacecolor = 'OrRd' # coolwarm
    linecolor = 'black'
    fig = plt.figure(figsize=(4, 3))
    
    # ax: 等高線グラフの描画
    ax = fig.add_subplot(111)
    c = ax.contour(x, y, z, levels=10, colors=linecolor, linewidths=0.5)
    c = ax.contourf(x, y, z, cmap=surfacecolor)
    fig.colorbar(c, ax=ax)
    ax.set_xlabel('$X$')
    ax.set_ylabel('$Y$')
    fig.tight_layout()

    plt.show()

③グラフのプロット
パラメータのmean_vecとcov_mtxの値を変えてみてください。
なお、パラメータの妥当性チェックをしていないので、誤った値を設定するとエラーになります。

■パターン１

# 2変量正規分布のパラメータ設定：平均ベクトル、分散共分散行列
mean_vec = [0, 0]
cov_mtx = [[30, 0], 
           [0, 30]]
# プロット
bi_norm_plot(mean_vec, cov_mtx)

インタラクティブモードの場合、つぎのような画面が表示されます。
マウスでグリグリ回しましょう！

回転の例です。
この状態でインタラクティブモードを終了すると、この状態で画像が残ります。

bi_norm_plot2(mean_vec, cov_mtx)
出力イメージ

■パターン２

# 2変量正規分布のパラメータ設定：平均ベクトル、分散共分散行列
mean_vec = [0, 0]
cov_mtx = [[30, 0], 
           [10, 30]]
# プロット
bi_norm_plot(mean_vec, cov_mtx)
出力イメージ
bi_norm_plot2(mean_vec, cov_mtx)
出力イメージ

■パターン３

# 2変量正規分布のパラメータ設定：平均ベクトル、分散共分散行列
mean_vec = [0, 0]
cov_mtx = [[10, 0], 
           [10, 30]]
# プロット
bi_norm_plot(mean_vec, cov_mtx)
出力イメージ
bi_norm_plot2(mean_vec, cov_mtx)
出力イメージ

■パターン４

# 2変量正規分布のパラメータ設定：平均ベクトル、分散共分散行列
mean_vec = [0, 0]
cov_mtx = [[1, 0], 
           [0, 30]]
# プロット
bi_norm_plot(mean_vec, cov_mtx)
出力イメージ
bi_norm_plot2(mean_vec, cov_mtx)
出力イメージ


Pythonサンプルファイルのダウンロード
こちらのリンクからJupyter Notebook形式のサンプルファイルをダウンロードできます。


おわりに

今回のグラフ作成はとても楽しく準備することができました。
１日程度のWeb調査でやりたいことを実現する方法をたくさん見つけられたからです。
特に３次元グラフと周辺分布の２次元グラフを一緒に表示する方法を知ったときは嬉しかったです。
新しい知見を獲得できること＝幸せです。

最後までお読みいただきまして、ありがとうございました。

のんびり統計シリーズの記事

次の記事

前の記事

目次

いいなと思ったら応援しよう！