#はじめに
千葉大学/Nospareの米倉です．今回はインポータンス・サンプリングについて解説していきたいと思います．
#モチベーション
統計学・機械学習ではしばしばある関数$f$に関して，その期待値$\int f(x)\mu(x)dx$を評価したいことがあります．ここで$\mu(x)$は密度関数だと理解してください．モンテカルロ法では$\mu(x)$から乱数$\xi$を例えば$N$個iidでサンプリングして，$\int f(x)\mu(x)dx$を$N^{-1}\sum_i f(\xi_i)$で近似します．この推定値は一致性や漸近正規性等良い性質があることは知られていて，実際良く使われています．
一方で，しばしば$\mu(x)$は複雑な形をしていて，それからのサンプリングが不可能な時があります．この様な状況の時に$\int f(x)\mu(x)dx$を評価するテクニックの一つがインポータンス・サンプリングです．
#基本的なインポータンス・サンプリング
インポータンス・サンプリングは次のような自明な関係を用いて導出されます．$$\int f(x)\mu(x)dx=\int f(x)\frac{\mu(x)}{\nu(x)}\nu(x)dx$$
ここで比$$\frac{\mu(x)}{\nu(x)}=w(x)$$
はしばしば重み関数と呼ばれ，$\nu(x)$はインポータンス密度と呼ばれます．重み関数がきちんと定義されるためには任意の$x$に対して，$\mu(x)>0$ならば$\nu(x)>0$が成立することが必要です．以上より，インポータンス・サンプリングを用いた$\int f(x)\mu(x)dx$の評価は，$\nu(x)$から乱数$\xi$を例えば$N$個iidでサンプリングして$N^{-1}\sum_i f(\xi_i)w(\xi_i)$で行うことが可能です．
$I:=\int f(x)\mu(x)dx$，$\hat{I}:=N^{-1}\sum_i f(\xi_i)w(\xi_i)$とすると，明らかに$\hat{I}$は$I$は不偏推定量となり，また強い意味での一致性があることも示せます．例えばGeweke(1989)を参照してください．
#最適なインポータンス密度
次に$\hat{I}$の分散を最小化する意味で，最適なインポータンス密度$\nu(x)$を導出します．まずコーシー・シュワルツの不等式より，$\int \mid f(x)\mid w(x)\nu(x)dx\leq(\int f(x)^2w(x)^2\nu(x)dx)^{1/2}$が成立するので，$\int\mid f(x)\mid\mu(x)dx\leq(\int f(x)^2w(x)^2\nu(x)dx)^{1/2}$を得ます．$\sigma^2 :=\int f(x)^2w(x)^2\nu(x)dx-I$と置くと，$(\int f(x)^2w(x)^2\nu(x)dx)^{1/2}=(\sigma^2+I^2)^{1/2}$と書き直せます．これより，$$(\int\mid f(x)\mid\mu(x)dx)^2-I^2\leq \sigma^2$$をえます．この時ある定数$c>0$が存在して，$\mid f(x)\mid\mu(x)=c\nu (x)$となれば，上の不等式は$0=\sigma^2$となり分散が最小化されます．密度関数の定義より，$\int\nu(x)dx=1$なので，結局$c=\int\mid f(x)\mid\mu(x)dx$となり，$\nu(x)=\frac{\mid f(x)\mid\mu(x)}{\int\mid f(x)\mid\mu(x)dx}$が最適となります．しかし問題設定（$\int f(x)\mu(x)dx$が評価できない）から明らかですが，最適なインポータンス密度は通常使えません．
#自己正規化されたインポータンス・サンプリング
統計学・機械学習の多くの問題では，残念ながら$\mu(x),\nu(x)$の正規化定数が不明なことが多いです．例えば事後分布を例にして考えると，尤度と事前分布の積からはサンプリングは行えるものの，その積分（正規化定数）は不明なことが多い等です．つまり，$Z_\mu,Z_\nu$をそれぞれの正規化定数とすると，$\mu(x)=\frac{\mu_u(x)}{Z_\mu},\nu(x)=\frac{\nu_u(x)}{Z_\nu}$のそれぞれの分子の項のみが計算可能な状況を指します．
この場合でも次の等式を用いてインポータンス・サンプリングを行うことが可能です．$$\int f(x)\mu(x)dx=\frac{\int f(x)\frac{\mu(x)}{\nu(x)}\nu(x)dx}{\int\frac{\mu(x)}{\nu(x)}\nu(x)dx}=\frac{\int f(x)\frac{\mu_u(x)}{\nu_u(x)}\nu(x)dx}{\int\frac{\mu_u(x)}{\nu_u(x)}\nu(x)dx}$$
つまり，乱数$\xi$を例えば$N$個iidでサンプリングして$W(\xi_i):=\frac{w(\xi_i)}{\sum_jw(\xi_j)}$と定義すれば，$\int f(x)\mu(x)dx$の評価は$N^{-1}\sum_i f(\xi_i)W(\xi_i)$で行うことが可能です．これをしばしば自己正規化されたインポータンス・サンプリングと呼びます．
不偏推定推定量の比は不偏推定量にならないので，自己正規化されたインポータンス・サンプリングは不偏性がありません．一方で強い意味での一致性はあります．
#数値実験
$f(x)=x^2$，$\mu(x)=\frac{1}{2}\exp(-\mid x\mid)$とし$\int f(x)\mu(x)dx$を計算することを考えます．これはラプラス分布の一種で，サンプリングを行うことが難しい分布の一つです．ここで$$\int x^2\frac{1}{2}\exp(-\mid x\mid)dx=\int x^2\frac{\frac{1}{2}\exp(-\mid x\mid)}{\frac{1}{\sqrt{8\pi}}\exp(-\frac{x^2}{8})}\frac{1}{\sqrt{8\pi}}\exp(-\frac{x^2}{8})dx$$と例えば書き直すことが出来ます．つまり$\nu(x)=N(0,4)$としました．これよりインポータンス・サンプリングを用いて$\int f(x)\mu(x)dx$を計算できます．この積分は解析的に行えて，答えは２です．下のコードはインポータンス・サンプリングを行うJuliaのコードです．
using Distributions
using Random
Random.seed!(1234)
N = 10000
x = rand(Normal(0,2),N)
weight = exp.(-abs.(x))./pdf(Normal(0,2),x)
Ihat = x.^2
Ihat = Ihat.*weight/2
mean(Ihat)
このコードで$\int f(x)\mu(x)dx$を計算すると，約2.02となり解析解とほぼ一致することが分かります．
#おわりに
株式会社Nospareでは統計学の様々な分野を専門とする研究者が所属しております．統計アドバイザリーやビジネスデータの分析につきましては株式会社Nospare までお問い合わせください．
