A Representer Theorem for Hawkes Processes via Penalized Least Squares Minimization
Authors: Hideaki Kim
Tomoharu Iwata
Venue: ICLR 2026 Oral (2026)
Keywords: Hawkes processes
kernel methods
representer theorem
point processes
least squares loss
Primary area: learning theory

Abstract: The representer theorem is a cornerstone of kernel methods, which aim to estimate latent functions in reproducing kernel Hilbert spaces (RKHSs) in a nonparametric manner. Its significance lies in converting inherently infinite-dimensional optimization problems into finite-dimensional ones over dual coefficients, thereby enabling practical and computationally tractable algorithms. In this paper, we address the problem of estimating the latent triggering kernels--functions that encode the interaction structure between events--for linear multivariate Hawkes processes based on observed event sequences within an RKHS framework. We show that, under the principle of penalized least squares minimization, a novel form of representer theorem emerges: a family of transformed kernels can be defined via a system of simultaneous integral equations, and the optimal estimator of each triggering kernel is expressed as a linear combination of these transformed kernels evaluated at the data points. Remarkably, the dual coefficients are all analytically fixed to unity, obviating the need to solve a costly optimization problem to obtain the dual coefficients. This leads to a highly efficient estimator capable of handling large-scale data more effectively than conventional nonparametric approaches. Empirical evaluations on synthetic and real-world datasets reveal that the proposed method achieves competitive predictive accuracy while substantially improving computational efficiency compared to state-of-the-art kernel method-based estimators.

Published as a conference paper at ICLR 2026
A REPRESENTER THEOREM FOR HAWKES PROCESSES
VIA PENALIZED LEAST SQUARES MINIMIZATION
Hideaki Kim & Tomoharu Iwata
NTT, Inc.
fhideaki. kin, tomoharu. (
ABSTRACT
The representer theorem is a cornerstone of kernel methods, which aim to esti- mate latent functions in reproducing kernel Hilbert spaces (RKHSs) in nonpara- metric manner. Its significance lies in converting inherently infinite-dimensional optimization problems into finite-dimensional ones over dual coefficients, thereby enabling practical and computationally tractable algorithms. In this paper, we ad- dress the problem of estimating the latent triggering cernels-funetion that encode the interaction structure between events-for linear multivariate Hawkes processes based on observed event sequences within RKHS framework. We show that, under the principle of penalized least squares minimization, a novel form of rep- resenter theorem emerges: family of transformed kernels can be defined via a system of simultaneous integral equations, and the optimal estimator of each trig- gering kernel expressed as a linear combination of these transformed kernels evaluated at the data points. Remarkably, the dual coefficients are all analytically fixed to unity, obviating the need to solve a costly optimization problem to ob- tain the dual coefficients. This leads to a highly efficient estimator capable of handling large-scale data more effectively than conventional nonparametric ap- proaches. Empirical evaluations synthetic and real-world datasets reveal that the proposed method achieves competitive predictive accuracy while substantially improving computational efficiency compared state-of ofte--oftie of-the-art kernel method- based estimators.
INTRODUCTION
Nonparametric estimation of latent functions remains central topic in both theoretical and ap- plied research, spanning domains such as signal and image processing (Lm 201;; Takeda eta (0077), system control (Luu et all 2018), geostatistics (Chiles Delfiner, 2012),, bioinfor- matics (Schölkopt et all 20TE), and clinical studies (Caller 2023). Among various nonparametric approaches, kernel methods stand out as one of the most powerful and mature frameworks. These methods enable flexible function approximation by embedding data into high-dimensional repro- ducing kernel Hilbert spaces (RKHSs) (Schölkopf & Smola, 2018; Shawe- & Cristianini] 2004)). In classical supervised settings with i.i.d. data, the representer theorem plays a pivotal role in kernel methods. It states that the solution to a broad class of infinit-ddimensons optimization problems RKHSs admits dimensiaill representation: the optimal function estimator can be expressed as a linear combination of kernel functions evaluated at the training points (Schölkopl 2C0I: Wahbal 1990). This linear form not only provides theoretical insight but also brings practical advantages in optimization and inference.
Recently, the kernel method literature has begun to address the nonparametric estimation of intensity functions in point process models. The problems are fundamentally more challenging than i.i.d. cases, primarily because the loss functions to minimize (e.g., negative log-likelihood functions) involve integrals of latent intensity functions over observation domains and violate independence assumptions, which renders classical representer theorems inapplicable. A seminal contribution by Flaxman et al (2017) demonstrated that representer theorem still hold for point process: specifically, they showed that if the square root of the intensity function lies an RKHS, then the solution to the penalized maximum likelihood estimation problem admits finite-dimensional

---

Published as a conference paper at ICLR 2026
representation. Interestingly, the optimal estimator is expressed not via standard RKHS kernels but via equivalent kemels-RKHS kernels transformed through a Fredholm integral equation. The result has since been extended to settings with covariate-dependerr intensity functions (Kimettal 2022) and survival point processes (Kim, 2023). providing a broader foundation for kernel method-based
learning in point process models.
More recently, Bonnet & Sangnerr (122 addressed a more intricate setting of multivariate Hawkes
processes (Bremaud & Massonlie 1996: Hawkes, 192T), which offer powerful framework for modeling self- and mutually-interactie
event dynamics in real-world applications such as finance (Bacry et all 2015). neuroscience (Gerhard al, 2012), social networks (Zhou all, 20TS), and seismology (Ogata 1988). By leveraging the approximations of both the log- ligllike and least- squares loss functions, they obtained representer theorem for the estimation of triggering kernels in an RKHS. To ensure the non-negativity of the intensity functions, the model employs non-linear link functions, allowing to capture both excitatory and inhibitory interactions. Although the method demonstrates strong empirical performance, it requires solving non-linear optimization problem over dual coefficients whose dimensionality scales with the data size, posing serious scalability issues for large scale of datasets often seen in multivariate Hawkes processes.
In this paper, we consider a kernel method-basee least squares loss formulation for estimating la- tent triggering kernels in linear multivariate Hawkes processes, where the identity link function is assumed. By leveraging Mercer's s theorem, we establish a novel representer theorem tailored to the functional optimization problem: the obtained estimator of each triggering kernel admits linear expansion in terms of equivalent kernels defined through a system of Fredholm integral equations. Notably, all dual coefficients are analytically fixed to unity, eliminating the need solve a costly coefficient optimization problem. To the best of our knowledge, this paper is the first to establish a representer theorem for the non-a -anaappooximteed penalized least squares formulation of linear Hawkes processes. Furthermore, we propose an efficient algorithm to solve the integral equations using the random feature map approximation of RKHS kernels (Rahimi & Rechi, ((M))), where all required integrals are obtained in closed form, in contrast to the Bonnet & Sangnier (2025) model that re- lies on Riemann approximation. Consequently, the proposed estimator consists solely of additive matrix operations and an inversion of a matrix whose size is independent of the data size. This yields a highly lightweight and scalable estimator that remains effective even on large-scale event data, offering practical and theoretically grounded solution for learning in multivariate Hawkes processes.
PROPOSED METHOD
2.1 PRELIMINARY: LINEAR HAWKES PROCESSES
We consider a multivariate linear Hawkes process (Bremaud & Massoulie 1996 Hawkes 1971) on a time domain R4, i.e., a U-dimensional counting process (N) (t) characterized by the following conditional intensity functions:
(1(1) = + El (o(1-s)dN,(6)) LER,, i iEU = [1, JELL"
(1)
where μ4 E R.+, denotes the baseline intensity for dimension i, and gtr(1 8) RR I is the triggering kernel quantifying the change in the dimension i's intensity at time t caused by the event of dimension j occurring at time
Let {(n,, un) E R4 x N(T) denote sequence of N(T) - Eiez (T) observed events over an interval []. where each pair (tn, un) indicates that the n-th event occurred time t, on dimension u,, In the literature on point processes, two common approaches have been used to estimate the intensity functions: one based on the negative log- likelihoodd function (Daley & Vere- lones, 2001),, and the other on the least squares contrast (Hansen et all, 2015), defined respectively as
LL = E[[ Ai(t)dt- Σ (oE(Ac(t.)). Lis = [[[ x1t Σ selta) (2) NELL
HEA:
SEU
REN,
'Codes available https:: //g/ thub. htt.bb..ooo..iiiikk http///Hiii//////e/eee/ee
2

---

Published as a conference paper at ICLR 2026
denotes subset of event indices associated with dimension Notably, the least squares contrast, Lts, arises from the principle of empirical risk minimization (van Geen. 2000,, and has recently attracted attention due its favorable computational properties in Hawkes process modeling (Bacry et all 2020 Ca 202E!).
where N₁ (,; (T)
2.2 A REPRESENTER THEOREM FOR LINEAR HAWKES PROCESSES
Let T R be positive semi-definite kernel on a dieedmmmnniinnll compact space C R. Then there exists unique reproducing kernel Hilbert space (RKHS) H₂ (Schölkopt & Smola, 2018; Shawe- Bhylor Cristiamini 2004) associated with RKHS kernel k(.,))
Given an observed sequence of events {(t,, u, [n) N(2) over an interval [0, we consider the fol- and baseline
lowing regularized optimization problem of triggering kernels, g intensities, μ (14)iuu in the linear Hawkes process (II)
L(g. 1 Σ 19 (k.))E?
P = arg min gEAd pERu
(3)
where L represents the loss functional, -1 represents the squared Hilbert space norm, and 7 E 14 represents the regularization hyperparameter. In this paper, we adopt the least squares contrast for point processes, denoted by (Ls(g.H) H)), a loss functional, which takes quadratic form in terms of the triggering kernels and baseline intensities,
Les(g,
Σ EEN
(++ Σ nEN
Iian (t. bn (octt-nan)))) )av
2
Σ 14+ +
Σ REN
gia (tw - tn loct ( tnn ).
N'EN,
where ( denotes the indicator, and W [m] In represents the whole indices of events observed.
N(T)
Here, finite support window E RR for the triggering kernels introduced, as commonly done in Hawkes process modeling to reduce computational cost (Bonnet ell 2023: Halp 2013). Theorem D establishes a novel representer theorem for the functional optimization problem defined in Equations (EH+). Notably, all dual coefficients are analytically fixed unity, eliminating the need to optimize them. We present a proof sketch based on Mercer's theorem (Mercen, 1909), with the complete proof deferred Appendix Bl For additional intuition,
also provide alternative derivation via the path integral representation (Kim, 2021) in Appendix D
Theorem 1. Given the estimation of the baseline intensity {iu heu, the solutions of the functional optimization problem (B-E),, denoted as (9ig (oss())),eeee involve the representer theorem under a set of equivalent kernels (h,(:-))) and their dual coefficients are equal to unity:
(ir(s) Σ hn(s,tt Ps (-" hy(s,t)dt, = SET, (i.j) U2, neN.
J.
(5)
where {a 11 denote the dual coefficients, and the equivalent kernels (hy(: )eeu are defined through system of Fredholm integral equations,
1 hn,(, + Σ
Vil(s, (t)
Σ nEN
k(s, tn))
tEll"
(6)
Vs((s,t)= Σ Σ
k(s,t + tw tn 1 Sat.
nEN, n'EM
Proof sketch. Through Mercer's theorem, the RKHS kernel kc(;)) can be expressed through its Mer- cer's expansion:
k(1,
em (1)mm (s), f
een(()) (1)mm (I)dl Tmm 6mmn"
Following Flaxman al (2)1)), call the transformed kernel functions where representer theorem holds the equivolent kernels.

---

Published as a conference paper at ICLR 2026
where fem(-))m=) and {m Jm denote the eigenfunctions and the eigenvalues, respectively, of the integral operator J k(t, s)ds. Accordingly, the triggering kemelss in the RKHS, {gi(() E He (i,)ezze, and their squared RKHS norms, |1 admit the representation
Pi (s) = Σ by (s),
(og))
(4,5) EUe
(8)
where b {b,] R] is the expansion coefficient. problem (L) can be reformulated as follows:
Using this representation, the optimization
b. i arg min Li(b,, + Σ 1
(i,j)E
where
2
LLs(b, RENL
Us + Σ
ba Cm (t tn to) dt
EEN m=1
(10)
2
b 6m (bw' bn -o)
E'EN.
REV m=1
Given the estimate of the baseline intensity jt, the optimal coefficient vector b must satisfy the equation obtained by setting the gradient of the objective with respect to b equal to zero:
DC
Las A)
Σ
Σ
= o, EU, m E (1,
(1I)
067
(KJ) ml
b-
Using the Mercer's expansion leads to Equations (5-66.
( and the kernel trick, k(t,, 8) Em Cm (t) (8), Equation
Note that for any positive semi-definite kernel k(:)),, )), the corresponding of equivalent kernels (hy(--)))eu defined in Equation (6)) exists and uniquely defined. See Appendix G for the proof.
While Bonnet & Sangmer (2025) has also explored representer theorems under the least squares contrast for linear Hawkes processes, their formulation relies on a discretized approximation (see Proposition 1 in (Bonnet & Sangmier, 2025)), which introduces additional optimization over dual coefficients and obscures the elegant mathematical properties of the least squares contrast. To the best of our knowledge, this work is the first establish representer theorem for the non-discretized penalized least squares formulation of linear Hawkes processes.
In Theorem MI the optimal estimators of baseline intensities, {ār jiw are treated given constants. Proposition D demonstrates that by substituting Equation (5) into Equation (L), Hiu can be
obtained closed form in terms of the equivalent kernels. The proof is provided in Appendix DI
The solutions, {patieu of the functional optimization problem (G+) have closed forms in terms of the equivalent kernels defined by Equation (GG) as follows:
Proposition 2.
WW1 -TEEN EWEN, 1o ha,, (t tni tua Doct-ts T T So So Ir.., (1- ln,
A4
i EEu,
(12)
where [Na] denotes the number of observed events associated with dimension i.
2.3 CONSTRUCTION OF EQUIVALENT KERNELS
In Section 2.2, showed that the optimal estimators of g and can be expressed in closed form using the equivalent kernels [hy (., (-.))))))54) yE4. However, obtaining (hy(:-))yeu in practice requires solving the coupled integral equations (E), which is generally non-trivial task. In Proposition we propose solution based on the degenerate kernel approximation methods (Attimsoni 2010:
Polyanin & Manzhirov, (198). The proof is provided in Appendix B
4

---

Published as a conference paper at ICLR 2026
Proposition 3. Let RKHS kernel k(.,-) have degenerate form with M feature maps {om (s)).
M
k(s, = Σ
(m()))0m())) ( (6(5), (6)))
(13)
mmml
where ((s)= (d) (8))T Then the solution of the simultaneous Fredholm integral
equations (0) can be obtained in closed form as follows:
hy(s, =<())) [GTTA
(MCC+=) +3) = 2(5')
jeu,
(14)
140-1MM:M
where Hos denotes the slice of matrix between the a-th row and the b-th one, IM E MUMU0U denotes the identity matrix, E == RMUUMMU defined as a symmetric block matrix whose (i,j)-th th block is given by an M-by-M submatrix,
E,j Σ
Σ
ma af) ( ( max((tast))
Φ(t fn En)((t -t,)) (15)
nEN: N'EN;
and ē(s) : T R MU denotes a concatenated vector function,
((s) (5)(s) (s) 1...) (b,(s)
Σ - tnno-tt) neN,
(16)
Proposition A shows that substituting the equivalent kernels (L4) into Equations (5) and (12). we can obtain the optimal estimators in terms of the feature maps. The proof is provided Appendix H Proposition 4. For degenerate form of RKHS kernel in (ZE)), the optimal estimators, و and i., are obtained in closed form in terms of the feature maps:
(s)(s) <bsss [2) IMU +3)
Σ x(tn) - iss 1' ((()at) neN,
1+G-1))M.jM
(17)
W& (Jo ((1)dt)
(on0))))((H)))) I + +3) = $(1.))
A4 T- (To ((t)dt) (, +3) - "Io" ((t)dt)
In this paper, we assume that RKHS kernels are shift-invariant, i.c., k(s, s') k(ss s'I), which includes popular kernels such as Gaussian, Matérn, and Laplace kernels. We employ the random Fourier feature method (Rahimi Recht, 2007)). approximating the shift-invarian RKHS kernel as sum of Fourier features sampled from the Fourier transform of the kernel, denoted by k((w),
2 (s))) = V M cos(wm + 0m); Wm
k(w) W
m 0n
m < T
M
(18)
M 2
M
m
To enhance the approximation accuracy of the random Fourier features, we employed the quasi- Monte Carlo feature maps (Yang all, 2014), and used M = 100 in Section EL Then the integral operations appeared in (1E) and (177 can be performed analytically follows:
h L. (s)dss V
2 1 M |sin(w min A (n) + 0) - -in E ERM,
(19)
[ Φ(t (n-(t)((tt ) dt b M
[( w 0 + s(w, -w. -0T] E MRM
where o denotes the Hadamard product, = (wi,..., w T. 0 = (01,...,0. and ((w,w',0,0 = cos [0 + () + /2)/2 + 0' wtn - w tw. sinc (bb a)( + w/// (20) Here, sinc(x) = sin(x)/x is the unnormalized sinc function. As result, the optimal estimators obtained closed form without requiring any discretization approximation the integral operators. It is worth noting that the number of feature maps, M, required to approximate RKHS kernel remains modest regardless of the data size, whereas the number of discretization nodes needed for accurate integral evaluations grows with the data size. See Section & for details.
dineinnenii
5

---

Published as a conference paper at ICLR 2026
2.4 COMPLEXITY ANALYSIS
The computational complexity of obtaining our estimators in Equation (LZ) is O +
M$U3), where N. - max((Ni),......)) (o:]]) the first term arises from the computation of and the second from the inversion of IMU + =). Its memory complexity is U( which stems of E. In contrast, the prior kernel method-based estimator (Bonnet & Sangnier, 20255) requires the computation of (NN where P denotes the number of iterations needed for convergence in an iterative optimization algorithm. Its memory complexity (N) Therefore, our ap- proach achieves significantly better scalability with respect to the data size compared to the previous method, making it well-suited for large-scale data scenarios. Moreover, our method requires only a single matrix inversion and avoids the need to carefully tune convergence criteria and learning rates, offering a more stable and practical solution, which is in contrast to the prior kernel method-based method that relies on iterative optimization.
3 RELATED WORK
Hawkes processes (Hawkes, 1971)), particularly in the multivariate setting, have been extensively studied due to their expressive power in modeling self- and mutuall--xcittin temporal dynamics on networks. One of the simplest approaches to learning the triggering kernels in Hawkes processes parametric modeling, where exponential kernels are particularly popular owing to their ability to encode interaction strength and temporal decay compactly. In the case of linear Hawkes processes (W. maximum likelihood estimation has been the gold standard (Bacry et all 2015: Ozaki. 197%
Zhou et all 2013). However, several alternative estimation strategies have been proposed, including approaches that exploit analytic tractability (Bacry et all, 2020), spectral methods (Adamopoulos (976), and moment-matching methods (Da Fonseca & Zaatour X014))
least
Most of the above models assume mutual excitation, i.e., non-negative triggering kernels. Recently. however, Bonnet all (2023) introduced a flexible non-linear Hawkes model (see Equation (2I)) with exponential triggering kernels, which enables us to estimate both excitatory and inhibitory in- teractions efficiently within exponential forms. While the complexity of model fitting scales linearly with the number of events, the model requires processing event times successively to evaluate the likelihood function and is therefore difficult to parallelize across multiple cores, limiting scalability.
In the nonparametric regime, wide variety of methods have been explored for linear Hawkes pro- cesses, which include piece-wise constant Soyynuu-BBurce et 2016) representations of triggering kernels, and the estimation method via the solution of Wiener-Hopf equations (Bacry Muzy, 2016). For non-linear Hawkes processes, nonparametric formulations such as those using Bernstein-type polynomials (Lemonnerr & Vayatıs, am))) and B- spline expansions caal, 2022)) have been proposed. Additionally, many neural network-based models have been developed to learn event dynamics directly from data, ranging from based approaches (Mei & Eisner, 2017) to Transformers (Zuo et all 2020). For a more comprehensive
anl 20) and Gaussian mixture (Xu
survey, we refer readers to (Bonnet & Sangnier, 2025; Bonnet all 2023).
While prior works using reproducing kernel Hilbert spaces (RKHSs) remain relatively underex-
plored in the context of Hawkes processes, two notable exceptions exist. Yang et al. (201 propose an online estimation method RKHS under regret minimization framework, which fundamentally differs from the batch learning setting considered in this paper. Bonnet Sangnier (2025), on the other hand, considered a non-linear multivariate Hawkes process,
xxtt) = Pi
I' gty(t - s))NN(s))) JEU
Σ
LER,
(21)
where ((x) = log(1 e//) is non-negative soft-plus function (w 100), and demonstrated a representer theorem for the problem by adopting an approximate likelihood function (Lemonnier Vayatis, M or an upper bound on the least squares loss (Lemonmer & Vayatıs. 2012) as the
objective. Specifically, they showed that in both cases, the optimal estimator of each triggering kernel admits linear representation in terms of RKHS kernels as follows:
k(,,t t, 1.) /'
0ug a Σ REN,
Σ nEN.
a DJ Σ
k(,,hn tw Moct (22)
R'EA

---

Published as a conference paper at ICLR 2026
where {ag } denote dual coefficients of dimension (T) + U). Rather than solving the associated dual optimization problem directly, they adopted the linear representation (222) as a semi-
parametric model and estimated the dual coefficients by maximizing the objective in (3), where the
integrals in the loss functionals ( were evaluated via discretization approximation. For a detailed
formulation of the resulting objective function, Section 3.3 of (Bonnet & Sangnier, (025)
While the approach ()) exhibits strong empirical performance, it involves solving non-linear optimization problem over dual coefficients, which demands the computation of U() for each evaluation of the objective function, resulting in significant scalability challenges for large-scale datasets, as is often the case multivariate Hawkes processes. What makes the situation worse is that the intensity functions in Hawkes processes are usually discontinuous at the observed event times, and to accurately evaluate the integrals of the intensity functions the loss functionals (2), dense set of discretization nodes required as the number of events increases, As result, the computational cost can grow significantly as the dataset becomes larger.
Here, have primarily focused on kernel method-based Hawkes processes. For discussion of their relationship to another major approach, i.c., deep neural network models, we refer the reader to Appendix II
4 EXPERIMENTS
We evaluated the validity of our proposed method (Ours) by comparing with prior parametric and non-parametric approaches. In accordance with Bonnet & Sangnier (2025), we adopted the following four approaches as baselines: Exp is the sfate-of-the-ard parametric approach based on non-linear Hawkes process with exponential triggering kernels (Bonnet all 2023); Gau is non non-parametric approach based on a linear Hawkes process ((u et allll 2016), where the triggering kernels represented as Gaussian mixtures; Ber is non-parametric approach based linear Hawkes process (Lemonnier& & Vayatis 20))),, where the triggering kernels are represented as Bernstein -type polynomials; Bonnet is the state-of-the-art kernel method-based approach (Bonnet & Sangnier, 20255), which assumes the triggering kernels lie in an RKHS. For Bonnet, the node size of discretization approximation was set at max (100, 2|A1 I........ l following (Bonnet Sangmer, ((255). For Ours and Bonnet, a Gaussian RKHS kernel was employed: k(s, = <1)? where 3 is the inverse scale hyperparameter. For Gau and Ber, the number of basis functions was set at 50. For the models except Exp, the support window A was set at 5.
The four baselines were implemented using the Python code in Bonnet & Sangmerr (2025) (MIT License), and our model using TensorFlow-2.. (Ahadi et all, 2015). All the experiments were executed on a MacBook Pro equipped with 12-core CPU (Apple M2 Max), with the GPU disabled.
4.1 SYNTHETIC DATA
Except for Exp, the models have hyperparameters to optimize. Gau and Ber have the regulariza- tion hyperparameter y for quadratic penalty the coefficients of mixture models, and Ours and Bonnet have the inverse scale hyperparameter 3 in addition to 7, We optimized the hyperparame- ters the grids of 7 E (0.1, 0.5, 1.0) and 3 E {0.5, 1.0, 1.5), based the negative log-likelihood (Gau, Ber, Bonnet) and the least squares loss (Curs) minimization. Specifically, for a sequence of events observed in an interval [o, T]. each model was fitted with the events in [0. 0.8T), evaluated the negative log- log-likelihood/least squares loss for the rest of the data in [0.8T, T), and the hyperpa- rameters minimize the criteria were chosen.
Predictive performance was assessed using the integrated squared error 3 defined as follows:
A
E= NEU JEM
(s)(()) (s)(s)
(23)
where A A= 5, and 9i (s) and (s)(s) denote the true and estimated triggering kernels, respectively. Efficiency was evaluated based on the CPU time, denoted by cpu, required to execute the model fitting given the optimized hyperparameters.
7

---

Published as a conference paper at ICLR 2026
Table 1: Results of Exp (Bonnel all 2023), Gau (Ku all 2016). Ber (Lemonnier & Vayatıs, 2014), Bonnet (Bonnet & Sangnier, 2025), and Ours on mutually-excitngg scenario dataset across 10 trials with standard errors in brackets. N is the average data size per trial. opu is the CPU time in seconds. The performances not significantly (p 2 0.01) different from the best one under the Mann- Mann-Whitney U test (Holm, 1979) are shown in bold.
Exp
Gau
Ber
Bonnet
Ours
T
N
A²
quu
A²
cpu
A²
cpu
A²
cpu
A²
cpu
0.29
46.1
0.20
3.42
0.51
4.06
0.21
135
0.38
3.16
2000
1318
(0.01)
(23.3)
(0.02)
(1.33)
(0.16)
(2.03)
(0.05)
(144)
(0.15)
(2.42)
0.29
74.7
0.19
4.94
0.69
6.61
0.19
272
0.27
4.76
3000
2055
(0.00)
(41.3)
(0.02)
(3.00)
(0.50)
(4.71)
(0.04)
(298)
(0.06)
(3.97)
0.28
134
0.18
12.2
0.30
19.0
0.15
1358
0.20
10.7
5000
4081
(0.00)
(64.9)
(0.01)
(8.91)
(0.06)
(15.4)
(0.02)
(1270)
(0.06)
(7.71)
0.28
180.4
0.18
17.3
0.27
29.0
0.14
2070
0.16
13.1
7000
5380
(0.00)
(56.6)
(0.02)
(8.00)
(0.04)
(16.1)
(0.04)
(1210)
(0.04)
(5.55)
4.1.1
MUTUALLY-EXCITING SCENARIO
We consider synthetic data generated from 3-variate linear Hawkes process (I) with baseline in- tensities, 14 = 0.01 for E U, and mutually-exciting triggering kernels (g:; > 0) defined as follows:
911 (8) = 0.5e
(12(s) = 0.5e -10ss-1
(1a($) = 0.5e -20e-20(e-33)
921 (8) 2-5s-1
(22 0.3e
g2al (s) = 0.5e-20(e-
(24)
9si (s)(s) (s) 0.2e -3.2--((s22)
(s2(s) 0.25(1 + cTs())
933 (s) 0.5e
of which setting is a modification of the one that appeared in Bonnet & Sangmer (2025)) We simu-
lated 10 trial sequences of events over the interval (0, T and performed the estimation of triggering kernels 10 times using the compared methods. To clarify the model efficiency regarding data size, set the horizon at 7 E [2000, 3000, 7000).
Table W displays the predictive error and computational efficiency on the mutually-exciting scenario dataset. Some estimation results are displayed in Appendix Al The results demonstrate that, given sufficiently large amount of data, our proposal achieved comparable predictive accuracy to Bonnet, the SOTA kernel method-based approach, while requiring several orders of magnitude less compu- tation time for model fitting. In small data regimes, Bonnet tended to achieve higher accuracy than Ours, which may be attributed the theoretical advantage of negative log likeeliodood over least squares loss in reducing estimation biases (Bacry et al, 2016). Gau consistently achieved high pre- dictive accuracy with small computation time, because for mutually-exciting interactions, aligned with the underlying process. Exp performed well only for the exponential triggering kernels (ga(.) (-eeu (see Figure M).
the only baseline specifically designed
We additionally conducted experiments investigate the relationship between the estimation accu- racy of the latent triggering kernels (measured by A and the predictive performance of the point process model (measured by the negative aigliiklihood). The results provided in Appendix H
4.1.2 REFRACTORY SCENARIO
Refractory phenomena arise in point processes exhibiting short-ters self-inhibition, where the oc- currence of an event temporarily suppresses the likelihood of subsequent events. Such behavior is observed in neuronal spike trains (Berry & Meister, (1997) and in sequences of mainshock events
(Rotondi & Varini, 2019). Here, we consider synthetic data generated from a 3-variate non-linear Hawkes process (4 l) with short-term self-inhibition adopted in (Bonnet & Sangneer (025).
(s) = (8s² - 1)1.50.5 + 2.5(s 1. 1.00.5.
(25)
(22(s) (s) (852 1)1, 10.5 (s-00.5
1>00.5
and various non-inhibitory inter-interactions.
912(s) = 0.6e - gia(5) 0.8e
18e-2005 20(-
y2l (8) = 0.6 -5e
(26)
(23(s) = 0.8e -0(
gai (s) = 0
(32(s) = 0
8

---

Published as a conference paper at ICLR 2026
Table 2: Results on refractory scenario dataset across 10 trials. Notations follow Table II
Exp
Gau
Ber
Bonnet
Ours
T
A
cpu
A
2
Gpu
A²
cpu
A²
cpu
2.19
124
1,51
7.03
1.23
10,8
0.63
413
0.95
5.04
2000
2050
(0.03)
(61.3)
(0.02)
(3.58)
(0.31)
(6.19)
(0.19)
(291)
(0.24)
(3.92)
2.19
183
1,50
10.8
0.96
18.5
0.56
927
0.85
7.66
3000
2956
(0.02)
(47.9)
(0.02)
(4.11)
(0.14)
(7.58)
(0.22)
(460)
(0.19)
(4.74)
2.20
355
1.47
24.7
0.82
44.0
0.44
3197
0.59
14.9
5000
5222
(0.02)
(74.0)
(0.01)
(8.41)
(0.16)
(15.8)
(0.18)
(1323)
(0.13)
(5.48)
2.20
503
1,47
37.3
0.71
74.7
0.41
5884
0.50
16.1
7000
6887
(0.02)
(108)
(0.01)
(14.5)
(0.09)
(29.0)
(0.19)
(2972)
(0.07)
(5.26)
Table 3: Performance of Ours regarding the support window mutually-exciting scenario dataset (T = 5000) across 10 trials with standard errors in brackets.
A
A
2
A
10
A
20
A
cpu
A
cpu
cpu
A
cpu
A
cpu
0.51
1.75
0.24
3.89
0.20
10.7
0.29
22.5
0.33
39.1
(0.04)
(1.24)
(0.03)
(2.74)
(0.06)
(7.71)
(0.11)
(16.9)
(0.07)
32.4
The remaining experimental conditions follow those of the mutually-excttnng scenario. We adopted the soft-plus function as the link function, which is consistent with the assumption in Bonnet.
Table @ presents the predictive error and computational efficiency on the refractory scenario dataset. Some estimation results are displayed in Figure @ As shown, Ours was outperformed by Bonnet terms of accuracy on small datasets, but the gap became less significant as the dataset size creased. In contrast to the mutually-excitin scenario, Gau performed poorly here due to its inability to model inhibitory interaction. While Ber succeeded in reconstructing the inhibitory interactions. it still fell short of the kernel method-based approaches in terms of accuracy. Increasing the com- ponent number may improve Ber's performance, but at the cost of higher computational time. Our proposed method was the fastest, achieving speed-up of several hundred times compared to the SOTA Bonnet. Note that Ber exhibited unstable behavior in some trials, where A exceeded 1000. Such outlier samples
excluded from the results in Table D
Our proposed method, based on a linear Hawkes process, as well as the approach based on a non- linear one (Bonnet & Sangnier, (02E), estimate inhibitory effects, i.e., triggering kernels that take negative values (see the estimated triggering kernels on the diagonal of Figure 2). However, when using the estimated triggering kernels to predict intensity function X(-), our method requires post-hoc clipping, such as applying max(}(-), )) for small to ensure the non-negativity of the intensity function. Note that the computational overhead this post-processing is negligible.
4.1.3 EFFECTS OF SUPPORT WINDOW ON PERFORMANCE
The support window A can be regarded as hyperparameter that controls the shape of the triggering kernel and can therefore be selected from data, for example, via cross-validatie validdaii...... In general, if A too small, estimation methods cannot capture the true shape of the underlying triggering kernel. Conversely, if A is too large, non-negligible values remain in regions where the true triggering kernel essentially zero, thereby increasing the estimation error. It is worth noting that choosing A too small is substantially more detrimental than choosing A too large; therefore, it is preferable to set A on the larger side.
To verify this behavior, we evaluated the predictive performance of the proposed method Ours for the support window A {1, 10, 20), using data from the mutually-exciting scenario (I' 5000). Here, predictive performance was assessed using the integrated squared emror (A²) defined follows: 2 Eieu Ejeu 20 (9x]s)) - ds. Table B summarizes the results, suggesting that optimal value of A exists.

---

Published as a conference paper at ICLR 2026
Table 4: Results on financial dataset across 10 trials with standard errors in brackets. nll is the negative oggllikelhood on test data. N is the data size per trial. cpu is the CPU time in seconds. The other notations follow Table
Exp
Gau
Ber
Bonnet
Durs
vell
qpu
nll
cpu
rill
cpu
all
cpu
vell
Gpu
8086
38.2
2899
7.17
1487
13.3
1829
550
1721
3.4
3319
(720)
(15.8)
(1015)
(0.69)
(1291)
(1.40)
(1147)
(71.3)
(896)
(1.43)
4.2 REAL-WORLD DATA
We conducted an experiment on financial dataset that is widely used for evaluating Hawkes process models etll 2016). This dataset contains transaction records of single stock over one day, with two event types (U = 2): "buy" and "sell". The event sequence is further partitioned by timestamps. From the 100 sequences available on the GitHub repository of Zuo al (2020) each of which size is 3319, we randomly constructed 10 pairs of sequences; for a pair, one was used for model training, and the other was used to evaluate the negative log- legliikeihood as the predictive error (the lower, the better). For our proposed method, we applied post-hoc clipping, max max(ā(t), 10 to the estimated intensity function. For the models except Exp, the support window A was set at 3, and the hyperparameters were optimized on the grids of E {0.01, 0.101,01,117 1.0) and B E (0.5, 1.0, 10), based on the negative log-likelihood minimization. All other procedures followed those described in Section
Table ! summarizes the results: compared to the baseline methods, approach achieves robust performance while remaining computationally efficient.
CONCLUSION
We have proposed penalized least squares loss formulation for estimating triggering kernels in multivariate Hawkes processes that reside in RKHS. We demonstrated that novel represen- ter theorem holds for the optimization problem and derived a feasible estimator based on kernel methods. We evaluated the proposed estimator synthetic and real-world data, confirming that achieved comparable predictive accuracy while being substantially faster than the state-of-the-ar kernel method estimator. To the best of our knowledge, this is the first example in the long history of kernel methods where representer theorem holds under RKHS kernels defined by a system of integral equations. We therefore view our result as unique and potentially impactful for the kernel methods community.
Limitations:
Our proposed method is based on a linear Hawkes process. which does not guarantee the non-negativity of the intensity function. As result, when using the estimated triggering kernel to predict future intensity values, post-hoc clippings such as applying max (ax({x(), c) for a small are required. it is known priori that the underlying triggering kernels are excitatory, it is more appropriate to enforce non-negativity directly on the estimated triggering kernels.
The computational complexity of our method scales cubically with the dimensionality U of the Hawkes process, making it less suitable for high-dimensional processes, while it empirically works robustly for moderate (up to a few hundred) dimensions (see Appendix H
A rigorous theoretical analysis of the learning guarantees is beyond the scope of this study and
remains an important direction for future research.
10

---

Published as a conference paper at ICLR 2026
REPRODUCIBILITY STATEMENT
We provide detailed implementation instructions and reproducibility guidelines in Section EL and the Python code is provided
THE USE LARGE LANGUAGE MODELS (LLMs)
We used LLMs to polish writing and correct typos.
REFERENCES
Martín Abadi, Ashish Agarwal, Paul Barham, Eugene Brevdo, Zhifeng Chen, Craig Citro, Greg S. Corrado, Andy Davis, Jeffrey Dean, Matthieu Devin, Sanjay Ghemawat, Ian Goodfellow, Andrew Harp, Geoffrey Irving, Michael Isard, Yangqing Jia, Rafal Jozefowicz, Lukasz Kaiser, Manjunath Kudlur, Josh Levenberg, Dandelion Mané, Rajat Monga, Sherry Moore, Derek Murray, Chris Olah, Mike Schuster, Jonathon Shlens, Benoit Steiner, llya Sutskever, Kunal Talwar, Paul Tucker, Vincent Vanhoucke, Vijay Vasudevan, Fernanda Viégas, Oriol Vinyals, Pete Warden, Martin Wat- tenberg. Martin Wicke, Yuan Yu, and Xiaoqiang Zheng. TensorFlow: Large-scale machine learn- ing on heterogeneous systems, 2015. URL nttps:/ tensortion org/. Software available from tensorflow.org.
L. Adamopoulos. Cluster models for earthquakes: Regional comparisons. Journal of the Interna- tional Association for Mathematical Geology, 8:463-475, 1976.
Kendall Atkinson. A personal perspective on the history of the numerical analysis of Fredholm integral equations of the second kind. In The Birth of Numerical Analysis, pp. 53-72. World Scientific, 2010.
Emmanuel Bacry and Jean-François Muzy. First-and second-order statistics characterization of Hawkes processes and non-parametric estimation. JEEE Transactions on Information Theory, 62 (4):2184-220)) 2016.
Emmanuel Bacry, lacopo Mastromatteo, and Jean-François Muzy. Hawkes processes in finance. Market Microstructure and Liquidity, 1001):155005, 2015.
Emmanuel Bacry, Stéphane Gaiffas, lacopo Mastromatteo, and Jean-François Muzy. Mean-field inference of Hawkes point processes. Journal of Physics A: Mathematical and Theoretical, 49 (17):174006, 2016.
Emmanuel Bacry, Martin Bompaire, Stéphane Gaiffas, and Jean-François Muzy. Sparse and low- rank multivariate Hawkes processes. Journal of Machine Learning Research, 21(50)11-32, 2020.
Michael Berry and Markus Meister. Refractoriness and neural precision. In Advances in Neural Information Processing Systems 10, 1997.
Anna Bonnet and Maxime Sangnier. Nonparametric estimation of Hawkes processes with RKHSs. In Artificial Intelligence and Statistics, pp. 3574-3582. PMLR, 2025.
Anna Bonnet, Miguel Martinez Herrera, and Maxime Sangmier. Inference of multivariate exponen- tial Hawkes processes with inhibition and application to neuronal activity. Statistics and Comput- ing, 33(4):91, 2023.
Pierre Brémaud and Laurent Massoulié. Stability of nonlinear Hawkes processes. The Annals of Probability, pp. 1163-1588, 1996.
Biao Cai, Jingfei Zhang, and Yongtao Guan. Latent network structure learning from high- dimensional multivariate point processes. Journal of the American Statistical Association, 119 (545):95-108, 2024.
Jean-Paul Chiles and Pierre Delfiner. Geostatistics: Modeling Spatial Uncertainty, volume 713. John Wiley & Sons, 2012.
David Collett. Modelling Survival Data in Medical Research. Chapman and Hall/CRC, 2023.

---

Scandinavian Journal of Published as a conference paper at ICLR 2026
José Da Fonseca and Riadh Zaatour. Hawkes process: Fast calibration, application trade cluster- ing, and diffusive limit. Journal of Futures Markets, 34(6):548- 0708885772200 2014.
Daryl J. Daley and David Vere-Jones. An Introduction to the Theory of Point Processes: Volume I: Elementary Theory and Methods. Springer Science & Business Media, 2006.
Nan Du, Hanjun Dai, Rakshit Trivedi, Utkarsh Upadhyay, Manuel Gomez-Rodriigue,, and Le Song- Recurrent marked temporal point processes: Embedding event history to vector. In ACM SIGKDD International Conference on Knowledge Discovery and Data Mining, pp. 1555-1564, 2016.
Seth Flaxman, Yee Whye Teh, and Dino Sejdinovic. Poisson intensity estimation with reproducing kemells. In Artificial Intelligence and Statistics, pp. 270-279. PMLR, 2017.
Felipe Gerhard, Moritz Deger, and Wilson Truccolo. On the stability and dynamics of stochastic spiking neuron models: Nonlinear Hawkes process and point process GLMs. PLoS Computa- tional Biology. (5())ee005390 2017.
Peter F. Halpin. A scalable EM algorithm for Hawkes processes. In New Developments in Quantita- tive Psychology: Presentations the 77th Annual Psychometric Society Meeting, pp. 403-414. Springer, 2013.
Niels Richard Hansen, Patricia Syynaud-Boured and Vincent Rivoirard. Lasso and probabilistic inequalities for multivariate point processes. Bernoulli, 21(1):83-143,, 2015.
Alan G. Hawkes. Spectra of some self-exciting and mutually exciting point processes. Biometrika, 58(1):83-90, 1971.
Sture Holm. A simple sequentially rejective multiple test procedure. Statistics, pp. 65-70, 1979.
Hideaki Kim. Fast Bayesian inference for Gaussian Cox processes via path integral formulation. In Advances in Neural Information Processing Systems 34, 2021.
Hideaki Kim. Survival permanental processes for survival analysis with time-varying covariates. In Advances in Neural Information Processing Systems 36, 2023.
Hideaki Kim, Taichi Asami, and Hiroyuki Toda. Fast Bayesian estimation of point process intensity as function of covariates. In Advances in Neural Information Processing Systems 35, 2022.
Remi Lemonnier and Nicolas Vayatis. Nonparametric Markovian learning of triggering kernels for mutually exciting and mutually inhibiting multivariate Hawkes processes. In Machine Learning and Knowledge Discovery in Databases, pp. 161-176. Springer, 2014.
Miao Liu, Girish Chowdhary, Bruno Castra Da Silva, Shih-Yuan Liu, and Jonathan P. How. Gaussian processes for learning and control: A tutorial with examples. IEEE Control Systems Magazine, 38(5):53-B8) 86, 2018.
Weifeng Liu, Jose C. Principe, and Simon Haykin. Kernel Adaptive Filtering: Comprehensive Introduction. John Wiley & Sons, 2011.
Hongyuan Mei and Jason M. Eisner. The neural Hawkes process: A neurally self-modulating mul- tivariate point process. In Advances in Neural Information Processing Systems 30, 2017.
James Mercer. XVI. Functions of positive and negative type, and their connection with the theory of integral equations. Philosophical Transactions of the Royal Society of London, Series A, 209: 415-446, 1909.
Yosihiko Ogata. The asymptotic behaviour of maximum likelihood estimators for stationary point processes. Annals of the Institute Statistical Mathematics, 30:243-266 1978.
Yosihiko Ogata. Statistical models for carthquake occurrences and residual analysis for point pro- cesses. Journal of the American Statistical Association, (3(401):9-27, 1988.
Tohru Ozaki. Maximum likelihood estimation of Hawkes' self-exciting point processes. Annals of the Institute of Statistical Mathematics, 31:145-155, 1979.
12

---

Published as a conference paper at ICLR 2026
Andrei D. Polyanin and Alexander V. Manzhirov. Handbook of Integral Equations. CRC press, 1998.
Ali Rahimi and Benjamin Recht. Random features for large-scale kernel machines. In Advances in Neural Information Processing Systems 20, 2007.
Patricia Reynaud- Syynaud-Bouret, Vincent Rivoirard, Franck Grammont, and Christine Tuleau-Mable Goodness-offfi tests and nonparametric adaptive estimation for spike train analysis. The Journal of Mathematical Neuroscience, 4:1-41, 2014.
Renata Rotondi and Elisa Varini. Failure models driven by a self-correcting point process in earth- quake occurrence modeling. Stochastic Environmental Research and Risk Assessment, 33(3): 709-724, 2019.
Bernhard Schölkopf and Alexander J. Smola. Learning with Kernels: Support Vector Machines, Regularization, Optimization, and Beyond. MIT press, 2018.
Bernhard Schölkopf, Ralf Herbrich, and Alex Smola. A generalized representer theorem. In International Conference Computational Learning Theory, pp. 416-426. Springer, 2001.
Bernhard Schölkopf, Koji Tsuda, and Jean- Cean-Philippe Vert. Kernel Methods in Computational Biology. MIT press, 2004.
John Shawe-Tayee and Nello Cristianini. Kernel Methods for Pattern Analysis. Cambridge Univer- sity Press, 2004.
Hiroyuki Takeda, Sina Farsiu, and Peyman Milanfar. Kernel regression for image processing and reconstruction. IEEE Transactions on Image Processing, 16(2):349-366 2007.
Sara de Geer. Empirical Processes in M-estimation, volume 6. Cambridge University Press, 2000.
Grace Wahba. Spline Models for Observational Data, volume 59. SIAM, 1990.
Andrew Gordon Wilson, Zhiting Hu, Ruslan Salakhutdinov, and Eric P Xing. Deep kernel learning. In Artificial Intelligence and Statistics, pp. 370-378. PMLR, 2016.
Hongteng Xu, Mehrdad Farajtabar, and Hongyuan Zha. Learning Granger causality for Hawkes processes. In International Conference Machine Learning, pp. 1717-1726. PMLR, 2016.
Jiyan Yang, Vikas Sindhwani, Haim Avron, and Michael Mahoney. Quasi-Monte Carlo feature maps for shift-invariant kernels. In International Conference on Machine Learning, pp. 485-493.
PMLR, 2014.
Yingxiang Yang, Jalal Etesami, Niao He, and Negar Kiyavash. Online learning for multivariate Hawkes processes. In Advances in Neural Information Processing Systems 30, 2017.
Ke Zhou, Hongyuan Zha, and Le Song. Learning social infectivity in sparse low-rank networks using multi-dimensional Hawkes processes. Artificial Intelligence and Statistics, pp. 641-649. PMLR, 2013.
Simiao Zuo, Haoming Jiang, Zichong Li, Tuo Zhao, and Hongyuan Zha. Transformer Hawkes process. In International Conference on Machine Learning, pp. 11692-11700 1702. PMLR, 2020.
13

---

Published as a conference paper at ICLR 2026
A
EXAMPLES OF THE ESTIMATED TRIGGERING KERNELS ON SYNTHETIC
DATA
0.6
True 91
True 912
Ours
Exp
True 913
0.4
Bonnet
Ber
0.2
Gau
0.0
0.6
True 921
True 922
True 923
0.4
0.0
0.6
True 931
True 932
True 93s
0.2
0.0
3
0
Figure 1: Examples of the estimated triggering kernels in the mutually-exciine scenario. Dashed lines represent the true triggering kernels.
True 9
True 9
True 9,3
True 921
True 922
True 923
0
1
True 931
True 9yp
True 9
0
0
2
5
3
5
Figure Examples of the estimated triggering kernels in the refractory scenario. Dashed lines represent the true triggering kernels.
14

---

Published as a conference paper at ICLR 2026
B PROOF OF THEOREM II VIA MERCER'S THEOREM
Proof. Through Mercer's theorem, the RKHS kernel k);,,) can be expressed through its Mercer's expansion:
em (m)))mm()) J, e (1)()) (t)dt 1/m Omm':
k(L,8)
where f-m(-))m and {m Jm denote the eigenfunctions and the eigenvalues, respectively, of the integral operator J kitt, s)ds. Accordingly, the triggering kernels in the RKHS, {gi(()) Hy (b.j)ee,, and their squared RKHS norms, |19h) admit the representation
96 = E by (s),
1195 = (og)), (4,3) e
where b (b]) E R is the expansion coefficient. Using this representation, the optimization problem (L)-E) be reformulated as follows:
b. P = arg min
Ls(b.
1 7 (ii)) m=1
Σ
Σ
(6m
b. 1
where
> iEll
0
Lbs
(1
bm em (t tn 10 dt
REN
E
biu em (tn' tn (n)loct
n'EM
Given the estimate of the baseline intensity P. the optimal coefficient vector b must satisfy the equation obtained by setting the gradient of the objective with respect to b equal to zero:
2C
(bs A)
Σ
Σ
mm
o, (1,j) EU, m E (1,2....
(27)
Obm
(KJ)EU m²==1
b=6
Equation (27)) can be written explicitly
a
Tns(b. (1) + Σ
Σ
2 abm
(l,))
m/=1
b=b
T
(in+
8mn' enn bn loct-
Σ
em bn lost (28)
w E.N m'=1
neN,
Σ
Σ
Em (m - tn)t
jm =
E'EN; nEN,
Operating E Em (s) on both sides of Equation (2S) yields the following system Fredholm integral equations of the second kind:
/'
A4 Σ mEN,
k(s, l loct + Σ
Vsi(s, ())a(t)d
(29)
Σ Σ
1 k(s, tn 1000 + (oy(s) ( 0,
-t
n' R'EN nEN;
where used the relation, (s) Σm= im em))) and the kernel trick, k(1, Em= em (t)e here, Vj((s, t) is defined Equation (8). Equation (2 indicates that the opti- mal estimator (15()) admits a linear representation in terms of a set of transformed kernel functions, (1.)n.e as
(s)(s) Σ tig(s,tt) - Po. ''" hig(s,t)dt, s SE T, (i.j) E u2 neA,
15

---

Published as a conference paper at ICLR 2026
where the transformed kernel functions are defined by system of simultaneous Fredholm integral equations,
[nn)),, Σ I. Vuls ()alt,) s')dt LEU
Σ neN,
h(s,8' -ocs''
(30)
Since the coefficients in Equation (6D) are independent of the index its solution (ey().) is also independent of allowing to write as hig(t-) (.)).. This completes the proof.
C PROOF OF THEOREM II VIA PATH INTEGRAL REPRESENTATION
Proof. Let Let (8) Jy k(s, t) be the integral operator with RKHS kernel k(:-). and KTt((s) JT. (t(s, t)dt be its inverse operator. Then, through the path integral representation (Kim, 2021), the squared RKHS norm term can be represented in functional form,
Σ (5)E
Σ
[I (*(s, (5) TXT
(K))E
Using the representation, the functional derivatives of the least squares term and the penalization term in (L), with respect to 9.z(+), can be written as follows:
8T1 (g)
=2 (
E n'EN
Iiun (t - tne Loct Σ nEN,
8(s (t
ö(s (tm'-tn)) Loct) -InS:
n'EN
nEN,
119ss 122e. = 8g4 1
11
(t') 2 / kt (s,
8gts
1.
where 8()) denotes the Dirac delta function. The optimal estimator (.)(-) should solve the equation where the functional derivative the penalized least squares contrast equal zero:
1 Σ
Lis A) +
19955
o,
SET, E..
(5)
1==
Then applying operator K. to the equation leads to the following simultaneous Fredholm integral equations of the second kind:
1 [)s((s) +
Σ f
Vsi(s, (f)
Ield
(31)
[" k(s, sad
k(s,tnr - eis,ttr-tty Σ Σ
to P4 Σ nEN,
N'EN nEN;
where V8((s,t) t) is defined in (6), the second term the left-hand side of Equation (3II) is derived using the following relation,
I'
Σ Σ
Iin., (t tw w to)
N'EN nEN,
ΣΣΣ
by CanMssttHt-aaaao-----sAAeeeeeee
(t- tw + t)
IEU REN, n'EN
Σ
ga(t) Σ
Σ
k(s, + t, by (Ct,, dt,
1EW1
n'Ei
16

---

Published as a conference paper at ICLR 2026
and the relation, (KC = estimator ())()) admits linear representation in terms of a set of transformed kernel functions, as
-0(s - t)dt, used. Equation (+l) indicates that the optimal
(+) Σ hig(8, (n) Ps ! his(s, t)dt, SE,, (i,5) EuP, nEN.
where the transformed kernel functions are defined by system of simultaneous Fredholm integral equations,
I''
hiy(s,
Σ tEt
s')ddt = E s's,'' tt, to)
(32)
nEN,
Since the coefficients in Equation (32) are independent of the index its solution (+y())) is also independent of i, allowing us to write as he)(;,)) = hy(;,)) This completes the proof.
D PROOF OF PROPOSITION
Proof. Given the estimated baseline intensities {Ar hizu, the optimal estimators of the triggering kernels, (ux())),, are determined by the equations,
8
Lssgg (1)
Σ (i))eu
11
=0, SET, (1,5) EPP.
(33)
(s)
9-
These equations are equivalently expressed as:
L (pr (
Diu - tx to))locc-
Σ
8(s - (t tn) tt))
WEN
nEN,
(34)
ΣΣ 8(s (n)) W'EN nEN,
kt (1), =0.
Applying the operator ST (8)) to both sides of Equations (34) yields the following representa- tion of the RKHS regularization term under the estimated triggering kernels:
19:5s12ee = = //
1t (s., (s) (t)dsct
TXT
Σ
Σ
(tr (Fa - to lost SASA
(35)
W'EN nEN;
(1s Σ
Giv.. (t txe) (te lost Σ
Gis(t l lost
W'EN
nEN,
Substituting this representation into the objective function in (B) leads to: Les ()) +
Σ
10
(L))
= || (iit Σ
Biun (t (t-tt tu) )
2
dt + E
Σ
giun (tn - tnt Jost
Fell
neN
W'EN neN
Σ
(Ri
Σ
Dium (ta tn)o-tt
WeV.
WEN
(36)
( Σ
Dim.. (1 lm lostt
Σ neN
giv, (1- (n) Sadt
Σ iEll
The + Pis
Σ
Jrun (t to) 2//
nEN
Σ
Σ
give (twr tn to)
WEN WEN
17

---

Published as a conference paper at ICLR 2026
By invoking the representer theorem (C), the objective can be rewritten function of the estimated baseline intensities:
Z LLs(0(1).)))))
Σ
(
(a.))
T
Σ
Nun (t bys sn
IEU1
EE
(37)
[' hun (t tn tw Juct-
- iui
2//
ΣΣ
RENWEN N'NN
ΣΣ NEN E'EN
/ haw (tw - tn t):ttt/ ti ] +C.
where C the constant term. 7(t)) can be more simplified as
% hap (t tn, 8)
Z(p) =
pu T Σ
iEM
EEN
(38)
20: (Nail
Σ
Σ
[' has tnstn (oct
] +C,
NEA W'E,
where the identity
ΣΣ nEN W'EN
hun, (t - bn )oct-
(39)
ΣΣΣ DEN N'N,
/- Mu, (bw bys toott to - sadl,
used (for proof, see Appendix H). Finally, the optimal estimators Úú should solve the equation where the functional derivative of Z(ū) regarding μæ is equal to zero (dZ/djjz = 0), which leads to the following representation in terms of the equivalent kernels:
LNG EnEN EmEN; So Iu., (I - bns ba loct
A
i E U.
(40)
T EnEN So So Ney (1 Ins 1 S sAdlds
This completes the proof.
E PROOF OF PROPOSITION 3
Proof. For an RKHS kernel (s((.) with the degenerate form given in Equation (L3), the coefficient
functions. (Vx(s,t) (.nn)t) appearing in the system of integral equations (Ki) can also be expressed in degenerate forms as follows:
M
(s) Σ
Om (1)
(41)
vill (1) Σ
Σ
tm + twy 1 Stnn,,,
neN, n'EM
18

---

Published as a conference paper at ICLR 2026
Substituting Equation (+II into Equation (G), we find that the solutions, (')(s,s))) admit de- generate forms as
h5( 7 Σ k(s. - -to) ~
Vals,, t)),t,,)
neN,
M
% 1ell
?
(s)
Σ
émls' (ools'-nn)lnppp -on)
Σ
w (t)(tt sf s')dt
(42)
m=1
LnEN
M
(m) (s'),
m=
where (6') ( (m..)eeDMMM are unknown coefficient functions. By substituting Equations (
and (II) into Equation (E), we obtain the following linear system that the coefficient functions must satisty:
1 Eem(s
(5))
Σ IEZ
/ (t)(t)
Eo (t)cme) (5))dt - Σ dml - to)s = 0, nEN,
2d (s) + ΣΣ dmee (8) S (((1)mm (1)dt
Σ
Φm (8 -o (43)
7
nEN,
for (m,j)= E [1 x u. Let us define the MU-dimensional stacked vector of coefficient functions
ē(s) (ct ...... [c (,, (1
(44)
Then, the linear system can be written compactly as
C Ju + E(s ()ss)
(45)
where x(s) and =\t are defined in Equations (155) and (IG), respectively. Substitting Equation into Equation (22) yields the solution the system of integral equations (6) as,
M h5(8,8)= = Σ dm (s)
m=1
+ss) T [c(s') 110--)MM
(46)
+ss) (=IMU + =)) *5(s))
140-1)M:MM
This completes the proof.
F PROOF OF PROPOSITION
Proof. Substituting Equation (44] into Equation (5) yields the expression for the estimated triggering
kernels in terms of the feature maps:
(1)(5) d(s)T (+ + +
=
Σ ((tn) Po 1-¹ ((t))t)]
(47)
neA,
1+(--))MMMM
19

---

Published as a conference paper at ICLR 2026
By using Equation (I=), the double integral in the denominator of Equation (l2) be rewritten using the feature maps as follows:
L.
Σ MEN
Nun (t in 8) s)1occttn
Φ(t tn) sa) ( IMU +8) -1 f-" ((s)ds =
Σ
REN
1+(un M3M M
/ /'' Φit tn) (n) "[ + 3 -1 Lo o(s)ds
leu "ENT
1-(1-1)M:M
EL (((t)dtt
[GIMO +0)
[¹" 2 o(s)ds)
leu
1+(1-1)MMM
= (1," ((0aaa)
(H
" x(0)))
(48)
Similarly, the integral term in the numerator of Equation becomes,
! hu,, (1 bus li Doctt
ΣE MEN WEN,
Σ nEN
Φ(1 1o) [
[TT IMU + =)
:(bxr))
n'EN;
1+(u, 11)M M
EV (i(O))))
INsU + +))
Σ
T(tw)))
leu
N'EN:
1+ 111
= -(. a(t)dt) ) (+)
=
Σ
S(t )ct
(49)
W'EN
This completes the proof.
Furthermore, from Equation (99) and the following identity:
['
Σ Σ BEN BEN n'EN,
Ne.,, (lwr In -
@ctn, - tn In) [G IMU + +3) +=)"'((" (ce)at) It(un 1)M8:a. M
Σ
Σ
EEN '''NN
[Grm=++ +3) +=) '1/" )(t)at)
Σ
Σ
Σ
p(tn' tu Lo) - Lo)
W'EN: mENi
1+0-1)MMMM
[+ +3) + ((t)att)
Σ di(tw)
Tell '
1+(1-1)MaM
(H +
((( o(t)dt)
Σ olt )dt
n'ENs
( ((())t)
(+ +3)
Σ
)(t:))dt)
n'EN
where the final equality holds because Im-'IMU + =) is symmetric, we obtain the relation in Equa- tion (39),, which holds for any M I oo and feature map Φ(s).
20

---

Published as a conference paper at ICLR 2026
Table 5: Average CPU time in seconds across 10 trials. denotes the average data size per trial.
Bonnet
Ours cpu
T
N
cpu
10000
8248
9250 26406
19.1
15000
12748
29.5
G EXISTENCE OF EQUIVALENT KERNELS
(+, :) Jee defined by the system of Fredholm inte- gral equations (6) exist and are uniquely determined for any positive semi cemi-definite kernel k(-,)).
We demonstrate that the equivalent kernels, {hy(-
A positive semi-definite kernel can be represented as the inner product of an M M-dimensional feature vector, (b.) (<h1 Φ (.)) for some M Dc. Under this representation, the equivalent
kernels, (hn jeu, can be expressed by Equation (14). The block matrix E = [.y] E TR. ML MMU
Equation (1El) has set submatrices defined in Equation (LS), which can be rewritten as follows:
min min(T..Att..Atte))
Bi = Σ
Σ NEN, n N'E,
Φ 4n)((t (x)) de
(emnart....ncmm........... A+t..)
max((ln)e))
(50)
/
Σ Σ NeN; R'E,,
(t))) - tn) tn)) Tdt. dt.
where (t) = 1,,, cr (...)) Therefore, for all c the following inequality holds:
RM
c =
[' (L)) mm) t )dl,
ΣΣΣ
Σ
Cim Gym' Σ
Σ
VEll JEUl m= m'=1
neN, R'EN
Cim Σ
(n(t) - tmx)
Σ LiEll m/=1
Cjmm E WeN,
(t)( two) dt (51)
Liel m=1
neV:
M Σ
Σ
Cim
(a(l)mml - kn)
di
LIEN m=1
nEN:
2 0.
This relation shows that ( T - in Equation (144 is positive definite and invertible for >> 0.
Therefore, the equivalent kernels defined by Equations (6) and (1E) exist and are uniquely deter- mined for any positive semi-definite kernels k(-,.).
H ADDITIONAL EXPERIMENTS
H.1 SCALABILITY ON LARGER DATA SIZE
In Section we discussed the scalability of the proposed method on the data size. To confirm it, we conducted experiment in the refractory scenario with T E (10000, 15000), and evaluated the CPU times of Ours and Bonnet on these larger datasets. The results Table (5 demonstrate that Ours remains scalable for the larger data sizes.
H.2 SCALABILITY ON LARGER DIMENSIONALITY
The computational cost our method (Ours) scales cubically with the dimensionality U, which disadvantage compared the quadratic scaling of the prior kernel method (Bonnet). We
21

---

Published as a conference paper at ICLR 2026
Table 6: Average CPU time in seconds across 5 trials.
Exp cpu
Gau
3er
Bonnet quu
Ours
U
qpu
cpu
cpu
3
124
7.03
10.8
413
5.04
1369
300
254
8513
29.9
Table CPU time in seconds of Ours regarding the dimensionality across 5 trials with standard emrors in brackets.
U 10
U 50
U = 100
U 200
=03
U 500
cpu
cpu
cpu
cpu
cpu
cpu
160
478
(0.26)
(0.82)
(0.67)
(1.91)
(6,70)
(15.5)
conducted experiment under refractory scenario with T = 2000 and U = 15 to examine this block-diagonal matrix obtained by
issue. The triggering kernel matrix was constructed as a U placing copies of the 3 x 3 triggering kernel matrix, g(s) [gin (s)]aj. used in Section 1.2 along the diagonal. We fixed the hyperparameters to (7,, (3) = (1, 1) and evaluated only the computation
time.
The results in Table & (note that the U = case is identical to that reported in Table for 2000) show that all methods exhibited increase in computation time U grows. However, the increase for Ours is more moderate compared to the conventional methods (Exp, Gau, Ber, and Bonnet). Although this trend may contradict the complexity analysis presented in Section 24L it can be attributed to the fact that our method relies solely on matrix additions and matrix inversions (performed via Cholesky decomposition), which are highly amenable to parallelization across multiple CPU cores.
10'
Training [sec]
10
10
10
10
10a
Dimensionality
Figure 3: Scaling curve of the training time regarding the dimensionality. The blue line represents the training time of Ours with respect to the dimensionality U. The black line represents a reference curve of (UF)).
We further conducted an experiment to empirically characterize how the computational cost of our method scales with respect to U. The computational burden of solving the linear system is most evident settings where the number of events per dimension is small while U itself is large. Accordingly, we set μ4 0.1. Str(.) (.) 0, and T = 100, and generated synthetic datasets for U E {10, 50, 100, 200, 300, 500). We then measured the training time of the proposed method. Table @ and Figure B summarize the results: empirically, the training time grows sub-cubically in U up to around ~ 100, and the scaling curve gradually approaches O At least up to U 500, the observed scaling therefore better than O(U3). We attribute this behavior largely the highly
22

---

Published as a conference paper at ICLR 2026
Table Estimation accuracy of triggering kernels and predictive performance by Ours. Results are on the mutually exciting scenario dataset, averaged over 10 trials, with standard errors in brackets. all is the negative log-like oogliikeihood on the test data. A denotes the integrated squared error for triggering kernel estimation, of which results are reproduction of those presented in Table III
T = 2000
T
3000
T = 5000
nll A²
rll
A*
rll
A
3294
0.38
3013
0.27
2842
0.20
(481)
(0.15)
(226)
(0.06)
(194)
(0.16)
optimized implementation of Cholesky factorization in TensorFlow, which efficiently exploits the 12-core CPU architecture.
Additionally, we implemented naïve conjugate gradient (CG) solver compare its runtime with the Cholesky factorization (CF) approach, using very simple preconditioner given by the inverse of the diagonal of ( (H Tare + 2). Unfortunately, this configuration led to worse runtime than CF, suggesting that a substantially more suitable preconditioner is required to fully benefit from CG. Designing such preconditioner is nontrivial and beyond the scope of this study. Based on these findings, our current model should be regarded as being practically targeted at event data with up to few hundred dimensions.
H.3 EFFECTS OF HYPERPARAMETER GRID ON PERFORMANCE
For the proposed model Ours, we conducted an experiment under the refractory scenario with T 5000, where the grid was refined from
to 10 10. The resulting squared error A was 0.58 + 0.12, which represents only marginal improvement over the 3 x 3 grid (0.59 + 0.13). This result suggests that the performance Tables μНД especially the gap between Bonnet and Ours, could not be solely attributed to the hyperparameter tuning strategy. Since Bonnet: is based on the likelihood function, is expected achieve higher accuracy than Ours, which relies on the least squares loss. Note that maximum likelihood estimation is known to be statistically efficient asymptotically for Hawkes processes (see (Ogata, (1789))
H.4 EFFECTS OF ESTIMATION ACCURACY OF TRIGGERING KERNELS ON PREDICTIVE PERFORMANCE
We conducted an additional experiment to examine the relationship between the estimation accuracy of the latent triggering kernel and the predictive performance of the point process model. Based on the mutually-exciting scenario dataset in Section I. we estimated the triggering kernel for each trial of T {2000, 3000, 5000} by using the proposed model (Ours), and evaluated its predictive performance on T = 7000 data. Here, predictive performance was assessed using the negative log- likelihood (the lower, the better). Since the T 7000 data composed of 10 trials, we report the average predictive performance over these trials. In addition, because each training dataset also contains 10 trials, we repeated the above prediction experiment 10 times. We applied post-hoc clipping, max(ā(t), 10 10= to the estimated intensity function Å(-). The results, summarized in Table & show that the predictive performance of the point process model indeed improves as the estimation accuracy of the latent triggering kernel increases.
I DEEP NEURAL NETWORK MODELS VS. KERNEL METHODS
Deep neural network (DNN) models currently dominate much of the machine learning literature. Nevertheless, approaches based on kernel methods remain highly valuable for several important reasons.
Reliability of optimization: Kernel method-based models (including proposed method) typ- ically reduce to convex optimization problems, which enjoy unique global optimum and stable training behavior. As result, they are subject to the reproducibility issues often observed in temporary DNNs, where the final performance can vary significantly depending on the optimization algorithm, its tuning parameters, and the parameter initialization. Moreover, kemel method-based
23

---

Published as a conference paper at ICLR 2026
models involve only few regularization parameters and a small number of kernel parameters, sulting in significantly fewer hyperparameters than DNNs. This makes the models substantially easier to tune in practice.
Model interpretability: Many DNN-based Hawkes processes (or more generally, history- dependent temporal point processes) embed the past event history into a latent vector via DNN and directly output the conditional intensity the next event (Zuo et all 2020; While this approach is powerful, especially for marked point processes with unstructured marks and complex history dependence, the influence of past events on future events is effectively a black box, making the event-generation dynamics difficult interpret or explain. contrast, both con- ventional kernel method-based models and our proposed method explicitly learn triggering kernels, which provide transparent description of how past events affect future intensity functions, thereby offering high interpretability.
20TG).
Complementarity of kernel methods and deep models: We also view unproductive to posi- tion kernel methods and DNNs fundamentally competing paradigms. Hybrid approaches such as deep kernel learning (Wilsson all 20TG,, which combine the computational advantages of kernel methods with the representational power of DNNs, are promising directions for future research. Our current work contributes to this broader agenda by strengthening the theoretical and methodological foundations on the kernel side, which can subsequently be integrated with deep architectures.
24


Reviews:
[{"review_id": "UqcUE86sg5", "reviewer": "Reviewer_EXZy", "summary": "This paper proposes a new representer theorem for linear multivariate Hawkes processes trained with penalized least squares in an RKHS. The theorem says that each triggering kernel can be written as a sum of equivalent kernels** that solve a system of Fredholm integral equations. A key outcome is that the dual coefficients are all 1, so there is no separate optimization over them. The authors show an efficient construction using random Fourier features, which reduces learning to one matrix inversion of size $MU \\times MU$ (independent of the number of events). On synthetic datasets, the method keeps similar accuracy to strong baselines while being much faster. The paper also discusses limits: linear link (no guaranteed non-negativity) and cubic cost in the process dimension $U$.", "questions": "1. For large \\(U\\), do **conjugate-gradient** solvers or block preconditioners help? Please provide scaling curves vs. \\(U\\).  \n2. How sensitive are results to the **support window $A$**? Can $A$ be learned ?  \n3. Consider moving the long proof to the **appendix**, and keeping only a short, readable proof sketch in the main text.", "limitations": "", "rating": 6, "confidence": 3, "soundness": 3, "presentation": 2, "contribution": 3, "strengths": "- **Originality:** First representer theorem for the non-discretized penalized LS problem in linear Hawkes; fixing all dual coefficients to 1 removes a costly optimization step.  \n- **Quality:** Clear path from theorem to algorithm using degenerate kernels / random Fourier features, with closed-form estimators.  \n- **Clarity:** Assumptions and complexity are explicit; the narrative from theory to implementation is easy to follow.  \n- **Significance:** Very large speed-ups compared to strong baselines while keeping similar accuracy for larger datasets.", "weaknesses": "- **Evidence scope:** Only synthetic data are used; real-world robustness and data issues are not tested.  \n- **Proof placement:** The main-text sketch uses inverse operators and delta distributions and may feel less rigorous; the appendix version is better suited for details.  \n- **Modeling limit:** Linear link does not guarantee non-negativity of intensity; post-hoc clipping or other fixes may be needed.  \n- **Scalability in $U$:** Cost scales as $O(N^2 M^2 U^2 + M^3 U^3)$; very high-dimensional processes or very long sequences can still be heavy.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "uYsBdhyDD8", "reviewer": "Reviewer_MBXV", "summary": "This paper proposes a new representer-theorem–based method for estimating triggering kernels of linear multivariate Hawkes processes in an RKHS framework. By formulating a system of simultaneous integral equations, the optimal estimator is expressed as a linear combination of transformed kernels evaluated at data points, with all dual coefficients analytically fixed to one. As a result, the method avoids solving high-dimensional optimization problems and achieves substantial computational efficiency gains. Experiments on synthetic data show that the approach matches state-of-the-art predictive accuracy while being significantly faster than existing nonparametric kernel estimators.", "questions": "1. The analytical framework appears to rely on classical techniques, and at first glance the contribution seems incremental. Could you clarify precisely where the conceptual or technical novelty lies in your representer-theorem formulation for Hawkes processes?\n2. In the era of deep neural networks, kernel methods are sometimes viewed as outdated or less competitive. Why is a kernel-based approach appropriate and timely for ICLR, and what advantages does it offer over contemporary deep-learning-based models for point processes?\n3. The representer theorem is a well-established and classical concept in kernel methods. How does the representer theorem presented in this paper differ from existing versions, including variants used in prior work on point processes? Is the contribution more than a minor extension, and in what sense is it a meaningful and substantial generalization?", "limitations": "", "rating": 2, "confidence": 3, "soundness": 2, "presentation": 2, "contribution": 2, "strengths": "This paper establishes the first representer theorem for the non-approximated penalized least squares formulation of linear multivariate Hawkes processes, filling a clear theoretical gap in kernel-based point process modeling. The resulting estimator has all dual coefficients analytically fixed to one, eliminating the need for high-dimensional nonlinear optimization and dramatically improving scalability over prior methods. The proposed approach yields a closed-form solution of the integral equations via random feature approximations, avoiding Riemann-sum discretization used in earlier work and enabling efficient computation. The final estimator consists only of additive matrix operations and a matrix inversion whose size is independent of the data size, making the method lightweight, scalable, and well-suited for large-scale Hawkes process datasets. Overall, the paper offers both theoretical novelty and practical computational advantages, advancing kernel-based learning for multivariate Hawkes processes.", "weaknesses": "Although the proposed method offers strong computational advantages, it relies on the identity link function and therefore applies only to linear Hawkes processes; models requiring non-linear link functions to enforce non-negativity or capture inhibitory effects fall outside its scope. In addition, the representer theorem is developed under a penalized least-squares formulation rather than maximum likelihood, which may limit statistical efficiency in some settings. The approach also requires solving Fredholm integral equations whose accuracy depends on the random feature approximation, potentially introducing approximation error relative to exact kernel evaluations. Finally, the empirical evaluation is conducted on synthetic datasets, leaving open questions about robustness and performance on real-world event data.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "kNRGzdEy00", "reviewer": "Reviewer_wcyH", "summary": "The paper aims to provide a theoretically principled methodology for reducing the computational complexity of nonparametric estimation of the intensity functions in multivariate linear Hawkes process, where the triggering \"kernels\" $g_{ij}$ are assumed to belong to an RKHS $\\mathcal{H}_{k}$ with reproducing kernel $k$ and a Tikhonov regularization is employed. \n\nThe novel technical contribution of this paper and the workhorse behind this paper is a representer theorem for the sample based estimates $\\hat{g}_{ij}$ in terms of  linear combinations of Fredholm integral transformed versions of the base kernel function $k$, denoted by $h_j$'s. The representer theorem ensures that there is no need for dual variable (coefficients in terms of $h_ij$'s) optimization, leaving only the system of Fredholm integral equations for $h_j$'s to be solved. Using the machinery of Random Fourier Features (RFF), a practical algorithm is proposed with corresponding experiments based on synthetic data that achieves $O(N^{2}M^{2}U^{2} +M^{3}U^{3})$ computational complexity where $N$ is the maximum sampling frequency in each dimension, M is the number of Random features and U is the dimensionality of the multivariate linear Hawkes process, improving upon $O(N^{4}U^{2})$ computational complexity of the SOTA method (Bonnet and Sangnier, 2025), along with the usual memory complexity improvement expected from RFF. While computation time on synthetic experiments drops by around 2 orders of magnitude, inferential accuracy loss is reasonably low.", "questions": "Based on the discussion of the weaknesses, the paper will improve if\n1.  the authors can report the specific conditions (properties of kernels etc.) under which the Fredholm integral equations for $h_j$'s (Equation 6 under Theorem 1) admit a solution and \n2. demonstrate the implication of the small loss of inferential power for the latent trigerring \"kernels\" $g_{ij}$'s for actual predictive performance for the point process. \n\nAlthough this did not play a role in my assessment, but, as always, a real-life experimment will always help, especially regarding the evaluation of predictive performance.", "limitations": "", "rating": 8, "confidence": 4, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "The paper presents a very clear and precise solution to a significant computational problem in the context of linear Hawkes processes,and explains the utility behind analyzing the continuous variational formulation of the data fitting problem, along with a clean and complete derivation of the key technical result being included. The representer theorem is novel in the context of this particular problem, and shows the theoretical advantage of using the Tikhonov regualrization of RKHS norm in contrast to the Hakwes process log-likelihood based estimation approach. Synthetic experiments also show around 100 times relative improvement in computational complexity compared to SOTA with small loss of inferential power for the trigerring \"kernels\" $g$, in terms of integrated squared error.", "weaknesses": "While the RFF approach indeed allows a practical solution to the problem, from a theoretical perspective it would be helpful if the authors can report the specific conditions (properties of kernels etc.) under which the Fredholm integral equations for $h_j$'s (Equation 6 under Theorem 1) admit a solution. The authors themselves acknowledge the limitation of the paper being restricted to linear Hawkes processes and requires adhoc post-processing in empirical experiments, while approach of the SOTA method is more general and is applicable to nonlinear processes as well. Also, from a empirical perspective, the implication of the small loss of inferential power for the latent trigerring \"kernels\" $g_{ij}$'s for actual predictive performance for the point process is not immediately clear.\n\nMinor typos:\n1. Replace $x$'s by $s$ in Equation 22 (Lines 452 and 453)\n2. Replace \"quantifing\" by quantifying on Line 096\n3. Replace \"evet\" by event on Ine 134", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "LePWBRBXPg", "reviewer": "Reviewer_JiiZ", "summary": "The authors consider nonparametric estimation of the triggering kernels in a multivariate linear Hawkes process using an RKHS (reproducing kernel Hilbert space) framework. They formulate a penalised least-squares loss for the latent triggering functions, and derive a novel representer theorem showing that the optimal kernel estimators admits an expansion over a family of transformed kernels evaluated at data points, and remarkably, the dual coefficients are analytically fixed to unity, so one does not need to solve an optimization problem. These transformed kernels are defined by a system of Fredholm integral equations, which allows an efficient random-feature approximation algorithm, further enabling scalability since number of random features can be typically set much smaller than number of samples. \n\nThe method applies to linear multivariate Hawkes (identity link function, no inhibitory effects or non-linear link). So the generality is narrower than the Bonnet & Sangnier (2025) paper which handles inhibition and non-linear link functions. The trade-off is computational simplicity and analytic coefficient expressions.\n\nEmpirical results on synthetic toy-models show competitive accuracy and much better efficiency compared to existing kernel-based methods.", "questions": "Minor issues:\n- line 278 it should read \"number\"\n- line 132, definition of ${\\cal N}_i$ is missing, though intuivily clear what it should represent\n- line 147, would be good to explicitly write what dual coefficients equal to\n- line 155 personally I would prefer using $\\dagger$ for pseudo inverse sainted of $*$ (typically adjoint)", "limitations": "", "rating": 6, "confidence": 3, "soundness": 3, "presentation": 3, "contribution": 2, "strengths": "(1) Novel representer theorem: While prior works (Bonnet & Sangnier 2025) derived a general representer theorem for approximate loss formulations (and requiring optimization of dual coefficients), this paper focuses on linear Hawkes processes and provides one for the exact penalised least squares loss obtaining that dual coefficients equal unity w.r.t. transformed kernels characterised by Fredholm integral equations  That is a non-trivial theoretical advance.\n\n(2) Practical algorithm: To obtain practical algorithm, authors use random features to avoid discretization approximation of the integral operators,  and, in process, reduce the inversion to a matrix independent of data size. This results in practical algorithm that scales much better than previous kernel-based methods for Hawkes. Further, authors empirically show that you can maintain competitive predictive accuracy while significantly improving computational efficiency. This allows algorithm's deployment in large event-data applications.\n\n(3) Paper is generally well written and related works are, up to my knowledge, well covered.", "weaknesses": "(1) Limited depth: The paper is fundamentally methodological. New kernel method is backed by the representer theorem and coupled with random features to result in practical algorithm. However, any form of learning guarantees is lacking. \n\n(2) Scalability via random features: Using random features and reducing the inversion to a matrix independent of data size is golden standard in scaling kernel methods. Hence, there is hardly some big novelty in this. Even more so, since there is lack on any theory such as e.g. convergence w.r.t. number of random features.\n\n(3) Empirical evaluation: Since the paper is methodological,  testing it only on synthetic toy-models is a bit underwhelming. \n\nTo conclude, up to my knowledge, the paper delivers a clear meaningful advance: the theoretical representer result is strong (analytic coefficients) and the algorithmic scalability is obtained via random features. As such, it can make a solid contribution to kernel-based nonparametric estimation for Hawkes processes if the empirical study is improved, or some theoretical  learning guarantees included.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}]

Author rebuttal:

Decision: Accept (Oral)