I'm trying to get my head around online (recursive) maximum-likelihood parameter estimation in the language of stochastic processes and in the context of stochastic filtering, i.e. where we have a partially observed Ito diffusion process \begin{align} dX_t&=f(t,\theta_0,X_t)dt+g(t,\theta_0,X_t)dW_t, \\ dY_t&=h(t,\theta_0,X_t)dt+dV_t. \end{align} Here, everything is scalar, $X_t$ is the (unobserved) signal, $Y_t$ is the observation, $W_t$ and $V_t$ are independent standard Wiener processes, and $\theta_0$ is a parameter. Suppose we don't know $\theta_0$ and we want to estimate it recursively as observations $Y_t$ come in.

Moura & Mitter (1986) give a formula for the log likelihood function for parameter estimation based on only observations of $Y_t$, namely \begin{equation} L_t(\theta)=\frac{1}{r}\int_0^t\left ( \hat h_s dY_s-\frac{1}{2} \hat h^2_s ds \right), \quad \hat h_t=E \left[ h(t,\theta,X_t)|\mathcal{F}^Y_t \right]. \end{equation} First of all, it is not entirely clear why the above function would be the correct thing to look at, and the paper mentioned above does not really help me with a lot of detail. What I can guess so far is that they use a sequence of changes of measures, ending up with one (let's call it $\mathcal{Q}_{\theta}$) under which the innovations process $n_t$, defined as \begin{equation} dn_t=dY_t-\hat h_t dt, \end{equation} is a Wiener process. By Girsanov's theorem, the above log likelihood function is the log of the Radon-Nikodym (RN) derivative $d\mathcal{Q}_{\theta}/d\tilde{\mathcal{P}}_{\theta}|_{\mathcal{F}_t}$ of that measure with respect to another measure $\tilde{\mathcal{P}}_{\theta}$ under which $Y_t$ is a Wiener process (is this correct?). However, what I would want to maximize looks something like \begin{equation} L'_t(\theta)=\log \frac{d\mathcal{P}^Y_{\theta}}{d\mathcal{P}^Y_0}, \end{equation} i.e. the RN derivative of the measure of the process $Y_t$, with $X_t$ 'integrated' or 'marginalized' out, with respect to a parameter-independent reference measure $\mathcal{P}^Y_0$. So my first question is: why should $L_t$ be the log likelihood function instead of $L'_t$, or alternatively, why is it equivalent to maximize $L'_t$?

For my second question, let's pretend that we understand why $L_t$ is the log likelihood function, and let's maximize it. For simplicity, I assume that everything is linear, \begin{align} dX_t&=-X_tdt+dW_t, \\ dY_t&=\theta X_tdt+dV_t, \end{align} such that \begin{equation} L_t(\theta)=\frac{1}{r}\int_0^t\left ( \theta\mu_s dY_s-\frac{1}{2} \theta^2(\mu_s)^2 ds \right), \end{equation} where $\mu_t$ is the mean of the Kalman-Bucy filter \begin{align} d\mu_t&=-\mu_tdt+\theta \sigma^2_t(dY_t-\theta\mu_tdt), \\ \frac{d\sigma^2_t}{dt}&=1-2\sigma^2_t-\theta^2\sigma^4_t, \end{align} which is integrated using the parameter $\theta$, acquiring an implicit parameter dependence. To maximize $L_t$ for fixed $t$, we could use a steepest ascent algorithm, i.e. \begin{equation} \begin{split} \frac{d}{d\tau}\hat \theta_{\text{ML,offline}}(\tau)&=\eta \partial_{\theta}L_t(\theta)\Big|_{\theta=\hat \theta_{\text{ML,offline}}(\tau)}\\ &=\frac{\eta}{r}\int_0^t\left(dY_s-\theta\mu_sds\right)\left(\mu_s+\theta\nu_s\right)\Big|_{\theta=\hat \theta_{\text{ML,offline}}(\tau)}. \end{split} \end{equation} Here, $\tau$ is the parameter of the offline gradient ascent, i.e. while the process $(X,Y)$ is stopped (observations from 0 to $t$ have been collected and stored), the gradient ascent is initialized at $\tau=0$ with an initial guess $\theta_{\text{ML,offline}}(0)$, and then updated according to the differential equation above until it converges. We also introduced a learning rate $\eta>0$ and the derivative of the mean of the Kalman-Bucy filter, $\nu_t=\partial_{\theta}\mu_t$, for which we find \begin{align} d\nu_t&=-\nu_tdt+(\sigma^2_t+\theta\rho^2_t)(dY_t-\theta\mu_tdt)-\theta\sigma^2_t(\mu_t+\theta\nu_t)dt, \quad &\nu_0=0,\\ \frac{d\rho^2_t}{dt}&=-2\rho^2_t-2\theta^2\sigma^2_t\rho^2_t-2\theta\sigma^4_t\quad &\rho^2_0=0. \end{align} Finally, we could turn the above gradient-based update equation into a recursive update rule by getting rid of the integral: \begin{equation} d\hat\theta_t=\frac{\eta}{r}\left(dY_t-\hat\theta_t\mu_tdt\right)\left(\mu_t+\hat\theta_t\nu_t\right), \end{equation} a rule which works well enough in practice (simulations). The problem is reconciling this with another derivation one could come up with. Let's discretize the Ito SDEs above: \begin{align} X_i&=X_{i-1}-\delta t X_{i-1}+\sqrt{\delta t}\xi_i, \\ Y_i&=Y_{i-1}+\theta\delta t X_{i-1}+\sqrt{\delta t}\eta_i, \end{align} where $\xi_i,\eta_i$ are sampled i.i.d from $\mathcal{N}(0,1)$. This is the Euler-Maruyama discretization. We can write down a log likelihood function in discrete time: \begin{align} L^D_n(\theta)&=\log P(\mathbf{Y}|\theta), \\ P(\mathbf{Y}|\theta)&=\int d\mathbf{X}\, P(\mathbf{X},\mathbf{Y}|\theta), \\ P(\mathbf{X},\mathbf{Y}|\theta)&=P(X_1)P(Y_1)\prod_{i=2}^n P(X_i|X_{i-1}) P(Y_i|X_{i-1},Y_{i-1},\theta). \end{align} Taking the derivative, we obtain \begin{equation} \begin{split} \partial_{\theta}L^D_n(\theta)&=\frac{1}{P(\mathbf{Y}|\theta)}\int d\mathbf{X}\,\partial_{\theta} P(\mathbf{X},\mathbf{Y}|\theta) \\ &=\int d\mathbf{X}\, P(\mathbf{X}|\mathbf{Y},\theta)\partial_{\theta}\log P(\mathbf{X},\mathbf{Y}|\theta). \end{split} \end{equation} Now, using the structure of the joint distribution in the integrand, and the explicit form of the Gaussian distribution, one ends up with \begin{equation} \partial_{\theta}\log P(\mathbf{X},\mathbf{Y}|\theta) = \sum_{i=2}^n \left[X_{i-1}(Y_i-Y_{i-1})-\theta\delta t X^2_{i-1} \right], \end{equation} which, after taking the continuum limit, expectations over $P(\mathbf{X}|\mathbf{Y},\theta)$ and getting rid of the integral, seems to lead to an update rule \begin{equation} d\hat\theta_t\propto \left[\mu_t dY_t-\theta(\mu^2_t+\sigma^2_t)dt \right]. \end{equation} So the second question (sorry for being so long-winded) is: Why are the two rules different? Note that they are fairly similar, with the only differences being that the first rule has expectations under the square, and the second one expectation of the square, and the second one does not contain filter derivatives. Is there a flaw in one of these derivations? Or is it to be expected that we get different rules because the problem is not well-posed?

  • $\begingroup$ Concerning your first question, $L_t(\theta)$ is the R-N derivative for the law of $Y$ with respect to the reference measure where $Y$ is simply a Wiener process (independent of course of the parameter $\theta$). (In the old filtering literature, this is referred to as the "Kallianpur-Striebel formula"). Since the reference measure does not depend on $\theta$, this seems to me as a good measure of likelihood. $\endgroup$ – ofer zeitouni May 2 '15 at 12:38
  • $\begingroup$ As for your second question, I am not sure about your notation, but don't you run into the trouble of non-adapted integration? $\endgroup$ – ofer zeitouni May 2 '15 at 12:38
  • $\begingroup$ Thanks for the hint, it just crossed my mind that the innovation process $n_t$ defined above is a Wiener process under the law of $Y$. I haven't seen it called the Kallianpur-Striebel formula, I thought that was the formula for expressing expectations conditioned on observations of $Y_t$ as quotient of expectations under the law where $Y_t$ is a Wiener process (for which one can then derive the DMZ equation)? $\endgroup$ – S.Surace May 2 '15 at 20:12
  • $\begingroup$ Could you elaborate on your comment re: non-adapted integration? What is unclear about the notation? I think up to $\partial \log P...$ it should be fine, but the step following that (continuum limits etc.) is not clear to me either, and I don't know how to make it rigorous. The main issue seems to be somehow going from the conditional distribution $P(\mathbf{X}|\mathbf{Y},\theta)$ to something resembling a filter (i.e. running forward in time). $\endgroup$ – S.Surace May 2 '15 at 20:15
  • $\begingroup$ I am not sure what is $\tau$, but in any case, in the equation defining $\hat \theta_{ML,offline}$, doesn't $\hat \theta$ depend on the whole path and therefore the stochastic integral is non-adapted? $\endgroup$ – ofer zeitouni May 2 '15 at 21:46

The formulation of a likelihood function for a continuous-time stochastic process requires the choice of a reference measure. In the context of parameter estimation, the reference measure has to be parameter-independent. In the case of partially observed diffusion processes with additive observation noise, the innovation process is a Brownian motion adapted to the observation filtration, and consequently there is a reference measure under which $Y_t$ is a Brownian motion.

The big miracle is that the Radon-Nikodym derivative gives us an expression of the likelihood function in terms of the filter, which makes it perfect for online parameter estimation.

The full story, including numerical examples, is written up in the following article: Online Maximum Likelihood Estimation of the Parameters of Partially Observed Diffusion Processes


Your Answer

By clicking "Post Your Answer", you acknowledge that you have read our updated terms of service, privacy policy and cookie policy, and that your continued use of the website is subject to these policies.

Not the answer you're looking for? Browse other questions tagged or ask your own question.