Pith. sign in

REVIEW 2 major objections 3 minor 25 references

Markovian multivariate Hawkes population processes: Efficient evaluation of moments

T0 review · 2 major / 3 minor · reviewed 2026-08-07 · deepseek-v4-flash

Challenge this review Re-run · record.json Download PDF Read on arXiv ↗

Pith's one-line read For Markovian multivariate Hawkes processes with exponential decay, the joint transform of population and intensity has a closed form, and differentiating it yields exact transient and stationary cross-moments of any order.

desk verdict The joint transform is a solid new result, but the moment recursions factorize joint jump-size moments incorrectly and only hold for independent marks. read the letter →

arxiv 2506.08775 v2 pith:JKEQSJOV submitted 2025-06-10 math.PR

R.S. Karim , R.J.A. Laeven , M , M. Mandjes This is my paper
classification math.PR MSC 60G5560J2560K2560E10
keywords HawkesprocessesmutualexcitationMarkovtransformanalysismomentcomputationsinfinite-serverqueuespopulationnearlyunstableregime
verification ladder T0 review T1 audit T2 compute T3 formal

The pith

A machine-rendered reading of the paper's core claim, the machinery that carries it, and where it could break.

The reading

This paper proves that for Markovian multivariate Hawkes processes—point processes whose intensities decay exponentially and jump at every event—the joint transform of the population size and the current intensity is available in closed form as a product of exponentials whose arguments solve a small system of ordinary differential equations. Because the transform is bijective to the joint distribution, differentiating it yields exact analytic expressions for transient and stationary cross-moments of any order, together with auto- and cross-covariance functions, rather than simulation-based estimates. The paper organizes the moment relations into a nested sequence of block lower triangular matrices, which makes the computation fast and, for transient moments, insensitive to the time horizon. It also shows that in a symmetric parameterization the stationary intensity, scaled by the distance to instability, converges to a multivariate Gamma law. If these results are correct, moment-based inference and comparative statics for contagious-event models become practical in settings where finite differences and Monte Carlo simulation were the only options.

What carries the argument

The load-bearing object is the joint transform $\zeta_{t_0}(t,s,z)$, a combined z-transform in $Q(t)$ and Laplace transform in $\lambda(t)$ that is bijective to the joint law of the population and intensity processes. The derivation writes the Markov generator of $(Q,\lambda)$ as a first-order PDE, (58), and solves it by the method of characteristics; the characteristic equations are precisely the ODE system (9). The moment machinery is repeated differentiation of this transform, yielding the moment recursions (21) and (23). The computational shortcut is the nested sequence of block lower triangular matrices $A_n$ and $F_n$ built from the blocks $M^{(k,n-k)}$, whose invertibility turns the moment ODEs into closed-form matrix-exponential expressions (40) and (44).

What would settle it

Set the bivariate parameters as in Section 7.2, solve the ODE system (9) for a nontrivial pair $(s,z)$, and compare the value of $\zeta_0(t,s,z)$ from Eq. (8) with a Monte Carlo estimate of $\mathbb{E}[z_1^{Q_1(t)}z_2^{Q_2(t)}e^{-s_1\lambda_1(t)-s_2\lambda_2(t)}]$ using enough replications that the Monte Carlo standard error is negligible; a persistent gap beyond that error would refute the closed-form transform, and agreement would confirm it.

Watch

Extended reading notes

Core claim

The central claim is Theorem 1: for a $d$-component Hawkes process with exponential decay rates $\alpha_i$ and independent jump marks $B_{ij}$, and for the induced infinite-server population $Q(t)$, the conditional joint transform $\zeta_{t_0}(t,s,z)=\mathbb{E}[\prod_{i=1}^d z_i^{Q_i(t)}e^{-s_i\lambda_i(t)}]$ is given in closed form by (8), where the functions $b_{z_j}$ and $\tilde{s}_j$ solve the explicit ODE system (9). Repeated differentiation of this transform produces the linear ODE recursion (21) for transient reduced moments and the algebraic recursion (23) for stationary moments, both of total order at most the target order, so all cross-moments can be solved recursively. In the bivariate setting the recursions are reorganized as nested block lower triangular matrices, and Proposition 4 gives transient moments up to order $n$ as a single matrix-exponential formula whose run time does not grow with the horizon. Separately, Theorem 3 establishes that, under a symmetric parameterization and with $\theta = \frac{1}{\alpha}\sum_i\mathbb{E}[B_i]\uparrow 1$, the Laplace transform of the scaled stationary intensity obeys $\lim_{\theta\uparrow1}\mathcal{T}\{\lambda\}(s(1-\theta))=(\sigma/(\sigma+s))^{\sigma\lambda}$, a multivariate Gamma limit.

Load-bearing premise

The whole construction rests on the joint process $(N(t),\lambda(t))$ being Markov, which holds only when the excitation decay functions are exponential and the jump marks are independent and identically distributed; drop that and the closed form in Theorem 1 is no longer claimed.

Editorial extensions

If this is right

  • Transient and stationary cross-moments of any order, as well as auto- and cross-covariances, can be evaluated exactly by solving linear ODEs or linear algebraic systems, without simulation error.
  • The nested block-matrix construction makes the cost of transient moment evaluation essentially independent of the time horizon, because only a matrix exponential of fixed dimension is involved.
  • Moment-based estimation and comparative statics for multivariate Hawkes models become practical, since each evaluation of a collection of moments is near-instant and exact.
  • In the nearly unstable symmetric regime, the stationary intensity, scaled by $1-\theta$, converges to a Gamma law with covariance $\lambda/\sigma$ between every pair of components.
  • The two-time transform of Theorem 2 supplies the processes' autocovariance functions directly, enabling likelihood-free inference that uses temporal dependence information.

Reading between the lines

Editorial extensions of the paper, not claims the author makes directly.

  • An implication the paper leaves implicit is that the transform-differentiation recipe depends only on the generator being a polynomial in the state variables, so the same scheme should extend to other Markovian intensity models, such as Markov-modulated or affine jump intensities, without reworking the recursion.
  • For non-exponential decay kernels the Markov property fails and Theorem 1 should not be applied; a natural test is to approximate a general kernel by a phase-type (sum-of-exponentials) kernel and compare the moment formulas against a cluster-representation simulation of the true Hawkes process.
  • The nearly-unstable Gamma limit suggests a further, unproved consequence: in the same parameterization the scaled population process $Q(t)$ may also possess a tractable heavy-traffic limit, by analogy with infinite-server queues in overload.
  • Because Theorem 2 gives the full two-time transform, spectral or frequency-domain inference for Hawkes populations is a plausible next step, though the paper does not explore it.
Share X Bluesky LinkedIn Reddit HN

Editorial analysis

A structured set of objections, weighed in public.

Desk editor's note, referee report, and a circularity audit.

Referee Report

2 major / 3 minor

Summary. The paper studies Markovian multivariate Hawkes processes with exponential decay and their induced infinite-server population processes. It derives a closed-form joint transform for (Q(t), λ(t)) via the Markov generator and the method of characteristics (Theorem 1), extends this to two time points (Theorem 2), and then differentiates the transform to obtain linear ODE systems for transient moments (Eq. (21)) and algebraic systems for stationary moments (Eq. (23)). For the bivariate case it develops recursive and nested block-matrix algorithms for moments of arbitrary order, studies the nearly unstable regime under a symmetry assumption (Theorem 3, Corollary 2), and reports numerical experiments comparing the proposed moment computations with finite differences and Monte Carlo simulation.

Significance. If correct, the paper would provide a substantial computational toolkit: exact moment recursions for all transient and stationary cross-moments of Markovian multivariate Hawkes population processes, closed-form autocovariance evaluation via Theorem 2, and an explicit heavy-traffic-type limit for the intensity. The transform results in Theorems 1 and 2 are derived in detail and appear sound. The paper also provides reproducible code, extensive appendices, and careful numerical comparisons. However, the moment recursion claimed for general joint distributions of the intensity jump sizes contains a factorization error that invalidates Eqs. (21) and (23) as stated; this is a load-bearing issue for the paper's central claims, even though the underlying transform theorem remains correct.

major comments (2)
  1. [Appendix B.2, Eqs. (62)-(63) and Eq. (21)] The passage from Eq. (62) to Eq. (63) factorizes E[z_j ∏_i B_{ij}^{n_{λ,i}-m_i} λ_i(t)^{...} z^{Q_i(t)}] incorrectly. Since B_j is independent of (λ,Q), the term equals E[∏_i B_{ij}^{n_{λ,i}-m_i}] · E[z_j ∏_i λ_i^{...} z^{Q_i}], so the coefficient is the joint moment E[∏_i B_{ij}^{n_{λ,i}-m_i}], not the product ∏_i E[B_{ij}^{n_{λ,i}-m_i}] written in Eq. (63) and in Eq. (21). These are equal only when the components of B_j are independent across i. The paper explicitly allows general joint distributions of the intensity jump sizes (Eq. (2), Introduction), and Eq. (18) even assumes finiteness of the joint moment E[∏_i B_{ij}^{n_{λ,i}}], so the factorization is not justified. A concrete consequence: for d=2, nλ=(1,1), nQ=0, j=1, the offending term in Eq. (21) gives E[B11]E[B21] E[λ1] instead of the direct-differentiation result E[B11B21] E[λ1]. The same error propagates into the stationary system Eq. (23), the recursive algorithms of Section 4, and the nested block-matrix construction of Section 5. The numerical experiments in Section 7 use independent exponentially distributed marks, which hide the discrepancy. The fix is local in principle—replace the products of marginal B-moments with the joint mixed moments of B_j—but it affects all moment-vector and block-matrix formulas and must be carried through Sections 4, 5, and Appendices C-D.
  2. [Sections 4-5 and Appendix D] Because the recursive and nested block-matrix procedures are built directly on Eq. (21), their validity for general dependent marks is currently unsupported. The matrices in Section 5 and Appendix D.4 contain objects such as E[B11]E[B21] and E[B_{11}^2]E[B_{21}] where the correct coefficients are joint moments E[B11B21] and E[B_{11}^2 B_{21}]. The presentation should either state and prove the recursions under an explicit independence assumption on the B_{ij} across i for each j, or replace every occurrence of product marginal B-moments by the corresponding joint moments and update the block-matrix formulas accordingly.
minor comments (3)
  1. [Corollary 2, Section 6] The parameterization of the limiting Gamma distribution is inconsistent with Theorem 3. From Eq. (49), the Laplace transform is (σ/(σ+s))^{σλ}, which corresponds to a Gamma distribution with shape σλ and rate σ (or scale 1/σ), not scale σ as stated in Corollary 2 and in the sentence defining Γ(r, λ). The covariance formula in Eq. (50) is consistent with the common-component interpretation, but the corollary should also state explicitly that the limiting vector is degenerate in the sense that all components are equal to a single common Gamma random variable.
  2. [Appendix B.1, sixth term derivation] In the displayed computation for the sixth term after Eq. (56), the second sum is written with a factor k_j ν_j inside the Laplace transform; this should be ν_j, since the k_j μ_j part is already separated. This is a typographical issue that does not affect the result.
  3. [Definition 1, Eq. (2)] The definition of (B_ij(t))_t as 'a sequence of independent random variables' leaves unclear whether independence is assumed across i for fixed j and across j for fixed i. Given that the paper's generality claim depends on the joint law of B_j, it would help to state explicitly what is independent and what may be dependent.

Circularity Check

0 steps flagged · score 0.0 of 10

No significant circularity: the moment systems are derived by differentiating a transform obtained from the Markov generator, with no fitted inputs or load-bearing self-citations.

full rationale

The paper's central derivation is self-contained. Theorem 1 (Eqs. (8)-(9)) is obtained by deriving the generator PDE (58) from the Markov property (Appendix A), applying Laplace/z-transforms and the method of characteristics; Theorem 2 follows by conditioning and repeated application of Theorem 1. The moment recursions (21) and (23) are not fitted or calibrated: they are obtained by repeated differentiation of the PDE (19)/(58) with respect to s and z and setting s=0, z=1 (Appendix B.2), with all expectations of jump marks entering through the transform's beta_j(s) function. The numerical parameters in Section 7 are chosen, not estimated, and the comparisons to finite differences and Monte Carlo are benchmarking exercises rather than predictions. The self-citations to [17] (cluster representation for simulation) and [18] (univariate moment approach being generalized) are contextual: [18] is external co-authored work and the present appendix supplies the actual derivation, while [17] is used only for the Monte Carlo sampling mechanism in the alternative method, not for the main theorem. The concluding remarks explicitly delimit the Markovian exponential-decay setting, so no hidden assumption is smuggled in. The only caveat is that the numerical 'true values' in Section 7 are produced by the paper's own Proposition 4 rather than by an independent external source; this weakens the numerical validation but does not make the analytic derivation circular. Potential algebraic errors in passing from Eq. (62) to Eq. (63) would be correctness defects, not circularity.

Assumptions & free parameters 0 free parameters · 6 assumptions · 0 invented entities

No parameters are fitted to data. The model primitives (base rates, decay rates, service rates, and the distributions of the marks) are inputs, not outputs. The central formulas depend on all of them, and the finite-moment and stability assumptions are explicit. The symmetry assumption in Section 6 is the only ad hoc restriction, and it limits the near-critical theorem to a symmetric parameterization.

assumptions (6)
  • domain assumption Individual marks B_ij(t) are i.i.d. copies of B_ij, independent of the history, and service times are i.i.d. exponential, independent of arrivals.
    Definition 1 and Definition 2 define the model this way; Theorem 1 uses these independence assumptions in the Markov generator step.
  • domain assumption The pair (N(t), lambda(t)) is Markovian under exponential decay.
    Used throughout Section 3, Eq. (58), and the proof of Theorem 1; the paper states this holds only for exponential g_i in Section 2.
  • domain assumption Stability condition rho(H) < 1 with H_ij = E[B_ij]/alpha_i (Assumption 1).
    Required for the stationary moments and for the nearly-unstable limit; used in Section 3.2 and Section 6.
  • domain assumption Finite product moments of jump sizes: E[prod_i B_ij^{n_lambda_i}] < infinity for the orders considered (Eq. (18)).
    Needed to make the derivatives of the transform and the binomial sums in Eq. (21) well-defined.
  • ad hoc to paper Symmetry Assumption 2: alpha_i = alpha, B_1i = ... = B_di = B_i, lambda_i = lambda.
    Imposed for tractability in Section 6; the nearly-unstable Theorem 3 only holds under this parameterization.
  • standard math Interchange of differentiation with expectation and Laplace and z-transform operations in Section 3.2 and Appendix B.
    Standard regularity implicit in the finite-moment assumption; the paper does not state a separate dominated-convergence argument.

how reviews work

0 comments
Cite this review

Pith. "Pith review of Markovian multivariate Hawkes population processes: Efficient evaluation of moments." pith.science (2026). https://pith.science/paper/JKEQSJOV

@misc{pith2026250608775,
  author       = {Pith},
  title        = {Pith review of: Markovian multivariate Hawkes population processes: Efficient evaluation of moments},
  year         = {2026},
  howpublished = {\url{https://pith.science/paper/JKEQSJOV}},
  note         = {Machine review of arXiv:2506.08775}
}
read the original abstract

We provide probabilistic and computational results on Markovian multivariate Hawkes processes and induced population processes. By applying the Markov property, we characterize in closed form a joint transform, bijective to the probability distribution, of the population process and its underlying intensity process. We demonstrate a method that exploits the transform to obtain analytic expressions for transient and stationary multivariate moments of any order, as well as auto- and cross-covariances. We reveal a nested sequence of block matrices that yields the moments in explicit form and brings important computational advantages. We also establish the asymptotic behavior of the intensity of the multivariate Hawkes process in its nearly unstable regime, under a specific parameterization. In extensive numerical experiments, we analyze the computational complexity, accuracy and efficiency of the established results.

Discussion (0). Sign in to comment.

Reference graph

Works this paper leans on

25 extracted references · 22 canonical work pages

  1. [17]

    S., Laeven, R

    Karim, R. S., Laeven, R. J. A. & Mandjes, M. (2021). Exact and asymptotic analysis of gen- eral multivariate Hawkes processes and induced population processes. Preprint. Available at https: //arxiv.org/pdf/2106.03560.pdf

    arXiv 2021

  2. [1]

    A¨ıt-Sahalia, Y., Cacho-Diaz, J. A. & Laeven, R. J. A. (2015). Modeling financial contagion using mutually exciting jump processes. Journal of Financial Economics 117, 585-606

    work page 2015

  3. [2]

    A¨ıt-Sahalia, Y., Laeven, R. J. A. & Pelizzon, L. (2014). Mutual excitation in Eurozone sovereign CDS. Journal of Econometrics 183, 151-167. MULTIV ARIATE HA WKES PROCESSES 25

    work page 2014

  4. [3]

    & Massouli´e, L

    Br´emaud, P. & Massouli´e, L. (1996). Stability of nonlinear Hawkes processes. Annals of Probability 24, 1562-1588

    work page 1996

  5. [4]

    & Hong, G

    Chen, X. & Hong, G. (2023). Steady-state analysis and online learning for queues with Hawkes arrivals. Preprint. Available at https://arxiv.org/pdf/2311.02577.pdf

    arXiv 2023

  6. [5]

    W., Qu, Y., Surya, B

    Chen, Z., Dassios, A., Kuan, V., Lim, J. W., Qu, Y., Surya, B. & Zhao, H. (2020). A two-phase dynamic contagion model for COVID-19. Preprint. Available athttps://www.lse.ac.uk/Statistics/ Assets/Documents/Paper-Covid-19-DCP.pdf

    work page 2020

  7. [6]

    Cui, L., Hawkes, A. & Yi, H. (2020). An elementary derivation of moments of Hawkes processes. Advances in Applied Probability 52, 102-137

    work page 2020

  8. [7]

    & Zaatour, R

    Da Fonseca, J. & Zaatour, R. (2015). Clustering and mean reversion in a Hawkes microstructure model. Journal of Futures Markets 35, 813-838

    work page 2015

Show all 25 references
  1. [8]

    Daley, D. J. & Vere-Jones, D. (2003). An Introduction to the Theory of Point Processes. Volume 1: Elementary Theory and Methods . Second Edition. Springer Science and Business Media, New York

    2003

  2. [9]

    & Zhou, H

    Dassios, A. & Zhou, H. (2011). A dynamic contagion process. Advances in Applied Probability 43, 814-846

    2011

  3. [10]

    & Gruendlinger, L

    Daw, A., Castellanos, A., Yom-Tov, G., Pender, J. & Gruendlinger, L. (2025). The co- production of service: Modeling service times in contact centers using Hawkes processes. Management Science 71, 1865-1888

    2025

  4. [11]

    & Pender, J

    Daw, A. & Pender, J. (2023). Matrix calculations for moments of Markov processes. Advances in Applied Probability 55, 126-150

    2023

  5. [12]

    & Pender, J

    Daw, A. & Pender, J. (2018). Queues driven by Hawkes processes. Stochastic Systems 8, 192-229

    2018

  6. [13]

    & Luo, D

    Du, D. & Luo, D. (2019). The pricing of jump propagation: Evidence from spot and options markets. Management Science 65, 2360-2387

    2019

  7. [14]

    Hawkes, A. G. (1971). Point spectra of some mutually exciting point processes. Journal of the Royal Statistical Society. Series B (Methodological), 438-443

    1971

  8. [15]

    Hawkes, A. G. (1971). Spectra of some self-exciting and mutually exciting point processes. Biometrika 58, 83-90

    1971

  9. [16]

    Hawkes, A. G. & Oakes, D. (1974). A cluster representation of a self-exciting process. Journal of Applied Probability 11, 493-503

    1974

  10. [18]

    T., Saxena, M., Boxma, O

    Koops, D. T., Saxena, M., Boxma, O. J. & Mandjes, M. (2018). Infinite-server queues with Hawkes input. Journal of Applied Probability 55, 920-943

    2018

  11. [19]

    J., Taimre, T

    Laub, P. J., Taimre, T. & Pollett, P. K. (2015). Hawkes Processes. Preprint. Available at https: //arxiv.org/abs/1507.02822

    2015 arXiv

  12. [20]

    Liniger, T. J. (2009). Multivariate Hawkes Processes. PhD thesis, ETH Z¨ urich

    2009

  13. [21]

    Oakes, D. (1975). The Markovian self-exciting process. Journal of Applied Probability 12, 69-77

    1975

  14. [22]

    On Lewis’ simulation method for point processes.IEEE Transactions on Information Theory 27, 23-31

    Ogata, Y.(1981). On Lewis’ simulation method for point processes.IEEE Transactions on Information Theory 27, 23-31

    1981

  15. [23]

    Kong, Q., Carman, M

    Rizoiu, M.-A., Mishra, S. Kong, Q., Carman, M. & Xie, L. (2018). SIR-Hawkes: Linking epi- demic models and Hawkes processes to model diffusions in finite populations. Preprint. Available at https://arxiv.org/abs/1711.01679

    2018 arXiv

  16. [24]

    & Whinston, A

    Xu, L., Duan, J. & Whinston, A. (2014). Path to purchase: A mutually exciting point process model for online advertising and conversion. Management Science 60, 1392-1412

    2014

  17. [25]

    −α1 E[B12] E[B21] −α2 # Ψ(0,1) t +

    Zhu, L. (2013). Central limit theorem for nonlinear Hawkes processes. Journal of Applied Probability 50, 760-771. Appendix A. Proofs of Theorems 1 and 2 Proof of Theorem 1. The proof is comprised of a number of steps. First, we use the Markov property on the distribution funct...

    2013

Pith tools

Reviewed August 7, 2026 · model on record in the stance chip above.