問題タブ [stan]

すでに devtools を使用してパッケージのスケルトンを作成してから、一連の R コード、メタデータ、ドキュメントなどを追加しました。このパッケージ内で rstan を使用したいと考えています。これを容易にするために rstan::rstan.package.skeleton がパッケージ スケルトンを作成することを理解しています。では、そのパッケージから rstan を使用するために必要な構造で既存のパッケージを拡張するためのベスト プラクティスは何ですか? ありがとうございました。

STANでスイッチポイント解析をしようとしています。yガウス確率変数の 2 つの異なるシーケンスを持つデータ ベクトルがあります。目標は、シフトが発生した可能性がある時期の事後分布を見つけることです。実行に使用RStanしていますが、エラーは STAN 内にあります。

これは STAN コードです。

パーサー (Rstudio に付属) は次のエラーを返します。

キャストを行う変数の割り当てを処理できないのはなぜですか? STAN では、この種の分析に別のパターンが必要ですか。で整数変数を作成しようとしましたparametersが、STAN はランダム整数変数をサポートしていないようで、連続変数のみをサポートしています。

観測データ y1 と y2 を考えます。y1 は連続スケールで測定され、y2 はバイナリ スケールで測定されます。連続潜在変数 z は、y2 = I(z > 0) として y2 を生成すると仮定されます。(z が正規の場合、y2 はわずかに 2 進プロビットです)。さらに、y1 と z の間の依存関係をモデル化するためにコピュラが使用されます。このモデルは、次のように階層的に記述できます (表記の乱用があります)。

y2 = I(z > 0)
(y1, z) ~ C(F_y1( |w), F_z( |w) | ファイ)
w, ファイ ~ 事前確率

ここで、w は y1 と z の周辺パラメーターのベクトル、F_y1 と F_z はそれぞれ y1 と z の周辺累積分布関数、phi はコピュラ パラメーターです。

これをスタンでどのようにモデル化できますか? コピュラによって生成された二変量尤度から y1 と z をサンプリングするカスタム確率関数を作成しました。どうすればよいかわからないのは、潜在変数 z を説明する (生成する) ことと、y2 と z の関係を指定する方法です。

stan でデータ拡張を使用したプロビット回帰を既に見ましたが、モデルにコピュラがあるため、これは役に立たないようです。

編集：上記のリンクが役に立たないことについて誤解している可能性があります。私は次のコードを書きました。(理論的に)正しいように見えるかどうかについてコメントをいただければ幸いです。

parameter を推測するために返されたというstanfitオブジェクトがあります。たとえば、を使用して分析できるようになりました。fitrstan::stan(...)thetathetarstan::summary(fit, pars="theta")

の 2 乗についての推論にもっと関心があることに後で気付きましたtheta。出力にパラメーターとしてtransformed parameters含めるブロックを STAN モデルに含める必要がありました。theta_squared

あたかもブロックで計算されたかのように、変換されたパラメータtheta_squared <- theta^2を既存のオブジェクトに追加することは可能ですか?stanfittransformed parameters

MC-StanをSparkで使いたいのですが、関連するページが Google で検索されていないようです。

このアプローチがSparkでも可能かどうか疑問に思うので、誰かが私に知らせてくれれば幸いです.

さらに、Spark で MCMC を使用するために広く使用されているアプローチは何なのかも気になります。Scala が広く使われていると聞きましたが、MC-Stan のような適切な MCMC ライブラリを備えた言語が必要です。

時変共変量を可能にする JAGS で生存モデルを構築しようとしています。たとえば、生存率がワイブル分布に従うと仮定します (ただし、ハザードを変化させたいので、指数関数は単純すぎます)。したがって、これは基本的にパッケージで実行できることのベイジアン バージョンでありflexsurv、パラメトリック モデルで時変共変量を可能にします。

したがって、「カウントプロセス」形式でデータを入力できるようにしたいと考えています。各被験者には複数の行があり、それぞれが共変量が一定のままである時間間隔に対応しています (この pdfまたはここで説明されています。これはまたはパッケージが許可する配合(start, stop]。survivalflexurv

残念ながら、JAGS で生存分析を実行する方法に関するすべての説明は、被験者ごとに 1 つの行を想定しているようです。

この単純なアプローチを採用して、カウント プロセス形式に拡張しようとしましたが、モデルは分布を正しく推定しません。

失敗した試み:

例を次に示します。まず、いくつかのデータを生成します。

次に、データを被験者ごとに 2 つの観測値に分割します。各被験者を時間 = 300 で分割しているだけです (時間 = 300 に達しなかった場合を除き、観察結果は 1 つしか得られません)。

これで、JAGS モデルをセットアップできます。

分割点がどこにあるかに応じて、モデルは基になる分布の非常に異なるパラメーターを推定します。データがカウント プロセス フォームに分割されていない場合にのみ、パラメーターが正しく取得されます。これは、この種の問題のデータをフォーマットする方法ではないようです。

仮定が欠けていて、問題が JAGS とはあまり関係がなく、問題の定式化に関連している場合は、提案を歓迎します。時変共変量がパラメトリック生存モデルで使用できないことに絶望しているかもしれません (また、一定のハザードを仮定し、実際には基になる分布を推定しない Cox モデルのようなモデルでのみ使用できます)。ただし、上で述べたように、flexsurvregR のパッケージは(start, stop]パラメトリック モデルの定式化に対応しています。

このようなモデルを別の言語 (JAGS の代わりに STAN など) で構築する方法を知っている人がいれば、それもありがたいです。

編集：

Chris Jackson が電子メールで役立つアドバイスを提供しています。

ここでは、JAGS の切り捨てのための T() コンストラクトが必要だと思います。基本的に、人が生きているが共変量が一定である各期間 (t[i]、t[i+1]) では、生存時間は期間の開始時に左打ち切られ、場合によっては右打ち切られます。終わり。だからあなたは次のようなものを書くでしょうy[i] ~ dweib(shape, scale[i])T(t[i], )

この提案を次のように実装しようとしました。

残念ながら、これではうまくいきません。古いコードでは、モデルはほとんどスケール パラメータを正しく取得していましたが、形状パラメータでは非常に悪い仕事をしていました。この新しいコードを使用すると、正しい形状パラメーターに非常に近くなりますが、スケール パラメーターは一貫して過大評価されます。過大評価の程度は、分割点が来るのが遅いことに相関していることに気付きました。分割点が早い ( cens_point = 50) 場合、過大評価はありません。遅い ( cens_point = 350) 場合は、たくさんあります。

この問題は、観測の「二重カウント」に関連している可能性があると考えました。t=300 で打ち切られた観測が見られた場合、同じ人物から、t=400 で打ち切られていない観測が見られた場合、この人物がは、ワイブル パラメータに関する推論に 2 つのデータ ポイントを提供していますが、実際には 1 つのポイントだけを提供する必要があります。したがって、私は各人にランダム効果を組み込んでみました。ただし、これは完全に失敗し、nuパラメーターの推定値が非常に大きくなりました (50 ～ 90 の範囲)。それがなぜなのかはわかりませんが、おそらくそれは別の投稿の質問です。私は問題が関連しているかどうかわからないので、そのモデルの JAGS コードを含む、この記事全体のコードをここで見つけることができます。

私はスタンの新人です。私は、地震のモデル化に使用されるモデルである時間 ETAS モデルに取り組んでいます。地震発生時刻 t[i] での震度は、

ここで、t は時間、p、c、mu は 3 つのパラメーターです。Rstanを使用しています。モデル用に次のスタンコードを書きました。

時間をベクトルとして指定していないことはわかっています。モデルセクションに尤度を書くのを手伝ってくれませんか? 強度を書く問題に直面しています。R で時間 t[i] での強度を記述するために使用した方法は、STAN でこれを行うための書き込み方法ではないと思います。

A small part (containing 20 times only) of the data is as follows: dat=list(0.0000,310.1907,948.4677,1007.2617,1029.7996,1065.7343,1199.8650, 1234.6809,1298.0234,1316.0350,1381.8400,1413.4311,1546.2059,1591.1326, 1669.5084, 1738.9363,1745.5503,1797.9980,1895.6705,1936.3146)

したがって、Python でプロットした次のデータがあります。

データは、私が扱っている微分方程式系の強制項の入力です。したがって、連続関数をこのデータに適合させる必要があるため、段階関数の不連続性に伴う安定性の問題に対処する必要はありません。残念ながら、それはかなり大きなデータセットです。

私は、微分方程式をコーディングしている言語であるスタンに翻訳するのが可能であり、それほど面倒ではない適合関数に終わろうとしています。手動でコーディングできる部分。

polyfitfrom から始めましたnumpyが、あまり良くありませんでした。UnivariateSplinefromを使用するscipyと、適切な適合が得られましたが、Stan に変換するのに扱いやすいものにはなりませんでした。したがって、他の言語に簡単に翻訳できる関数を返す、他のフィットへの提案を探していましたか? データの形状を見て、有用な周期的スプライン フィットはありますか?