How to get standard deviation using PyKalman?

Question

Having one dimensional measurement data I would like to know state standard deviations at each point using Kalman filter. My procedure is as follows:

from pykalman import KalmanFilter
import numpy as np

measurements = np.asarray([2, 1, 3, 6, 3, 2, 7, 3, 4, 4, 5, 1, 10, 3, 1, 5])
kf = KalmanFilter(transition_matrices=[1],
                  observation_matrices=[1],
                  initial_state_mean=measurements[0],
                  initial_state_covariance=1,
                  observation_covariance=1,
                  transition_covariance=0.01)
state_means, state_covariances = kf.filter(measurements)
state_std = np.sqrt(state_covariances[:,0])
print state_std

This leads to the following strange result:

[[ 0.70710678]
 [ 0.5811612 ]
 [ 0.50795838]
 [ 0.4597499 ]
 [ 0.42573145]
 [ 0.40067908]
 [ 0.38170166]
 [ 0.36704314]
 [ 0.35556214]
 [ 0.34647811]
 [ 0.33923608]
 [ 0.33342945]
 [ 0.32875331]
 [ 0.32497478]
 [ 0.32191347]
 [ 0.31942809]]

I would expect increased variance for the last data points. What am I doing wrong?

Answer 1

Since all the covariance matrix (measurement, transition) you provided are small (which means that you don't expect that much uncertainties in observation), the state covariance does not reflect your asymptotically increasing observation dispersion, and as a result the kalman fitler output is very smooth. But if you believe that there are more uncertainties in measurments, transition etc., I think you can provide higher covariance and as a result you will get KF output not very smooth (almost following the measurements), but the asymptotic increase will get reflected in the KF output covariance too, as shown below.

from pykalman import KalmanFilter
import numpy as np

measurements = np.asarray([2, 1, 3, 6, 3, 2, 7, 3, 4, 4, 5, 1, 10, 3, 1, 5])
kf = KalmanFilter(transition_matrices=[1],
                  observation_matrices=[1],
                  initial_state_mean=measurements[0],
                  initial_state_covariance=1,
                  observation_covariance=5,
                  transition_covariance=9) #0.01)
state_means, state_covariances = kf.filter(measurements)
state_std = np.sqrt(state_covariances[:,0])
print state_std
print state_means   
print state_covariances
import matplotlib.pyplot as plt
plt.plot(measurements, '-r', label='measurment')
plt.plot(state_means, '-g', label='kalman-filter output')
plt.legend(loc='upper left')
plt.show()

measurement_std = [np.std(measurements[:i]) for i in range(len(measurements))]
plt.plot(measurement_std, '-r', label='measurment std')
plt.plot(state_std, '-g', label='kalman-filter output std')
plt.legend(loc='upper left')
plt.show()