Bayesian Estimation – How to Perform Bayesian Parameter Estimation of a Poisson Process with Change/No-Change Observations at Irregular Intervals

bayesianconjugate-priorestimationpoisson process

Consider a Poisson process with unknown parameter $\lambda$.
We perform a sequence of $n$ observations at intervals $\overline{t}=t_1,\,t_2,\,\dots,\,t_n$. Each observation is a binary variable $x_i$ equal to zero if no changes occurred during interval $t_i$ or equal to one if one or more changes occurred:
$\forall i,$
$\Pr[X_i = 0 | \lambda,\,t_i] = \exp(-\lambda t_i)$
$\Pr[X_i = 1 | \lambda,\,t_i] = 1 – \exp(-\lambda t_i)$

The likelihood function for all observations is:
$P[\overline{x}| \lambda,\overline{t}] = (\prod_{x_i=0}\exp(-\lambda t_i)) \cdot (\prod_{x_i=1}(1-\exp(-\lambda t_i)))$

Assuming that the observation intervals are chosen independently from anything else, I would like to estimate the parameter $\lambda$ under some reasonable prior
$\lambda^\star = argmax_\lambda P[\lambda | \overline{x},\,\overline{t}]$
$P[\lambda | \overline{x},\,\overline{t}] \propto P[\overline{x}| \lambda,\overline{t}] \cdot P[\lambda]$

If possible, I would like to use a conjugate prior in order to perform recursive estimation, but I'm not sure whether one exists.

My questions:
1. Does a conjugate prior for that likelihood exist?
2. If a conjugate prior does not exist, what estimator can I use? I'm interested in an estimator that can be updated incrementally for each observation, without keeping track of all the observation history.

Best Answer

Notation

First I'll re-express the log-likelihood:

$$\log L(\lambda)=\sum_i \bigg[ x_i \log(1-e^{- \lambda t_i}) - \lambda (1-x_i)t_i\bigg].$$

Let

  • $\tau_0 = \sum_i (1-x_i)t_i$,
  • $\tau_1 = \sum_i x_i t_i$,
  • $\tau_2 = \sum_i x_i t_i^2$,
  • $n_1 = \sum_i x_i$.

Approximation

According to Wolfram Alpha,

$$ \log(1-e^{-x}) \approx \log(x) - \frac{x}{2} + \frac{x^2}{24}.$$

As $x \searrow 0$, $\log(1-e^{-x}) \rightarrow -\infty$, and as $x \rightarrow \infty$, $ \log(1-e^{-x}) \nearrow 0$. On the other hand, the approximation contains an upward-facing parabola, so it will diverge to $+\infty$ as $x \rightarrow \infty$. So when we find the maximum of the log-likelihood, only the smaller root of the parabola will be in the region where the approximation is accurate.

Estimation

Up to constant terms, substituting in the above approximation results in the following approximate log-likelihood:

$$\log \tilde{L}(\lambda) = n_1 \log(\lambda) - \lambda \bigg(\tau_0 + \frac{\tau_1}{2}\bigg) + \lambda^2\bigg(\frac{\tau_2}{24}\bigg).$$

The approximate MLE, $\hat{\lambda}$, satisfies

$$ \begin{array}{rcl}\frac{{\rm {d}}}{{\rm {d}\lambda}}\log\tilde{L}(\lambda)\bigg|_{\lambda=\hat{\lambda}} & = & 0\\n_{1}-\hat{\lambda}\bigg(\tau_{0}+\frac{\tau_{1}}{2}\bigg)+\hat{\lambda}^{2}\bigg(\frac{\tau_{2}}{12}\bigg) & = & 0\\ \Rightarrow \hat{\lambda} & = & \frac{6\tau_{0}}{\tau_{2}}+\frac{3\tau_{1}}{\tau_{2}}-\frac{6}{\tau_{2}}\sqrt{\bigg(\tau_{0}+\frac{\tau_{1}}{2}\bigg)^{2}-\frac{n_{1}\tau_{2}}{3}}.\end{array} $$

Prior distribution

The gamma distribution isn't exactly conjugate but will certainly be easy to combine with the above approximate log-likelihood.

Related Question