私は、やや型にはまらない事前分布を持つ単純な二変量正規モデルを使用しています。私が抱えている主な問題は、実行ごとに事後分布が一貫していないことです。これは、連続したサンプル間の依存度が高い問題に関連していると推測しています。ここに私の具体的な質問があります。
N個の独立したサンプルを取得する最良の方法は何ですか? 現時点では、sample() を呼び出して大きなチェーン (たとえば、長さ 10,000) を取得し、1,000 から始まる 100 番目のサンプルごとに取得しています。しかし、パラメータの 1 つの自己相関プロファイルを見ると、少なくとも 500 番目のサンプルごとに取得する必要があるようです! (相互情報量を使用して、ラグ間の依存関係をよりよく理解することもできます。)
私は、pymc3 チュートリアルの確率的ボラティリティの例で説明されているフィッティング手順に従っています。特に、最初に MAP を見つけ、それを使用して NUTS() オブジェクトを生成し、次に短いサンプルを取得し、それを使用して別の NUTS() オブジェクトを生成し、gamma=0.25 (???) を使用して、最後に私の大きなサンプル。これが適切かどうか、またはガンマ = 0.25 が必要かどうかはわかりません。
また、同じ例では、指数分布のテスト値があります。これらが必要かどうかわかりません。(平均値のデフォルト使用の何が問題なのですか?)
実際に使っているモデルはこちら。
import pymc3 as pymc
import numpy as np
import theano.tensor as th
from pymc3.distributions.continuous import Gamma, Uniform, Normal, Bounded
from pymc3.distributions.multivariate import MvNormal
from pymc3.model import Deterministic
data = np.random.randn(3000, 2) / 300 # I have actual data!
with pymc.Model():
tau = Gamma('tau', alpha=2, beta=1 / 20000)
sigma = Deterministic('sigma', 1 / th.sqrt(tau))
corr = Uniform('corr', lower=0, upper=1)
alpha_sig = Deterministic('alpha_sig', sigma / 50)
alpha_post = Normal('alpha_post', mu=0, sd=alpha_sig)
alpha_pre = Bounded(
'alpha_pre', Normal, alpha_post, np.Inf, mu=0, sd=alpha_sig)
corr_inv = th.stack([th.stack([1, -corr]),
th.stack([-corr, 1])]) / (1 - th.sqr(corr))
MvNormal(
'data', mu=th.stack([alpha_post, alpha_pre]),
tau=tau * corr_inv, observed=data)
map_ = pymc.find_MAP()
step1 = pymc.NUTS(scaling=map_)
trace1 = pymc.sample(1000, step=step1)
step2 = pymc.NUTS(scaling=trace1[-1], gamma=0.25)
trace2 = pymc.sample(10000, step=step2, start=trace1[-1])