Flexible Linear Kalman Filter

Keywords
kalman, time, series, signal, filter
BSD-3-Clause
Install
``` pip install linkalman==0.11.5 ```

### Documentation

`linkalman` is a python package that solves linear structural time series models with Gaussian noises. Compared with some other popular Kalman filter packages written in python, linkalman has a combination of several advantages:

• Account for partially and fully incomplete measurements
• Flexible and convenient model structure
• Robust and efficient implementation
• Proper implementation for unknown priors
• Built-in numerical and EM algorithm
• Open-source with a comprehensive user manual
• Modular design with intuitive model specification

### Installation

`linkalman` requires the following packages to run:

• numpy
• pandas
• networkx
• scipy

To install `linkalman`, simply use the standard `pip` command:

`\$ pip install linkalman`

### Example

Here I will provide a simple example using `linkalman`. See here for more examples, and user's manual for technical details.

```import pandas as pd
import numpy as np
from scipy.optimize import minimize
from linkalman.models import BaseConstantModel as BCM
import matplotlib.pyplot as plt

# Get data
df['x'] = 1
df.set_index('Date', inplace=True)```

First we define the system dynamics of a Bayesian Structural Time Series (BSTS) model. Here I define a Stochastic linear trend model to extract the trend information from the time series (referring to the example section of user's manual for details)

```def my_f(theta):
sig1 = np.exp(theta[0])
sig2 = np.exp(theta[1])
sig3 = np.exp(theta[2])

F = np.array([[1, 1], [0, 1]])
Q = np.array([[sig1, 0], [0, sig2]])
R = np.array([[sig3]])
H = np.array([[1, 0]])
# Collect system matrices
M = {'F': F, 'Q': Q, 'H': H, 'R': R}

return M ```

Next we define a solver or optimizer, you can choose any solver you prefer. Here I just use `scipy.optimize.minimize`.

```def my_solver(param, obj_func, verbose=False, **kwargs):
obj_ = lambda x: -obj_func(x)
res = minimize(obj_, param, **kwargs)
theta_opt = np.array(res.x)
fval_opt = res.fun
return theta_opt, fval_opt```

Now we can fit the data. First we initialize the model and feed the system dynamics (`my_f`) and solver (`my_solver`). You may also pass the keyworded arguments to for `my_f` and `my_solver`.

```model = BCM()
model.set_f(my_f)
options={'xatol': 1e-8, 'disp': True, 'maxiter': 10000})
theta_init = np.random.rand(3)
model.fit(df, theta_init, y_col=['Births'], x_col=['x'],
method='LLY')
df_LLY = model.predict(df)```

That is it! If you want to do additional work, you can do the following to plot a confidence interval around your predictions.

```df_LLY['kf_ub'] = df_LLY.Births_filtered + 1.96 * np.sqrt(df_LLY.Births_fvar)
df_LLY['kf_lb'] = df_LLY.Births_filtered - 1.96 * np.sqrt(df_LLY.Births_fvar)
df_LLY = df_LLY[df_LLY.index > '1959-01-01']
df_LLY.index = pd.to_datetime(df_LLY.index)

# Define plot function
def simple_plot(df, col_est, col_actual, col_ub, col_lb, label_est,
label_actual, title, figsize=(12, 8)):
ax = plt.figure(figsize=figsize)
plt.plot(df.index, df[col_est], 'r', label=label_est)
plt.scatter(df_LLY.index, df[col_actual], s=20, c='b',
marker='o', label=label_actual)
plt.fill_between(df.index, df[col_ub], df[col_lb], color='g', alpha=0.2)
ax.legend(loc='right', fontsize=9)
plt.title(title, fontsize=22)
plt.show()
simple_plot(df_LLY, 'Births_filtered', 'Births', 'kf_ub', 'kf_lb',
'Prediction', 'Births', 'Filtered Births Data')```