General result for the N-point correlation of the Poisson process (and its derivative)?

Question

I am specifically interested in computing: $$\mathbb{E}[S_p S_q S_s S_t]$$ where $S_t=\frac{dN_t}{dt}$ and $N_t$ is a Poisson process (so $S_t$ is a "Poisson pulse train"): $$\mathbb{P}(N_t=n)=\frac{(\lambda t)^n}{n!}e^{-\lambda t} \; , \; t>0$$

Here is my attempt:

Assuming $t>s>p>q$, I tried to compute it from the process $N_t$ rewriting it with the independent random variables $N_t-N_s, N_s-N_p, N_p-N_q, N_p$, that we can rename $N_j$, $j=1,2,3,4$: $$\mathbb{E}[N_p N_q N_s N_t]=\mathbb{E}[\prod_{k=1}^4\sum_{j=1}^k N_j]$$

This unfolds a combinatorics which somewhat looks like the one we obtain for normally distributed variables using the Isserlis / Wick theorem.

Is there an analogous theorem about such correlation function which applies to the Poisson process?

This looks simpler for the pulse train since the only term which does not vanish when we take the derivative is: $$\mathbb{E}[S_p S_q S_s S_t]=\frac{d^4}{dpdqdsdt}\mathbb{E}[N_p N_q N_s N_t]=\frac{d^4}{dpdqdsdt}\big(\mathbb{E}[N_1]\mathbb{E}[N_2]\mathbb{E}[N_3]\mathbb{E}[N_4]\Big)=\lambda^4$$ However, this is without considering coincidences among $t,s,p,q$. From here, I see (p.10): $$\mathbb{E}[S_s S_t]=\lambda^2+\lambda \delta(s-t)$$

I deduce that we should have something like $$\mathbb{E}[S_p S_q S_s S_t]=\lambda^4+\lambda^3\sum_{i\neq j}\delta(t_i-t_j)+\lambda^2\sum_{i\neq j,\; k\neq l }\delta(t_i-t_j)\delta(t_k-t_l)+\lambda\sum_{i\neq j,\; k\neq l,\; m\neq n }\delta(t_i-t_j)\delta(t_k-t_l)\delta(t_m-t_n)$$

Is this correct (except the sloppily specified summation indices)? Do you know any related general result?

Answer 1

This question (and its generalization) is conventiently addressed by considering the moment-generating functional. Let $N_t$ be a counting process with (possibly stochastic and time-dependent) intensity $\lambda_t$ (you recover your case by setting $\lambda_t=\lambda$, where $\lambda$ is a constant).

Suppose that under $\mathbb{P}$, the compensated counting process $N_t-\int_0^t\lambda_s ds$ is independent of $\mathcal{F}^{\lambda}_t$, the natural filtration of $\lambda_t$. Then we have

$$Z[a]\doteq \mathbb{E}\Bigg[\exp\Big[\int_0^Ta(t)dN_t\Big]\Bigg]=\mathbb{E}\Bigg[\exp\Big[\int_0^T\big(e^{a(t)}-1\big)\lambda_t dt\Big]\Bigg], \quad (\star)$$

for e.g. continuous $a$ for which the RHS exists. Indeed, under the above assumption there is an equivalent measure $\mathbb{P}^{\ast}$ under which $N_t$ has unit intensity and the law of $(\lambda_t)_{t\geq 0}$ is the same as under the original measure $\mathbb{P}$. By using a Girsanov-type change-of-measure formula (see e.g. Theorems 26.8 and 26.9 in Kallenberg, O. (2006). Foundations of Modern Probability.) for the measure $\mathbb{P}^{\ast}$, we obtain

$$ \begin{split} \mathbb{E}\Bigg[\exp\Big[\int_0^Ta(t)dN_t\Big]\Bigg]&=\mathbb{E}^{\ast}\Bigg[\frac{d\mathbb{P}}{d\mathbb{P}^{\ast}}\exp\Big[\int_0^Ta(t)dN_t\Big]\Bigg] \\ &=\mathbb{E}^{\ast}\Bigg[\exp\Big[\int_0^T\log \lambda_tdN_t+\int_0^T\big(1-\lambda_t\big)dt\Big]\exp\Big[\int_0^Ta(t)dN_t\Big]\Bigg]\\ &=\mathbb{E}^{\ast}\Bigg[\exp\Big[\int_0^T\log \big( e^{a(t)}\lambda_t\big)dN_t+\int_0^T\big(1-e^{a(t)}\lambda_t+(e^{a(t)}-1)\lambda_t\big)dt\Big]\Bigg] \\ &=\mathbb{E}^{\ast}\Bigg[\exp\Big[\int_0^T\log \big( e^{a(t)}\lambda_t\big)dN_t+\int_0^T\big(1-e^{a(t)}\lambda_t\big)dt\Big] \\ &\qquad\quad\times\exp\Big[\int_0^T(e^{a(t)}-1)\lambda_tdt\Big] \Bigg]\\ &=\mathbb{E}'\Bigg[\exp\Big[\int_0^T\big(e^{a(t)}-1\big)\lambda_t dt\Big]\Bigg]. \end{split}$$ In the last line, the measure $\mathbb{P}'$ is such that $N_t$ has rate $a(t)\lambda_t$, but the law of $(\lambda_t)_{t\geq 0}$ is the same as under both $\mathbb{P}$ and $\mathbb{P}^{\ast}$. Since the integrand only depends on $(\lambda_t)_{t\geq 0}$, $(\star)$ follows.

---

From $(\star)$, we may extract arbitrary correlation functions using functional derivatives. We have

$$\frac{\delta Z[a]}{\delta a(t)}=\mathbb{E}\Bigg[e^{a(t)}\lambda_t\exp\Big[\int_0^T\big(e^{a(s)}-1\big)\lambda_s ds\Big]\Bigg],$$

on which we can iteratively apply additional functional derivatives. E.g.

$$\mathbb{E}[S_t]=\frac{\delta Z[a]}{\delta a(t)}\Bigg|_{a=0}=\mathbb{E}[\lambda_t],$$

$$\mathbb{E}[S_tS_{t'}]=\frac{\delta^2 Z[a]}{\delta a(t)\delta a(t')}\Bigg|_{a=0}=\delta(t-t')\mathbb{E}[\lambda_t]+\mathbb{E}[\lambda_t\lambda_{t'}].$$

By specializing to the case of fixed deterministic rate, you can check your guesses for higher-order moments of a Poisson pulse train.

EDIT: I had to add an assumption for $(\star)$ to hold.