如何使用pymc創建離散狀態馬爾可夫模型？

我想弄清楚如何正確使用pymc製作離散狀態馬爾可夫鏈模型。作爲一個例子（在nbviewer中查看），讓我們建立一個長度爲T = 10的鏈，其中馬爾可夫狀態是二進制的，初始狀態分佈是[0.2,0.8]，並且在狀態1中切換狀態的概率是0.01，而在狀態2中爲0.5如何使用pymc創建離散狀態馬爾可夫模型？

import numpy as np 
import pymc as pm 
T = 10 
prior0 = [0.2, 0.8] 
transMat = [[0.99, 0.01], [0.5, 0.5]]

爲了使模型，我使狀態變量的陣列和轉變概率的陣列依賴於狀態變量（使用pymc.Index功能）

states = np.empty(T, dtype=object) 
states[0] = pm.Categorical('state_0', prior0) 
transPs = np.empty(T, dtype=object) 
transPs[0] = pm.Index('trans_0', transMat, states[0]) 

for i in range(1, T): 
    states[i] = pm.Categorical('state_%i' % i, transPs[i-1]) 
    transPs[i] = pm.Index('trans_%i' %i, transMat, states[i])

取樣模型顯示th Ë狀態邊緣人都應該是什麼（用Matlab中與凱文·墨菲的BNT包構建模式相比）

model = pm.MCMC([states, transPs]) 
model.sample(10000, 5000) 
[np.mean(model.trace('state_%i' %i)[:]) for i in range(T)]

打印出：

我的問題

[-----------------100%-----------------] 10000 of 10000 complete in 7.5 sec 

[0.80020000000000002, 
0.39839999999999998, 
0.20319999999999999, 
0.1118, 
0.064199999999999993, 
0.044600000000000001, 
0.033000000000000002, 
0.026200000000000001, 
0.024199999999999999, 
0.023800000000000002]

的 - 這似乎不是最優雅用pymc建立馬爾可夫鏈的方法。有沒有一種更清晰的方式，不需要一系列確定性函數？

我的目標是爲更一般的動態貝葉斯網絡編寫一個基於pymc的包。

來源

2014-03-25 shpigi

據我所知，您必須將每個時間步的分佈編碼爲前一個時間步的確定性函數，因爲這就是它的原因 - 參數中沒有涉及任何隨機性，因爲您在問題中定義了它們建立。不過，我認爲你的問題可能更多的是尋找更直觀的方式來表示模型。另一種方法是將時間步轉換作爲前一時間步的函數直接編碼。

from pymc import Bernoulli, MCMC 

def generate_timesteps(N,p_init,p_trans): 
    timesteps=np.empty(N,dtype=object) 
    # A success denotes being in state 2, a failure being in state 1 
    timesteps[0]=Bernoulli('T0',p_init) 
    for i in xrange(1,N): 
     # probability of being in state 1 at time step `i` given time step `i-1` 
     p_i = p_trans[1]*timesteps[i-1]+p_trans[0]*(1-timesteps[i-1]) 
     timesteps[i] = Bernoulli('T%d'%i,p_i) 
    return timesteps 

timesteps = generate_timesteps(10,0.8,[0.001,0.5]) 
model = MCMC(timesteps) 
model.sample(10000) # no burn in necessary since we're sampling directly from the distribution 
[np.mean(model.trace(t).gettrace()) for t in timesteps]

來源

2015-03-25 13:17:46 Austen

在你想看看你的馬爾可夫鏈的長遠行爲的情況下，discreteMarkovChain包可能是有用的。這些例子通過定義一個轉換函數來顯示構建離散狀態馬爾可夫鏈的一些想法，該轉換函數告訴你每個狀態下一步可達狀態及其概率。您可以使用相同的功能來模擬過程。

來源

2016-04-09 18:35:29 Forzaa

如何使用pymc創建離散狀態馬爾可夫模型？

回答

相關問題