pymc3の推定パラメータからの予測の生成

Question

私はよくある問題にぶつかります。私はしばしば、訓練（すなわち、実際にパラメータの推論を実行する）と評価（すなわち、予測を生成するために推論されたパラメータを使用する）の2つのモードでpymc3を使いたいと思う。

一般的に、私は、ポイントベースの推定だけでなく、ベイジアンフレームワークのメリットの一部でもあります。トレーニングデータが固定されている場合、これは通常、類似した形式のシミュレート変数を観測変数に追加することによって行われます。たとえば、

from pymc3 import * 

with basic_model: 

    # Priors for unknown model parameters 
    alpha = Normal('alpha', mu=0, sd=10) 
    beta = Normal('beta', mu=0, sd=10, shape=2) 
    sigma = HalfNormal('sigma', sd=1) 

    # Expected value of outcome 
    mu = alpha + beta[0]*X1 + beta[1]*X2 

    # Likelihood (sampling distribution) of observations 
    Y_obs = Normal('Y_obs', mu=mu, sd=sigma, observed=Y) 
    Y_sim = Normal('Y_sim', mu=mu, sd=sigma, shape=len(X1)) 

    start = find_MAP() 
    step = NUTS(scaling=start) 
    trace = sample(2000, step, start=start) 

データが変更されたらどうなりますか？新しいデータに基づいて予測を生成したいが、推論を繰り返し実行する必要はないとします。理想的には、私はpredict_posterior(X1_new, X2_new, 'Y_sim', trace=trace)またはpredict_point(X1_new, X2_new, 'Y_sim', vals=trace[-1])のような関数を持っていて、単にtheano計算グラフを通して新しいデータを実行するでしょう。

私の質問の一部は、pymc3がtheano計算グラフをどのように実装するかと関係していると思います。私は関数model.Y_sim.evalが欲しいものに似ているように見えましたが、それは入力としてY_simを必要とし、あなたがそれを与えるものを返すようです。

私はこのプロセスが非常に一般的だと思いますが、私はそれを行う方法を見つけることができないようです。どんな助けでも大歓迎です。 （これはpymc2でこれを行うためのハックを持っているので、pymc3ではtheanoのために難しいです。）

Answer 1

注：この機能は、pymc.sample_ppcメソッドとしてコアコードに組み込まれました。詳細はthe docsをご覧ください。

これに基づいて、（2017年7月現在で死んでいる）がtwieckiから私に送信されました。私の問題を解決するためにいくつかのトリックがあります。 1つ目は、トレーニングデータを共有するtheano変数に入れることです。これにより、theano計算グラフを壊すことなく、後でデータを変更することができます。

X1_shared = theano.shared(X1) 
X2_shared = theano.shared(X2) 

次に、モデルを構築し、通常のように推論を実行しますが、共有変数を使用します。

with basic_model: 

    # Priors for unknown model parameters 
    alpha = Normal('alpha', mu=0, sd=10) 
    beta = Normal('beta', mu=0, sd=10, shape=2) 
    sigma = HalfNormal('sigma', sd=1) 

    # Expected value of outcome 
    mu = alpha + beta[0]*X1_shared + beta[1]*X2_shared 

    # Likelihood (sampling distribution) of observations 
    Y_obs = Normal('Y_obs', mu=mu, sd=sigma, observed=Y) 

    start = find_MAP() 
    step = NUTS(scaling=start) 
    trace = sample(2000, step, start=start) 

最後に、新しいデータのための事後を予測することができます開発中の機能は、（おそらく最終的にpymc3に追加されます）があります。

X1_shared.set_value(X1_new) 
X2_shared.set_value(X2_new) 

最後に、あなたは新しいデータのための事後予測サンプルを生成することができます

from collections import defaultdict 

def run_ppc(trace, samples=100, model=None): 
    """Generate Posterior Predictive samples from a model given a trace. 
    """ 
    if model is None: 
     model = pm.modelcontext(model) 

    ppc = defaultdict(list) 
    for idx in np.random.randint(0, len(trace), samples): 
     param = trace[idx] 
     for obs in model.observed_RVs: 
      ppc[obs.name].append(obs.distribution.random(point=param)) 

    return ppc 

次に、あなたが予測を実行する新しいデータを渡します。

ppc = run_ppc(trace, model=model, samples=200) 

変数ppcは、モデル内の各観測変数のキーを持つ辞書です。したがって、この場合、ppc['Y_obs']には配列のリストが含まれ、それぞれの配列はトレースからの単一のパラメータセットを使用して生成されます。

トレースから抽出したパラメータを変更することもできます。たとえば、GaussianRandomWalk変数を使用するモデルがあり、将来の予測を生成したいと考えていました。将来的にpymc3をサンプリングできるようにすることができます（つまり、ランダムウォーク変数を発散させることができます）。最後の推定値に対応する係数の固定値を使用したかっただけです。このロジックは、run_ppc機能で実装できます。

run_ppc機能が非常に遅いことにも言及する価値があります。実際の推論を実行するのに多くの時間がかかります。私はこれがtheanoがどのように使用されているかに関連したいくつかの非効率性と関係があると考えています。

EDIT：元々含まれていたリンクが死んでいるようです。