Table of Contents
標準正規分布に従う確率変数$X$を$\exp(X)$に変換したとき
Problem
$X\sim N(0, 1)$のとき, $Y = \exp(x)$という変数変換を考えたとき,$Y$の密度関数, 平均, 分散を求めよ
$Y = \exp(x)$は単調増加 & 微分可能な関数による変換であることに留意します.
密度関数の導出
$X, Y$の密度関数をそれぞれ$f(x), g(y)$とすると
\[\begin{align*} g(y) = f(\ln(y))\frac{1}{y} \end{align*}\]$X$は標準正規分布に従うので
\[g(y) = \frac{1}{\sqrt{2\pi}}\exp\bigg(-\frac{[\ln(y)]^2}{2}\bigg)\frac{1}{y}\]平均の導出
$Y$の定義域が$(0, \infty)$であることに留意すると
\[\begin{align*} \mathbb E[Y] &= \int_0^\infty y\frac{1}{\sqrt{2\pi}}\exp\bigg(-\frac{[\ln(y)]^2}{2}\bigg)\frac{1}{y} dy\\[3pt] &= \int_0^\infty \frac{1}{\sqrt{2\pi}}\exp\bigg(-\frac{[\ln(y)]^2}{2}\bigg) dy \end{align*}\]ここで$\ln(y) = t$と変数変換すると
\[\begin{align*} &\int_0^\infty \frac{1}{\sqrt{2\pi}}\exp\bigg(-\frac{[\ln(y)]^2}{2}\bigg) dy\\[3pt] &= \int_{-\infty}^\infty \frac{1}{\sqrt{2\pi}}\exp\bigg(-\frac{t^2}{2}\bigg) \frac{dy}{dt}dt\\[3pt] &= \int_{-\infty}^\infty \frac{1}{\sqrt{2\pi}}\exp\bigg(-\frac{t^2}{2}\bigg) \exp(t) dt\\[3pt] &= \int_{-\infty}^\infty \frac{1}{\sqrt{2\pi}}\exp\bigg(-\frac{t^2-2t}{2}\bigg) dt\\[3pt] &= \int_{-\infty}^\infty \frac{1}{\sqrt{2\pi}}\exp\bigg(-\frac{(t-1)^2-1}{2}\bigg) dt\\[3pt] &= \exp\bigg(\frac{1}{2}\bigg)\int_{-\infty}^\infty \frac{1}{\sqrt{2\pi}}\exp\bigg(-\frac{(t-1)^2}{2}\bigg) dt\\[3pt] &= \exp\bigg(\frac{1}{2}\bigg) \end{align*}\]分散の導出
\[V(Y) = \mathbb E[Y^2] - \mathbb E[Y]^2\]で定義されるので$\mathbb E[Y^2]$を求めればよい.
\[\begin{align*} \mathbb E[Y^2] &= \int_{-\infty}^\infty \frac{1}{\sqrt{2\pi}}\exp\bigg(-\frac{t^2}{2}\bigg) \exp(2t) dt\\[3pt] &= \int_{-\infty}^\infty \frac{1}{\sqrt{2\pi}}\exp\bigg(-\frac{t^2-4t}{2}\bigg)dt\\[3pt] &= \int_{-\infty}^\infty \frac{1}{\sqrt{2\pi}}\exp\bigg(-\frac{(t-2)^2-4}{2}\bigg)dt\\[3pt] &= \exp(2)\int_{-\infty}^\infty \frac{1}{\sqrt{2\pi}}\exp\bigg(-\frac{(t-2)^2}{2}\bigg)dt\\[3pt] &= \exp(2) \end{align*}\]したがって,
\[\begin{align*} V(Y) &= \mathbb E[Y^2] - \mathbb E[Y]^2\\ &= \exp(2) - \exp(1)\\ &= e(e-1) \end{align*}\]Pythonで確認してみる
変数変換後のpdfに基づいた分布検定は少しめんどくさいので, 平均とthe second-momentの分布を乱数シミュレーションで確認する.
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
#%%
import numpy as np
import plotly.express as px
np.random.seed(42)
N = 1000
LOOP = 10000
X = np.random.normal(0, 1, (N, LOOP))
Y = np.exp(X)
y_mean = np.mean(Y, axis=0)
y_second_mean = np.mean(Y**2, axis=0)
fig1 = px.histogram(y_mean,
histnorm='probability density',
title="""Histogram of the mean of exp(x), true-val: {:.2f}, sample-mean: {:.2f}
""".format(np.exp(0.5), np.mean(y_mean)))
fig1.add_vline(x=np.exp(0.5), line_dash="dash", annotation_text="mean")
fig1.show()
#%%
fig2 = px.histogram(y_second_mean,
histnorm='probability density',
title="""Histogram of the second moment of exp(x), true-val: {:.2f}, sample-mean: {:.2f}
""".format(np.exp(2), np.mean(y_second_mean)))
fig2.add_vline(x=np.exp(2), line_dash="dash", annotation_text="the second moment")
fig2.show()
平均の分布
the second-momentの分布
テイラー展開で変数変換後確率変数の期待値を近似的に求める
Problem
$X\sim N(0, 1)$のとき, $Y = \exp(x)$という変数変換を考えたとき,$Y$の平均をテイラー近似で求めよ
一般に, 確率変数 $X$ の期待値を$\theta$としたとき, 関数$g(X)$の近似的な期待値はデルタ法で求めることができます. $\theta$まわりでの$g(X)$の2次のテイラー展開
\[g(x)\approx g(\theta) + g(\theta)(x-\theta) + \frac{1}{2}g''(\theta)(x-\theta)^2\]の両辺について期待値を取ることで求めることができます.
k次近似を$g_k(x)$と表すとしてまず, 2次近似で考えてみると
\[\mathbb E[g_2(x)] = \exp(0) + \frac{1}{2} = 1.5\]$n$を正の整数とし, 今回は $X$ は標準正規分布に従っているので$\mathbb E[X^{2n+1}] = 0$で有ることに留意すると, $2n$近似のみを考えれば良い. $2n$近似についてAppendixより
\[\begin{align*} \mathbb E[g_{2n}(x)] &= 1 + \frac{1}{2} + \frac{1}{2}\frac{1}{4}+ \frac{1}{2}\frac{1}{4}\frac{1}{6} + \cdots + \frac{1}{2}\cdots\frac{1}{2(n-1)}\frac{1}{2n}\\[3pt] &= 1 + \sum_{j=1}^{n}\bigg(\prod_{k=1}^j \frac{1}{2k}\bigg) \end{align*}\]ここで, 右辺第二項についてWolfram alphaで計算してみると
\[\lim_{n\to\infty}\sum_{j=1}^{n}\bigg(\prod_{k=1}^j \frac{1}{2k}\bigg)=\exp\bigg(\frac{1}{2}\bigg)-1\]従って,
\[\lim_{k\to\infty}g_k(x) = \exp\bigg(\frac{1}{2}\bigg)\]となり上記で求めた答えと一致します.
近似誤差が$n$が大きくなるに従ってどのように変化するか？
近似誤差のオーダーについてPythonで確認すると以下のようになります
Python code
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
import math
import numpy as np
import plotly.express as px
def taylor_term(n):
return np.sum(list(map(lambda x: (1/2)**x * 1/math.gamma(x+1), np.arange(1, n+1))))
x = 1 + np.array(list(map(taylor_term, np.arange(1, 10))))
y = np.exp(1/2)
fig = px.line(x=np.arange(1, 20),
y=(y-np.r_[1, np.repeat(x, 2)])/y,
title='テイラー近似とTrue-valueのpercentage error: 6次近似で誤差は0.2%を下回る')
fig.update_xaxes(title_text="近似レベル")
fig.update_yaxes(title_text="percentage error")
fig.update_traces(hovertemplate='テイラー展開%{x}次近似<br>percentage error: %{y:.5%}')
Appendix
変数変換のルール
確率変数 $X$の密度関数を $f(x)$, $y = h(x)$を単調増加で, 微分可能な関数とする. このとき, $Y = h(X)$の密度関数 $g(y)$を表せ.
$X, Y$のcdfをそれぞれ$F(x), G(y)$とする.
\[\begin{align*} G(y) &= \Pr(Y \leq y)\\ &= \Pr(h(X) \leq y)\\ &= \Pr(X \leq h^{-1}(y))\\ &= F(h^{-1}(y)) \end{align*}\]両辺を$y$で微分するとこのとき
\[g(y) = f(h^{-1}(y))\frac{d}{dy}h^{-1}(y)\]これで$g(y)$が表された.
標準正規分布のn次モーメントについて
$l$を自然数として$n=2l$について
\[\begin{align*} \mathbb E[X^n] &= \frac{1}{\sqrt{2\pi}}\int^{\infty}_{-\infty}x^n\exp\bigg(-\frac{x^2}{2}\bigg)dx\\[3pt] &= \frac{2}{\sqrt{2\pi}}\int^{\infty}_{0}x^{2l}\exp\bigg(-\frac{x^2}{2}\bigg)dx\\[3pt] \end{align*}\]$x^2/2 = u$として置換積分を行うと
\(\begin{align*} \mathbb E[X^n] &= \frac{2^{l+1}}{\sqrt{2\pi}}\int^{\infty}_{0}u^{l}\exp(-u)\frac{1}{\sqrt{2u}}du\\[3pt] &= \frac{2^l}{\sqrt{\pi}}\int^{\infty}_{0}u^{l-\frac{1}{2}}\exp(-u)du\\[3pt] &= \frac{2^l}{\sqrt{\pi}}\int^{\infty}_{0}u^{l+\frac{1}{2}-1}\exp(-u)du\\[3pt] &= \frac{2^l}{\sqrt{\pi}}\Gamma\bigg(l+\frac{1}{2}\bigg) \end{align*}\)
$\Gamma(1/2) = \sqrt{\pi}$であることに留意すると
\[\mathbb E[X^n] = 1\times 3\times 5\times\cdots\times (2l-1)\]なお,
\[\frac{1}{2l!}\mathbb E[X^{2l}] = \prod^l_{k=1}\frac{1}{2k}\]となる.
References
(注意：GitHub Accountが必要となります）