# How to Derive SISO Kalman Filter Update Equation Using only Probability Density Functions

I'm trying to prove to myself that a single state/single measurement kalman update can be derived using bayes theorem (as proof of concept for a more complicated task) only using the PDF. I am able to arrive at the correct value for the state distribution mean after the update, but cant figure out how the variance is determined. Any help would be appreciated.

To be clear, I fully understand how to derive the KF in the traditional matrix math way, but I would like to do it interms of only the pdfs and their parameters.

State Model (prior):

$$x_{k|k-1} \thicksim \mathcal{N}(\mu_{k|k-1},\sigma^2_{k|k-1}) = p(x_{k|k-1}|y_{1:k-1})$$

Measurement Model:

$$y_k = h*x + \mathcal{N}(0,\sigma^2_R)$$

$$y_k \thicksim \mathcal{N}(hx,\sigma^2_R) = p(y_k|x_k)$$

Bayes Theorem/Posterior:

$$p(x_{k|k}|y_k) \propto p(y_k|x_k)p(x_{k|k-1}|y_{1:k-1})$$

$$J = \log p(y_k|x_k) + \log p(x_{k|k-1}|y_{1:k-1})$$

The Optimization Problem:

I THINK I want to maximize the expected value of J with respect to the the new posterior. That is:

$$\underset{\mu_{k|k},\space \space \sigma^2_{k|k}}{\operatorname{argmax}} \int J(x,y,\mu_{k|k-1},\sigma^2_{k|k-1})\mathcal{N}(x,\mu_{k|k},\sigma^2_{k|k})dx$$

The above integral evalutes to:

$$\frac{\sigma^{2}_{R} \left(- \mu_{k}^{2} + 2 \mu_{k} \mu_{k-1} - \mu_{k-1}^{2} - \sigma^{2}_{k}\right) + \sigma^{2}_{k-1} \left(- h^{2} \mu_{k}^{2} - h^{2} \sigma^{2}_{k} + 2 h \mu_{k} y - y^{2}\right) - \log{\left (\left(4 \pi^{2} \sigma^{2}_{R} \sigma^{2}_{k-1}\right)^{\sigma^{2}_{R} \sigma^{2}_{k-1}} \right )}}{2 \sigma^{2}_{R} \sigma^{2}_{k-1}}$$

Results:

Mean:

$$\mu_{k|k} = \frac{h \sigma^{2}_{k|k-1} y + \mu_{k|k-1} \sigma^{2}_{R}}{h^{2} \sigma^{2}_{k|k-1} + \sigma^{2}_{R}}$$

Variance: The above optimization problem cant be optimized with respect to the posterior variance because the derivative with respect to it is a constant, which probably means i'm solving the wrong problem.

Any insights would be greatly appreciated!

Some python:

import sympy
sympy.init_printing()

x      = sympy.Symbol('x',real=True)
y      = sympy.Symbol('y',real=True)
h      = sympy.Symbol('h',real=True);
mu_km1 = sympy.Symbol('mu_k-1',real=True)
mu_k   = sympy.Symbol('mu_k',real=True)
varR   = sympy.Symbol('sigma^2_R',positive=True,real=True);
varXk  = sympy.Symbol('sigma^2_k',positive=True,real=True);
varXkm1  = sympy.Symbol('sigma^2_k-1',positive=True,real=True);

def pdf(x,mu,vari):
den = sympy.sqrt( 2*sympy.pi*vari  )
num = sympy.exp( -(1/(2*vari)) * (x-mu)**2  )
return num/den;

lik   = pdf(y,h*x,varR);
prior = pdf(x,mu_km1,varXkm1);
J     = ( (sympy.log( lik*prior  ))  ).expand();
JPost = (J*pdf(x,mu_k,varXk)).simplify();
integ = sympy.integrate(JPost,(x,-sympy.oo,sympy.oo)).simplify()
new_mu = sympy.solve(integ.diff(mu_k),mu_k)[0]


$$\int_{-\infty}^{\infty} \frac{\sqrt{2} \left(\sigma^{2}_{R} \sigma^{2}_{k} \left(\mu_{k-1} - x\right)^{2} - \sigma^{2}_{R} \sigma^{2}_{k-1} \left(\mu_{k} - x\right)^{2} + \sigma^{2}_{k} \sigma^{2}_{k-1} \left(h x - y\right)^{2} + \log{\left (\left(\frac{2 \pi \sigma^{2}_{R} \sigma^{2}_{k-1}}{\sigma^{2}_{k}}\right)^{\sigma^{2}_{R} \sigma^{2}_{k} \sigma^{2}_{k-1}} \right )}\right) e^{- \frac{\left(\mu_{k} - x\right)^{2}}{2 \sigma^{2}_{k}}}}{4 \sqrt{\pi} \sigma^{2}_{R} \left(\sigma^{2}_{k}\right)^{\frac{3}{2}} \sigma^{2}_{k-1}}\, dx$$