wwrechard / pydlm

A python library for Bayesian time series modeling
BSD 3-Clause "New" or "Revised" License
475 stars 98 forks source link

R and Python different results #49

Open SmitRohan opened 3 years ago

SmitRohan commented 3 years ago

I am trying to fit the following model (local level model) using pydlm package in Python-

Yt = θt + vt, vt ∼ N(0, Vt), θt = θt−1 + wt, wt ∼ N(0,Wt) …….(1)

I have simulated the above model in python with the following values for Vt and Wt-

sigma2_v = 0.5 ………..(2) sigma2_w = 0.25 ……….(3) x0 = 0 t = 500 t_burn = 100 t_tot = t + t_burn y_tot = np.zeros(t_tot) x_tot = np.zeros(t_tot)

v = np.random.normal(0, sqrt(sigma2_v), t_tot) w = np.random.normal(0, sqrt(sigma2_w), t_tot)

x_tot[0] = x0 + w[0] y_tot[0] = x_tot[0] + v[0]

for i in range(1,t_tot): x_tot[i] = x_tot[i - 1] + w[i] y_tot[i] = x_tot[i] + v[i]

I am using the following argument to fit model (1) in pydlm

myDLM1 = dlm(y_tot) trend1 = trend(degree=0, discount=0.99, name = 'trend') myDLM1 = myDLM1 + trend1 myDLM1.tune(maxit = 100) myDLM1.fit()

Below are my values for Vt and Wt- var_v = myDLM1.getVar(filterType='forwardFilter') var_v 0.996619 0.771301 1.183684 0.910351 0.774972 .. ... 1.003247 1.003131 1.001465 0.999822 0.999866

var_w = myDLM1.getVar(filterType='forwardFilter', name = 'trend') var_w

0.496571 0.285149 0.370662 0.257496 0.205303 .. ... 0.231519 0.231492 0.231107 0.230728 0.230738

Questions-

  1. If I fit model (1) with dlm package in R with same y series, I get the following estimates of Vt and Wt 0.55 0.23

Which are very close to the original sigma2_v and sigma2_w from which the y series has been simulated. I am unable to recover the similar values using pydlm.

Can you please help me out with this problem?

  1. Why does pydlm give Vt and Wt across all time points given they are assumed to be constant ?