標準正規分布と1/(1+Z^2)の期待値

はじめに

標準正規分布は数理統計学における最も基本的な確率分布の1つであり，また統計検定でもおなじみの確率分布である。

確率分布において興味があるものに，期待値が挙げられる。
標準正規分布に関する期待値は，久保川「現代数理統計学の基礎」の演習問題に多数紹介されているが，本記事においては，標準正規分布 $Z \sim \mathcal{N}(0, 1)$ にしたがう確率変数に対して，次のような期待値を考える(ただし $m$ は正の整数)。

$\begin{align} A_m = E \left[ \frac{1}{(1+Z^2)^m} \right] \end{align}$

一見するとただの積分であるが，計算方法がかなり興味深かったので，本記事において紹介する。

はじめに
- 目次
数理統計学における位置付け
m=1の場合
- 基本方針 : ラプラス変換
- 期待値の計算
m=2の場合
- 基本方針 : スタインの等式
- 期待値の計算
  - g(z)の選定
  - スタインの等式の適用
m≧3の場合
- 基本方針 : スタインの等式
- 期待値の計算
  - g(z)の選定
Pythonによる数値計算
まとめと感想

数理統計学における位置付け

期待値 $A_m = E[ (1+Z^2)^{-m} ]$ は，「多項式の逆数」のモーメントを評価する問題になる。また， $Z \sim \mathcal{N}(0, 1)$ から $Z^2$ は自由度1のカイ二乗分布 $\chi^2(1)$ にしたがうため，標準正規分布だけでなく，カイ二乗分布にも関連する問題である。

特に $m=1, 2$ の場合は，比較的シンプルな結果が得られる。

$\begin{align} &A_1 = E \left[ \frac{1}{1+Z^2} \right] = \frac{1 - \Phi(1)}{\phi(1)} \\ \\ &A_2 = E \left[ \frac{1}{(1+Z^2)^2} \right] = \frac{1}{2} \\ \\ \end{align}$

ただし $\Phi(z), \phi(z)$ はそれぞれ，標準正規分布の累積分布関数と確率密度関数である。

m=1の場合

基本方針 : ラプラス変換

$m=1$ ，すなわち

$\begin{align} A_1 = E \left[ \frac{1}{1+Z^2} \right] \\ \\ \end{align}$

の場合を考える。

この式は分母が多項式になっている形になるため，そのままでは積分が扱いづらい。したがって，以下のようなラプラス変換を考える。

$\begin{align} \frac{1}{1+Z^2} = \int_0^{\infty} e^{-t(1+z^2)} dt \\ \\ \end{align}$

ややトリッキーな変換ではあるが，以下のようなメリットがある。

$1+z^2$ が分母から指数に移動するため，分母に多項式がある状態を回避できる。
指数関数が登場するので，積分の中に出てくる正規分布(またはカイ二乗分布)と相性がよい。

期待値の計算

確率変数の変換

この後の計算のことを考慮し， $Z^2 =Y$ という変数変換を行なう。この変換により $Y$ は自由度1のカイ2乗分布分布にしたがう。

求める期待値を $Y$ の関数で表すことを考える。確率変数 $Z$ による期待値であることを明示的に表すために，期待値を $E_Z[ \cdot ]$ のように表す。
被積分関数は偶関数であるため，

$\begin{align} E_Z \left[ \frac{1}{1+Z^2} \right] &= \int_{\infty}^{\infty} \frac{1}{1+z^2} \cdot \frac{1}{\sqrt{2 \pi}} e^{-z^2/2} dz \\ \\ &= \int_{0}^{\infty} \frac{1}{1+z^2} \cdot \frac{2}{\sqrt{2 \pi}} e^{-z^2/2} dz \\ \\ \end{align}$

となる。

ここで， $Z^2 = Y$ という変数変換をすると， $dz = dy/2 \sqrt{y}$ となる。また積分範囲は $0 \leq y \lt \infty$ なので，

$\begin{align} E_Z \left[ \frac{1}{1+Z^2} \right] &= \int_{0}^{\infty} \frac{1}{1+y} \cdot \frac{2}{\sqrt{2 \pi}} \frac{1}{2 \sqrt{y}}e^{-y/2} dy \\ \\ &= \int_{0}^{\infty} \frac{1}{1+y} \cdot \frac{1}{\sqrt{2 \pi}} y^{-1/2} e^{-y/2} dy \\ \\ \therefore E_Z \left[ \frac{1}{1+Z^2} \right] &= E_Y \left[ \frac{1}{1+Y} \right] \\ \\ \end{align}$

となる。

積分順序の入替え

ラプラス変換と確率変数の変換を用いて，期待値の計算を行なう。

$\begin{align} E_Z \left[ \frac{1}{1+Z^2} \right] &= E_Y \left[ \frac{1}{1+Y} \right] \\ \\ &= E_Y \left[ \int_0^{\infty} e^{-t(1+y)} dt \right] \\ \\ &= \int_{0}^{\infty} dy \left( \int_0^{\infty} dt \; e^{-t(1+y)} \right) \cdot \frac{1}{\sqrt{2 \pi}} y^{-1/2} e^{-y/2} \\ \\ \end{align}$

ここで，被積分関数は非負であるため，積分順序の入れ替えが可能になる(※清水「統計学への確率論，その先へ」のP63 (2.2.3項)参照)。したがって，

$\begin{align} E_Z \left[ \frac{1}{1+Z^2} \right] &= \int_0^{\infty} dt \; e^{-t} \left( \int_{0}^{\infty} dy \; e^{-ty} \cdot \frac{1}{\sqrt{2 \pi}} y^{-1/2} e^{-y/2} \right) \\ \\ &= \int_0^{\infty} e^{-t} E[e^{-ty} ] dt \\ \\ \end{align}$

となる。

モーメント母関数を用いた期待値計算

$Y \sim \chi^2(1)$ であり，このモーメント母関数は $E[e^{ty} ] = M_Y(t) = (1-2t)^{-1/2}$ であるため，

$\begin{align} E[e^{-ty} ] = E[e^{(-t)y} ] = (1-2(-t))^{-1/2} = (1+2t)^{-1/2} \\ \\ \end{align}$

となる。よって，

$\begin{align} E_Z \left[ \frac{1}{1+Z^2} \right] = \int_0^{\infty} \frac{e^{-t}}{\sqrt{1+2t}} dt \\ \\ \end{align}$

となる。

$\sqrt{1+2t} = u$ とすると， $t = (u^2-1) /2$ であり， $dt = udu$ ，また $1 \leq u \lt \infty$ となるので，

$\begin{align} E_Z \left[ \frac{1}{1+Z^2} \right] &= \int_1^{\infty} e^{1/2} \cdot \frac{1}{u} e^{-u^2/2} \cdot u du \\ \\ &= e^{1/2} \int_1^{\infty} e^{-u^2/2} du \\ \\ &= \sqrt{2 \pi e} \times \int_1^{\infty} \frac{1}{\sqrt{2 \pi}}e^{-u^2/2} du \\ \\ &= \sqrt{2 \pi e}(1 - \Phi(1)) \\ \\ \end{align}$

となる。

標準正規分布の確率密度関数 $\phi(z)$ について， $\phi(1) = 1/\sqrt{2 \pi e}$ が成り立つので，最終的に

$\begin{align} E_Z \left[ \frac{1}{1+Z^2} \right] = \frac{1 - \Phi(1)}{\phi(1)} \\ \\ \end{align}$

が得られた。

なおこのような，生存関数と確率密度関数の比をミルズ比(Mills ratio)と呼ぶ。

m=2の場合

基本方針 : スタインの等式

$m=2$ ，すなわち

$\begin{align} A_2 = E \left[ \frac{1}{(1+Z^2)^2} \right] \\ \\ \end{align}$

の場合を考える。

$m=1$ の時と比べて，分母の次数が大きくなっている。そのため，次数を変えることを考える。正規分布にしたがう確率変数の期待値について，次数にかかわる公式として，スタインの等式(Stein identity)が挙げられる。

久保川「現代数理統計学の基礎」のP48(3.3.1項)を参照すると，標準正規分布版のスタインの等式は以下のようになる。

$Z \sim \mathcal{N}(0, 1)$ であり， $g(\cdot)$ が微分可能で $E[ \lvert g'(Z) \rvert ] \lt \infty$ のとき，次の等式が成り立つ。

$\begin{align} E[ Z g(Z) ] = E[ g'(Z) ] \end{align}$

この等式を用いると，被積分関数の次数を1つ下げて，微分で置き換えるということが可能になる。

期待値の計算

g(z)の選定

スタインの等式に当てはめる $g(z)$ を考える。
求めたい期待値 $A_2$ の形を見ると，分母に $(1+z^2)^2$ がある。 $1/(1+z^2)$ を微分すると，分母に $(1+z^2)^2$ が登場し，求めたい期待値 $A_2$ に近づくと考えられるので，ためしに，

$\begin{align} g(z) = \frac{1}{1+z^2} \\ \\ \end{align}$

と置いてみる。すると，

$\begin{align} g’(z) = -\frac{2z}{(1+z^2)^2} \\ \\ \end{align}$

となる。これでは分子に $z$ の1次式が出てしまい，求めたい $A_2$ の形に近づけることが難しくなる。

次に，

$\begin{align} g(z) = \frac{z}{1+z^2} \\ \\ \end{align}$

と置いてみる。すると，

$\begin{align} g’(z) = \frac{1-z^2}{(1+z^2)^2} &= \frac{1}{(1+z^2)^2} - \frac{z^2}{(1+z^2)^2} \\ \\ &= \frac{1}{(1+z^2)^2} - \frac{1 + z^2 -1 }{(1+z^2)^2} \\ \\ &= \frac{2}{(1+z^2)^2} - \frac{1}{1+z^2} \\ \\ \end{align}$

となり，分子に $z$ の1次式がなくなり，求めたい $A_2$ に近づいている。

スタインの等式の適用

$g(z)$ の選定にめどが立ったので，スタインの等式

$\begin{align} E[ Z g(Z) ] = E[ g'(Z) ] \\ \\ \end{align}$

を適用してみる。

まず左辺の被積分関数は，

$\begin{align} zg(z) = \frac{z^2}{1+z^2} = 1- \frac{1}{1+z^2} \\ \\ \end{align}$

となる。

次に右辺の被積分関数は，

$\begin{align} g’(z) = \frac{2}{(1+z^2)^2} - \frac{1}{1+z^2} \\ \\ \end{align}$

となる。

スタインの等式 $E[ Z g(Z) ] = E[ g'(Z) ]$ より，

$\begin{align} &1- E \left[ \frac{1}{1+z^2} \right] = 2E \left [ \frac{1}{(1+z^2)^2} \right ] - E \left [ \frac{1}{1+z^2} \right ] \\ \\ \Leftrightarrow & 1 = 2E \left [ \frac{1}{(1+z^2)^2} \right ] \\ \\ \therefore & A_2 = E \left [ \frac{1}{(1+z^2)^2} \right ] = \frac{1}{2} \\ \\ \end{align}$

が得られた。

m≧3の場合

基本方針 : スタインの等式

$m \geq 3$ ，の場合

$\begin{align} A_m = E \left[ \frac{1}{(1+Z^2)^m} \right] \\ \\ \end{align}$

においても， $m=2$ のときと同様に，スタインの等式を用いる。

期待値の計算

g(z)の選定

$m=2$ のときにならい，

$\begin{align} g(z) = \frac{z}{(1+z^2)^m} \\ \\ \end{align}$

と置いてみる。

すると，

$\begin{align} g’(z) = \frac{1-2m}{(1+z^2)^m} + \frac{2m}{(1+z^2)^{m+1}} \\ \\ \end{align}$

となる。また，

$\begin{align} zg(z) = \frac{z^2}{(1+z^2)^m} = \frac{1}{(1+z^2)^{m-1}} - \frac{1}{(1+z^2)^m} \\ \\ \end{align}$

となる。

スタインの等式 $E[ Z g(Z) ] = E[ g'(Z) ]$ より，

$\begin{align} & E \left[ \frac{1}{(1+z^2)^{m-1}} \right] - E \left[ \frac{1}{(1+z^2)^m} \right] = (1-2m) E \left[ \frac{1}{(1+z^2)^m} \right] + 2m E \left[ \frac{1}{(1+z^2)^{m+1}} \right] \\ \\ \therefore \; &A_{m-1} - A_m = (1-2m) A_m + 2m A_{m+1} \\ \\ \therefore \; &A_{m+1} = \frac{2(m-1)A_m + A_{m-1}}{2m} \\ \\ \end{align}$

という3項間漸化式が得られる。 $A_1, A_2$ は求められているので， $A_3, A_4, ...$ は順次求めることができる。

Pythonによる数値計算

$A_1, A_2$ について，期待値を数値計算した結果と，解析解を比較するPythonコードを作成した。

サンプルコード

クリックで展開

import scipy.integrate as integrate
import scipy.stats as stats
import numpy as np

# 標準正規分布の確率密度関数 phi(z)
def phi(z):
    return stats.norm.pdf(z)

# 標準正規分布の累積分布関数 Phi(Z)
def Phi(z):
    return stats.norm.cdf(z)

# 期待値 E[1/(1+Z^2)] を計算するための関数 f(z) = 1/(1+z^2) * phi(z)
def integrand(z):
    return (1 / (1 + z**2)) * phi(z)

# 期待値 E[1/(1+Z^2)^2] を計算するための関数 f(z) = 1/(1+z^2)^2 * phi(z)
def integrand2(z):
    return (1 / (1 + z**2)**2) * phi(z)


# m=1の場合 ###############
# 数値積分により期待値を計算
# 標準正規分布は無限の範囲を持つため，十分広い範囲で積分する
expected_value1, error1 = integrate.quad(integrand, -np.inf, np.inf)

# A_1 = (1 - Phi(1)) / phi(1) を計算
A1 = (1 - Phi(1)) / phi(1)

# 出力
print()
print(f"m=1 : ")
print(f"期待値 E[1/(1+Z^2)] の数値積分結果: {expected_value1:.6f}")
print(f"A_1 の解析解の結果: {A1:.6f}")

# m=2の場合 ###############
expected_value2, error2 = integrate.quad(integrand2, -np.inf, np.inf)

# A_2 = 1/2 を計算
A2 = 1/2

# 出力
print()
print(f"m=2 : ")
print(f"期待値 E[1/(1+Z^2)^2] の数値積分結果: {expected_value2:.6f}")
print(f"A_2 の解析解の結果: {A2:.6f}")

出力

m=1 : 
期待値 E[1/(1+Z^2)] の数値積分結果: 0.655680
A_1 の解析解の結果: 0.655680

m=2 : 
期待値 E[1/(1+Z^2)^2] の数値積分結果: 0.500000
A_2 の解析解の結果: 0.500000

以上のように，数値積分の結果と解析解の結果がほぼ一致することが確認できた。

まとめと感想

本記事で取り上げた $A_m = E[ (1+Z^2)^{-m} ]$ という期待値計算は，一見すると「積分計算の練習問題」に思えるが，実際には統計学の広い領域とつながっている題材であることが確認できた。

まず，標準正規分布から導かれる $Z^2 \sim \chi^2(1)$ という関係を通じて，分布の変換やモーメント母関数を扱う力が試される。またスタインの等式によって漸化式が得られる点は，「統計的推定の基礎にある道具が，積分評価にも応用できる」ことを示していて興味深い。

統計検定1級の学習を通じて強く感じたのは，「計算技術」だけでなく「分布の性質をどう解釈するか」が問われるという点である。今回の問題もその典型で，ただ積分を解くだけでなく，その過程や結果が統計学における他の話題とどう関係するか，ということを理解することが重要であると考える。統計を学ぶうえで，このような積分問題に挑戦し続けることは，自分の力を確実に底上げしてくれると感じた。