# Fit a linear State Space Model

I’would like to fit a state space model described as:

``````N: number of observations.
for t in 1......N :
y(t) = x(t) * z(t) +N(0,sigma_obs) : observed equation
x(t+1) = phi*x(t) + N(0,sigma_state) : state equation
x(0) ~ N(A0,P0) ; initial state prior
phi ~ N( mu = 0.89, sd = 0..89);
sigma_obs ~ InvGamma(a0,b0) ,;
sigma_state ~ InvGamma(c0,d0) ;
``````

the parameters to be estimated are `theta = {x(0), phi , sigma_obs, sigma_state }`
I’ve simulated data for a given values of the vector parameters theta, and coded the model in pymc3 but i get some issues.
Please if any body has any suggestions i’m great-full
the following code and simulate data :

``````def simulate_data(phi , nobs  , sigma_obs , sigma_state  , x0 ):
zt = np.random.normal(loc=0.0, scale=0.7, size=(nobs, 1))
v = np.random.normal(loc=0.0, scale=sigma_obs, size=(nobs, 1))
u = np.random.normal(loc=0.0, scale=sigma_state, size=(nobs , 1))
x = np.empty((nobs+1,1))
x[0,0] = x0
x[1:,] = u
D =  v
for i in range(0,nobs) :
x[i+1,0] += phi*x[i,0]
D[i,0]  += zt[i,0] * x[i+1,0]

data = {"phi": phi, "sigma_obs": sigma_obs, "sigma_state": sigma_state, 'zt': zt, 'y': D,
"x" : x}
return data
# true parameters
phi = 0.75 ;  nobs = 366  ; sigma_obs = 2.1 ; sigma_state = 1.5 ; x0 = 10
# simulation of Data
data = simulate_data(phi , nobs  , sigma_obs , sigma_state  , x0 )
y = data['y']; zt = data['zt'];  x =  data['x']
# Simulated data Viz
f, axarr  = plt.subplots( 3 ,1)
plt.subplots_adjust(wspace = 0.35 , hspace = 0.35 )

axarr[0].plot(x)
axarr[0].set_title('state  Evolution over time')
axarr[0].set_ylabel("latent state")

axarr[1].plot(y)
axarr[1].set_title('Demand Variation Evolution over time')
axarr[1].set_ylabel("Demand Variation")

axarr[2].plot(zt)
axarr[2].set_title('price  Evolution over time')
axarr[2].set_ylabel("price Variation")

# Model Definition

y_shared  = theano.shared(y, 'y')
zt_shared = theano.shared(zt, 'zt')
nobs  = y.shape[0]
ssm_model = Model()
with ssm_model:
# Specify prior distribution
a0 = 20
b0 = 20
c0 = 40
d0 = 20

sigma_obs = InverseGamma("sigma_obs", alpha=a0, beta=b0)
sigma_state = InverseGamma("sigma_state", alpha=c0, beta=d0)
phi = Normal('phi',mu = 0, tau = 0.0005)

x_0 = Normal('x0', mu=0, sd =200)  # Initioal Prior State at time t = 0
x = np.empty(nobs+1)
d = np.empty(nobs)
for i in range(0, nobs) :
#     # temp[i] = Deterministic("mu{}".format(i), x[i+1]*phi)
x[i+1] = Normal("x{}".format(i+1) , mu = x[i]*phi, sd = sigma_state)
d[i] = Normal("d{}".format(i), mu = zt_shared[i]*x[i+1], sd = sigma_obs, observed=y_shared[i])
``````

Issue : I got this error will compiling the ssm_model model

`````` File "C:\ProgramData\Anaconda3\lib\site-packages\IPython\core\interactiveshell.py", line 2881, in run_code
exec(code_obj, self.user_global_ns, self.user_ns)
File "<ipython-input-123-5f2b696cbea7>", line 22, in <module>
x[i+1] = Normal("x{}".format(i+1) , mu = x[i]*phi, sd = sigma_state)
ValueError: setting an array element with a sequence.
``````

Hi @ayoubsbitti, the way you are coding the state space model is very inefficient. You should use the time series distribution in PyMC3. For a starter, please have a look at the following 2 examples:
http://docs.pymc.io/en/latest/notebooks/stochastic_volatility.html
http://docs.pymc.io/en/latest/notebooks/survival_analysis.html