Functional central limit theorems for stationary Hawkes processes and application to infinite-server queues

Abstract

A univariate Hawkes process is a simple point process that is self-exciting and has a clustering effect. The intensity of this point process is given by the sum of a baseline intensity and another term that depends on the entire past history of the point process. Hawkes processes have wide applications in finance, neuroscience, social networks, criminology, seismology, and many other fields. In this paper, we prove a functional central limit theorem for stationary Hawkes processes in the asymptotic regime where the baseline intensity is large. The limit is a non-Markovian Gaussian process with dependent increments. We use the resulting approximation to study an infinite-server queue with high-volume Hawkes traffic. We show that the queue length process can be approximated by a Gaussian process, for which we compute explicitly the covariance function and the steady-state distribution. We also extend our results to multivariate stationary Hawkes processes and establish limit theorems for infinite-server queues with multivariate Hawkes traffic.

Appendices

Proofs of results in Section 3

This section collects the proofs of results in Sect. 3.

1.1 Proof of Theorem 2

Proof of Theorem 2

The proof relies on Hahn’s theorem (see Theorem 2 in [26] or Theorem 7.2.1. in [53]), and delicate estimates of moments of stationary Hawkes processes.

For the sake of simplicity, we first let $\mu$ be a positive integer. By the immigration-birth representation, we can decompose $N^{\mu }$ as the sum of $\mu$ independent and identically distributed (i.i.d.) Hawkes processes $N_i^1, i =1, 2,\ldots , \mu$, each distributed as a stationary Hawkes process with baseline intensity 1 (the superscript 1 in $N_i^1$) and the exciting function $h(\cdot )$. For notational simplicity, we use $N_i(\cdot )$ for $N_i^1 (\cdot )$. Therefore, we have

$$\begin{aligned} \frac{N^{\mu }(t)-{{\bar{\lambda }}} t}{\sqrt{\mu }} =\frac{1}{\sqrt{\mu }}\sum _{i=1}^{\mu }\left[ N_i ({t})-\frac{t}{1-\Vert h\Vert _{L^{1}}}\right] . \end{aligned}$$

Let ${\tilde{N}}_i ({t}):=N_i({t})-\frac{t}{1-\Vert h\Vert _{L^{1}}}$. Then, ${\tilde{N}}_i$ are i.i.d. random elements of $D([0,\infty ),{\mathbb {R}})$ with ${\mathbb {E}}[{\tilde{N}}_i(t)]=0$ (see, for example, Eq. (9) in [27]) and ${\mathbb {E}}[({\tilde{N}}_i (t))^{2}]<\infty$ for any t (see, for example, Lemma 2 in [58]).²

By Hahn’s theorem, since ${\tilde{N}}_{i}$ are i.i.d., as $\mu \rightarrow \infty ,$ we have

$$\begin{aligned} \frac{1}{\sqrt{\mu }}\sum _{i=1}^{\mu }\left[ N_i ({t})-\frac{t}{1-\Vert h\Vert _{L^{1}}}\right] =\frac{1}{\sqrt{\mu }}\sum _{i=1}^{\mu }{\tilde{N}}_{i}(t) \Rightarrow G(t), \end{aligned}$$

weakly in $(D([0,\infty ),{\mathbb {R}}),J_{1})$, where G is a mean-zero almost surely continuous Gaussian process with the covariance function of ${\tilde{N}}_1$, provided that the following condition is satisfied: For every $0<T<\infty$, there exist continuous non-decreasing real-valued functions g and f on [0, T] and numbers $\alpha >1/2$ and $\beta >1$ such that

$$\begin{aligned} {\mathbb {E}}\left[ \left( {\tilde{N}}_1 ({u}) -{\tilde{N}}_1 ({s}) \right) ^{2}\right] \le (g(u)-g(s))^{\alpha }, \end{aligned}$$

(A.1)

and

$$\begin{aligned} {\mathbb {E}}\left[ \left( {\tilde{N}}_1 ({u}) - {\tilde{N}}_1 ({t})\right) ^{2} \left( {\tilde{N}}_1 ({t}) -{\tilde{N}}_1 ({s})\right) ^{2}\right] \le (f(u)-f(s))^{\beta }, \end{aligned}$$

(A.2)

for all $0\le s\le t\le u\le T$ with $u-s<1$.

Let us prove (A.1) and (A.2). For notational simplicity, we use $N_1(a,b]$ to stand for $N_1((a,b])$ (equivalently, $N_1(b)-N_1(a)$), which records the number of points of the process $N_1$ in the interval (a, b]. We also use $\lambda _1$ to denote the intensity process of the stationary Hawkes process $N_1$ with baseline intensity 1. We now present a lemma which is the key to the proofs of (A.1) and (A.2). The proof is given at the end of this section.

Lemma 15

We have

$$\begin{aligned} {\mathbb {E}}[(\lambda _1(0))^4] < \infty . \end{aligned}$$

(A.3)

As a result, for all $0\le s\le t\le u\le T$ and $u-s<1$, there are some constants $c, C>0$ independent of s, t, u,  such that

$$\begin{aligned} {\mathbb {E}}\left[ (N_1(s,u])^{2}\right]\le & {} C \cdot (u-s), \end{aligned}$$

(A.4)

$$\begin{aligned} {\mathbb {E}}\left[ (N_1(t,u])^{2}(N_1(s,t])^{2}\right]\le & {} c \cdot (u-s)^2, \end{aligned}$$

(A.5)

With Lemma 15, we are ready to prove (A.1) and (A.2). First, let us prove (A.1). It is clear from the definition of ${{\tilde{N}}}_1$ that

$$\begin{aligned} {\mathbb {E}}\left[ \left( {\tilde{N}}_1 ({u}) -{\tilde{N}}_1 ({s}) \right) ^{2}\right]&={\mathbb {E}}\left[ \left( N_1(s,u]-\frac{u-s}{1-\Vert h\Vert _{L^{1}}}\right) ^{2}\right] \\&\le {\mathbb {E}}\left[ \left( N_1(s,u]\right) ^{2}\right] +\left( \frac{u-s}{1-\Vert h\Vert _{L^{1}}}\right) ^{2}. \end{aligned}$$

Using (A.4) in Lemma 15 and the fact that $0 \le u-s<1$, we immediately obtain that (A.1) is satisfied with $g(x)=\left( C + \frac{1}{(1- \Vert h\Vert _{L^{1}})^2} \right) x$ and $\alpha =1.$

Next, let us prove (A.2). Note that

$$\begin{aligned}&{\mathbb {E}}\left[ \left( {\tilde{N}}_1 ({u}) - {\tilde{N}}_1 ({t})\right) ^{2} \left( {\tilde{N}}_1 ({t}) -{\tilde{N}}_1 ({s})\right) ^{2}\right] \nonumber \\&\quad ={\mathbb {E}}\left[ \left( N_1(t,u]-\frac{u-t}{1-\Vert h\Vert _{L^{1}}}\right) ^{2} \left( N_1[s,t]-\frac{t-s}{1-\Vert h\Vert _{L^{1}}}\right) ^{2}\right] \nonumber \\&\quad \le {\mathbb {E}}\left[ \left( (N_1 (t,u])^{2}+\left( \frac{u-t}{1-\Vert h\Vert _{L^{1}}}\right) ^{2}\right) \left( (N_1 (s,t])^{2}+\left( \frac{t-s}{1-\Vert h\Vert _{L^{1}}}\right) ^{2}\right) \right] \nonumber \\&\quad = \left( \frac{u-t}{1-\Vert h\Vert _{L^{1}}}\right) ^{2}{\mathbb {E}}\left[ (N_1(s,t])^{2}\right] +\left( \frac{t-s}{1-\Vert h\Vert _{L^{1}}}\right) ^{2}{\mathbb {E}}\left[ (N_1(t,u])^{2}\right] \nonumber \\&\qquad +\left( \frac{u-t}{1-\Vert h\Vert _{L^{1}}}\right) ^{2}\left( \frac{t-s}{1-\Vert h\Vert _{L^{1}}}\right) ^{2} +{\mathbb {E}}\left[ (N_1 (s,t])^{2}(N_1(t,u])^{2}\right] . \end{aligned}$$

(A.6)

Since $0\le s\le t\le u\le T$ and $u-s<1$, we can then infer from Lemma 15 and (A.6) that (A.2) is satisfied with $f(x)=C'x$ for some positive constant $C'$ (independent of u, s, t) and $\beta =2$.

Now we have proved Theorem 2 by assuming $\mu$ is a positive integer in our discussions. The same result holds when $\mu \in (0, \infty )$ for $\mu \rightarrow \infty$. Note that for $\mu \in (0, \infty )$, by the immigration-birth representation, the process $N^{\mu }(\cdot )$ can be decomposed as the sum of two independent stationary Hawkes processes $N^{\lfloor \mu \rfloor }(\cdot )$ and $N^{\mu -\lfloor \mu \rfloor }(\cdot )$, where the superscripts $\lfloor \mu \rfloor , \mu -\lfloor \mu \rfloor$ represent the baseline intensities, respectively. Hence, to show the result holds for $\mu \in (0, \infty )$ with $\mu \rightarrow \infty$, it suffices to show that for any $T>0$,

$$\begin{aligned} \sup _{0\le t\le T}\frac{N^{\mu -\lfloor \mu \rfloor }(t)-\frac{(\mu -\lfloor \mu \rfloor )t}{1-\Vert h\Vert _{L^{1}}}}{\sqrt{\mu }}\rightarrow 0, \end{aligned}$$

in probability as $\mu \rightarrow \infty$. This can be easily verified since, for any $\epsilon >0$, for sufficiently large $\mu$, we have $\sup _{0\le t\le T}\frac{1}{\sqrt{\mu }} \frac{(\mu -\lfloor \mu \rfloor )t}{1-\Vert h\Vert _{L^{1}}} \le \frac{1}{\sqrt{\mu }}\frac{T}{1-\Vert h\Vert _{L^{1}}}\le \frac{\epsilon }{2}$, and

$$\begin{aligned} {\mathbb {P}}\left( \left| \sup _{0\le t\le T}\frac{N^{\mu -\lfloor \mu \rfloor }(t)-\frac{(\mu -\lfloor \mu \rfloor )t}{1-\Vert h\Vert _{L^{1}}}}{\sqrt{\mu }}\right| \ge \epsilon \right)&\le {\mathbb {P}}\left( \left| \sup _{0\le t\le T}\frac{N^{\mu -\lfloor \mu \rfloor }(t)}{\sqrt{\mu }}\right| \ge \frac{\epsilon }{2}\right) \nonumber \\&={\mathbb {P}}\left( N^{\mu -\lfloor \mu \rfloor }(T)\ge \frac{1}{2}\sqrt{\mu }\epsilon \right) \nonumber \\&\le {\mathbb {P}}\left( N^{1}(T)\ge \frac{1}{2}\sqrt{\mu }\epsilon \right) \\&\le \frac{2}{\sqrt{\mu }\epsilon }T \cdot {\mathbb {E}}[N^{1}(1)]\rightarrow 0, \end{aligned}$$

as $\mu \rightarrow \infty$. Here $N^1$ denotes a stationary Hawkes process with a baseline intensity one and an exciting function h.

Finally, let us compute the covariance function of G, or equivalently (from Hahn’s Theorem), the covariance function of ${\tilde{N}}_1(\cdot )$. Since ${\tilde{N}}_1 (t)$ and ${\tilde{N}}_1 (s)$ have mean zero, we can compute that, for any $t>s$,

$$\begin{aligned} \text{ Cov }({\tilde{N}}_1(t),{\tilde{N}}_1(s))&={\mathbb {E}}\left[ {\tilde{N}}_1(t) {\tilde{N}}_1(s)\right] \nonumber \\&={\mathbb {E}}\left[ \left( N_1(t)-\frac{t}{1-\Vert h\Vert _{L^{1}}}\right) \left( N_1(s)-\frac{s}{1-\Vert h\Vert _{L^{1}}}\right) \right] \nonumber \\&={\mathbb {E}}[N_1 (t) N_1(s)] -\frac{ts}{(1-\Vert h\Vert _{L^{1}})^{2}} \nonumber \\&={\mathbb {E}}[(N_1(t)-N_1(s)) N_1(s)] +{\mathbb {E}}[(N_1(s))^{2}]-\frac{ts}{(1-\Vert h\Vert _{L^{1}})^{2}}. \end{aligned}$$

(A.7)

It is clear that

$$\begin{aligned} {\mathbb {E}}[(N_1(s))^{2}] =\text{ Var }(N_1(s))+ \frac{s^{2}}{(1-\Vert h\Vert _{L^{1}})^{2}} =K(s)+\frac{s^{2}}{(1-\Vert h\Vert _{L^{1}})^{2}}. \end{aligned}$$

In addition, we can verify that

$$\begin{aligned} {\mathbb {E}}[(N_1(t)-N_1(s))N_1(s)]= & {} {\mathbb {E}}\left[ \int _{s+}^{t}N_1(\mathrm{d}u)\int _{0}^{s}N_1(\mathrm{d}v)\right] \nonumber \\= & {} \int _{s+}^{t}\int _{0}^{s}{\mathbb {E}}[N_1(\mathrm{d}v) N_1(\mathrm{d}u)] \nonumber \\= & {} \int _{s}^{t}\int _{0}^{s}\left[ \phi (u-v)+\frac{1}{(1-\Vert h\Vert _{L^{1}})^{2}}\right] \mathrm{d}v\mathrm{d}u \nonumber \\= & {} \int _{s}^{t}\int _{0}^{s}\phi (u-v)\mathrm{d}v\mathrm{d}u+\frac{s(t-s)}{(1-\Vert h\Vert _{L^{1}})^{2}}. \end{aligned}$$

Hence, we get

$$\begin{aligned}&\text{ Cov }(G(t),G(s)) = \text{ Cov }({\tilde{N}}_1(t),{\tilde{N}}_1(s))\nonumber \\&\quad =\int _{s}^{t}\int _{0}^{s}\phi (u-v)\mathrm{d}u\mathrm{d}v+\frac{s(t-s)}{(1-\Vert h\Vert _{L^{1}})^{2}} +K(s)\nonumber \\&\qquad +\frac{s^{2}}{(1-\Vert h\Vert _{L^{1}})^{2}}-\frac{ts}{(1-\Vert h\Vert _{L^{1}})^{2}} \nonumber \\&\quad =\int _{s}^{t}\int _{0}^{s}\phi (u-v)\mathrm{d}u\mathrm{d}v+K(s). \end{aligned}$$

(A.8)

The proof is therefore complete. $\square$

Proof of Lemma 15

We first prove (A.3). Using the definition of the intensity $\lambda _1$ and the simple inequality $\left( \frac{x+y}{2}\right) ^4 \le \frac{x^4 + y^4}{2}$, we obtain that, for sufficiently small $\delta >0$,

$$\begin{aligned}&{\mathbb {E}}\left[ (\lambda _{1}(0))^{4}\right] ={\mathbb {E}}\left[ \left( 1+\int _{-\infty }^{0}h(-s)N_1(\mathrm{d}s)\right) ^{4}\right] \nonumber \\&\quad \le 8+8{\mathbb {E}}\left[ \left( \int _{-\infty }^{0}h(-s)N_1(\mathrm{d}s)\right) ^{4}\right] \nonumber \\&\quad \le 8+8{\mathbb {E}}\left[ \left( \sum _{i=0}^{\infty }\max _{t\in [-(i+1)\delta ,-i\delta ]}h(-t)\cdot N_1[-(i+1)\delta ,-i\delta ]\right) ^{4}\right] \nonumber \\&\quad = 8+8\left( \sum _{i=0}^{\infty }\max _{t\in [-(i+1)\delta ,-i\delta ]}h(-t)\right) ^{4} \nonumber \\&\qquad \times {\mathbb {E}}\left[ \left( \sum _{i=0}^{\infty }\frac{\max _{t\in [-(i+1)\delta ,-i\delta ]}h(t)}{\sum _{i=0}^{\infty }\max _{t\in [-(i+1)\delta ,-i\delta ]}h(t)}N_1[-(i+1)\delta ,-i\delta ]\right) ^{4}\right] . \qquad \quad \end{aligned}$$

(A.9)

Note here that, under Assumption 1, we know $h(\cdot )$ is locally bounded and Riemann integrable, hence, for sufficiently small $\delta >0$,

$$\begin{aligned} \sum _{i=0}^{\infty }\max _{t\in [-(i+1)\delta ,-i\delta ]}h(-t) < \infty . \end{aligned}$$

(A.10)

Applying the Jensen’s inequality to (A.9), we get

$$\begin{aligned}&{\mathbb {E}}\left[ (\lambda _{1}(0))^{4}\right] \nonumber \\&\quad \le 8+8\left( \sum _{i=0}^{\infty }\max _{t\in [-(i+1)\delta ,-i\delta ]}h(-t)\right) ^{4} \nonumber \\&\qquad \times {\mathbb {E}}\left[ \sum _{i=0}^{\infty }\frac{\max _{t\in [-(i+1)\delta ,-i\delta ]}h(t)}{\sum _{i=0}^{\infty }\max _{t\in [-(i+1)\delta ,-i\delta ]}h(t)}\left( N_1[-(i+1)\delta ,-i\delta ]\right) ^{4}\right] \nonumber \\&\quad = 8+8\left( \sum _{i=0}^{\infty }\max _{t\in [-(i+1)\delta ,-i\delta ]}h(-t)\right) ^{4} {\mathbb {E}}\left[ \left( N_1[0,\delta ]\right) ^{4}\right] <\infty , \end{aligned}$$

(A.11)

where we have used the stationarity of $N_1$, (A.10) and the fact that ${\mathbb {E}}[e^{\theta N_1[0,1]}]<\infty$ for sufficiently small $\theta >0$; see, for example,

[58].

We next prove (A.4). We can directly compute that

$$\begin{aligned} {\mathbb {E}}\left[ \left( N_1(s,u]\right) ^{2}\right]&\le 2{\mathbb {E}}\left[ \left( N_{1}(s,u]-\int _{s}^{u}\lambda _{1}(v)\mathrm{d}v\right) ^{2}\right] +2{\mathbb {E}}\left[ \left( \int _{s}^{u}\lambda _{1}(v)\mathrm{d}v\right) ^{2}\right] \nonumber \\&= 2{\mathbb {E}}\left[ \int _{s}^{u}\lambda _{1}(v)\mathrm{d}v\right] +2{\mathbb {E}}\left[ \left( \int _{s}^{u}\lambda _{1}(v)\mathrm{d}v\right) ^{2}\right] \nonumber \\&=2\frac{(u-s)}{1-\Vert h\Vert _{L^{1}}} +2{\mathbb {E}}\left[ \left( \int _{s}^{u}\lambda _{1}(v)\mathrm{d}v\right) ^{2}\right] \nonumber \\&\le 2\frac{(u-s)}{1-\Vert h\Vert _{L^{1}}} +2(u-s){\mathbb {E}}\left[ \left( \int _{s}^{u}(\lambda _{1}(v))^{2}\mathrm{d}v\right) \right] \nonumber \\&=2\frac{(u-s)}{1-\Vert h\Vert _{L^{1}}} +2(u-s)^{2}{\mathbb {E}}[(\lambda _{1}(0))^{2}]\nonumber \\&\le C (u-s), \end{aligned}$$

(A.12)

for some positive constant C. Here, the second line follows from the martingale property, see Sect. 2; the third line follows from the stationarity of the intensity process and ${\mathbb {E}}[\lambda _{1}(0)]=\frac{1}{1-\Vert h\Vert _{L^{1}}}$; the fourth line follows from the Cauchy–Schwarz inequality, so that $\left( \int _{s}^{u}\lambda _1(v)\mathrm{d}v\right) ^{2}\le \int _{s}^{u}1\mathrm{d}v\cdot \int _{s}^{u}(\lambda _{1}(v))^{2}\mathrm{d}v$; the fifth line is due to the fact that $\lambda _1$ is a stationary process; and the last inequality is due to (A.3) and the fact that $0<u-s<1$.

Finally, we prove (A.5). We can directly compute that

$$\begin{aligned}&{\mathbb {E}}\left[ (N_1(t,u])^{2}(N_1(s,t])^{2}\right] \nonumber \\&\quad \le 2{\mathbb {E}}\left[ \left( N_1(t,u]-\int _{t}^{u}\lambda _1 (v)\mathrm{d}v\right) ^{2}(N_1(s,t])^{2}\right] \nonumber \\&\qquad +2{\mathbb {E}}\left[ \left( \int _{t}^{u}\lambda _1(v)\mathrm{d}v\right) ^{2}(N_1(s,t])^{2}\right] \nonumber \\&\quad \le 2{\mathbb {E}}\left[ \int _{t}^{u}\lambda _1 (v)\mathrm{d}v\cdot (N_1 (s,t])^{2}\right] +2(u-t){\mathbb {E}}\left[ \int _{t}^{u}(\lambda _1 (v))^{2}\mathrm{d}v\cdot (N_1 (s,t])^{2}\right] , \end{aligned}$$

(A.13)

where we use the fact that $N_1(t,u]-\int _{t}^{u}\lambda _1 (v)\mathrm{d}v$ is ${\mathcal {F}}_{t}$-measurable and is a martingale with predictable quadratic variation $\int _{t}^{u}\lambda _1 (v)\mathrm{d}v$, so that $\left( N_1(t,u]-\int _{t}^{u}\lambda _1 (v)\mathrm{d}v\right) ^{2}-\int _{t}^{u}\lambda _1 (v)\mathrm{d}v$ is also a martingale (see Sect. 2), and the Cauchy–Schwarz inequality, so that $\left( \int _{t}^{u}\lambda _1(v)\mathrm{d}v\right) ^{2}\le \int _{t}^{u}1\mathrm{d}v\cdot \int _{t}^{u}(\lambda _{1}(v))^{2}\mathrm{d}v$. From here, we can further estimate that

$$\begin{aligned}&{\mathbb {E}}\left[ (N_1(t,u])^{2}(N_1(s,t])^{2}\right] \nonumber \\&\quad \le 2\left( {\mathbb {E}}\left[ \left( \int _{t}^{u}\lambda _1 (v) \mathrm{d}v\right) ^{2}\right] \right) ^{1/2} \left( {\mathbb {E}}\left[ (N_1(s,t])^{4}\right] \right) ^{1/2} \nonumber \\&\qquad +(u-t){\mathbb {E}}\left[ \int _{t}^{u}\left( (\lambda _1 (v))^{4}+(N_1 (s,t])^{4}\right) \mathrm{d}v\right] , \end{aligned}$$

(A.14)

where we applied the Cauchy–Schwarz inequality to the first term in the last line in (A.13), so that ${\mathbb {E}}\left[ \int _{t}^{u}\lambda _1 (v)\mathrm{d}v\cdot (N_1 (s,t])^{2}\right] \le \left( {\mathbb {E}}\left[ \left( \int _{t}^{u}\lambda _1 (v) \mathrm{d}v\right) ^{2}\right] \right) ^{1/2} \left( {\mathbb {E}}\left[ (N_1(s,t])^{4}\right] \right) ^{1/2}$, and the simple inequality $2(\lambda _1 (v))^{2}(N_1 (s,t])^{2}\le (\lambda _1 (v))^{4}+(N_1 (s,t])^{4}$ to the second term in the last line in (A.13). From here, we can continue that

$$\begin{aligned}&{\mathbb {E}}\left[ (N_1(t,u])^{2}(N_1(s,t])^{2}\right] \nonumber \\&\quad \le 2(u-t)^{1/2}\left( \left[ \int _{t}^{u}{\mathbb {E}}[(\lambda _1(v))^{2}]\mathrm{d}v\right] \right) ^{1/2} \left( {\mathbb {E}}\left[ (N_1 (s,t])^{4}\right] \right) ^{1/2} \nonumber \\&\qquad +(u-t)\left[ \int _{t}^{u}\left( {\mathbb {E}}(\lambda _1(v))^{4}+{\mathbb {E}}(N_1 (s,t])^{4}\right) \mathrm{d}v\right] \nonumber \\&\quad = 2(u-t)\left( {\mathbb {E}}[(\lambda _1(0))^{2}]\right) ^{1/2} \left( {\mathbb {E}}\left[ (N_1(s,t])^{4}\right] \right) ^{1/2} \nonumber \\&\qquad +(u-t)^{2}{\mathbb {E}}\left[ (\lambda _1(0))^{4}\right] +(u-t)^{2}{\mathbb {E}}\left[ (N_1(s,t])^{4}\right] , \end{aligned}$$

(A.15)

where we used the Cauchy–Schwartz inequality $\left( \int _{t}^{u}\lambda _1(v)\mathrm{d}v\right) ^{2}\le \int _{t}^{u}1\mathrm{d}v\cdot \int _{t}^{u}(\lambda _{1}(v))^{2}\mathrm{d}v$, and the stationarity of the Hawkes processes; see Sect. 2.

Hence, to bound ${\mathbb {E}}\left[ (N_1(t,u])^{2}(N_1(s,t])^{2}\right]$, we need to estimate ${\mathbb {E}}\left[ (N_1 (s,t])^{4}\right]$. We can compute that (explanations follow below)

$$\begin{aligned} {\mathbb {E}}\left[ (N_1(s,t])^{4}\right]&\le 8{\mathbb {E}}\left[ \left( N_1 (s,t]-\int _{s}^{t}\lambda _1 (v) \mathrm{d}v\right) ^{4}\right] +8{\mathbb {E}}\left[ \left( \int _{s}^{t}\lambda _1(v)\mathrm{d}v\right) ^{4}\right] \nonumber \\&\le 8{\bar{C}}{\mathbb {E}}\left[ \left( \int _{s}^{t}\lambda _1(v)\mathrm{d}v\right) ^{2}\right] +8{\mathbb {E}}\left[ \left( \int _{s}^{t}\lambda _1(v) \mathrm{d}v\right) ^{4}\right] \nonumber \\&\le 8{\bar{C}}(t-s){\mathbb {E}}\left[ \int _{s}^{t}(\lambda _1(v))^{2}\mathrm{d}v\right] +8(t-s)^{3}{\mathbb {E}}\left[ \int _{s}^{t}(\lambda _1(v))^{4}\mathrm{d}v\right] \nonumber \\&= 8{\bar{C}}(t-s)^{2}{\mathbb {E}}\left[ (\lambda _1(0))^{2}\right] +8(t-s)^{4}{\mathbb {E}}\left[ (\lambda _1(0))^{4}\right] . \end{aligned}$$

(A.16)

The first inequality in (A.16) uses the inequality $(\frac{x+y}{2})^{4}\le \frac{x^{4}+y^{4}}{2}$. The second inequality in (A.16) uses the martingality of $N_1 [s,t]-\int _{s}^{t}\lambda _1 (v) \mathrm{d}v$ with the predictable quadratic variation $\int _{s}^{t}\lambda _1 (v) \mathrm{d}v$ and the Burkholder–Davis–Gundy inequality, where ${\bar{C}}>0$ is a constant from the Burkholder–Davis–Gundy inequality.³ The third inequality in (A.16) uses Jensen’s inequality, so that

$$\begin{aligned} \left( \frac{1}{t-s}\int _{s}^{t}\lambda _1(v) \mathrm{d}v\right) ^{2}\le \frac{1}{t-s}\int _{s}^{t}(\lambda _1(v))^{2}\mathrm{d}v, \end{aligned}$$

and

$$\begin{aligned} \left( \frac{1}{t-s}\int _{s}^{t}\lambda _1(v) \mathrm{d}v\right) ^{4}\le \frac{1}{t-s}\int _{s}^{t}(\lambda _1(v))^{4}\mathrm{d}v. \end{aligned}$$

Finally, the last equality in (A.16) is due to stationarity of the intensity process; see Sect. 2.

Combining (A.15), (A.16) and (A.3), we deduce that (A.5) holds. The proof is thus complete. $\square$

1.2 Proof of Proposition 3

Proof of Proposition 3

We first prove Part (b), then prove Parts (a) and (c).

To prove Part (b), we recall that

$$\begin{aligned} \phi (t)=\frac{h(t)}{1-\Vert h\Vert _{L^{1}}} +\int _{0}^{\infty }h(t+v)\phi (v)\mathrm{d}v+\int _{0}^{t}h(t-v)\phi (v)\mathrm{d}v. \end{aligned}$$

(A.17)

Denote $H(t):=\int _t^{\infty }h(s)\mathrm{d}s$. By integrating Eq. (A.17) on both sides, we get

$$\begin{aligned} \Vert \phi \Vert _{L^{1}}&=\frac{\Vert h\Vert _{L^{1}}}{1-\Vert h\Vert _{L^{1}}} +\int _{0}^{\infty }H(v)\phi (v)\mathrm{d}v +\Vert h\Vert _{L^{1}}\Vert \phi \Vert _{L^{1}}\\&=\frac{\Vert h\Vert _{L^{1}}}{1-\Vert h\Vert _{L^{1}}} +\int _{0}^{M}H(v)\phi (v)\mathrm{d}v +\int _{M}^{\infty }H(v)\phi (v)\mathrm{d}v +\Vert h\Vert _{L^{1}}\Vert \phi \Vert _{L^{1}}. \end{aligned}$$

Since $\Vert h\Vert _{L^{1}}<1$, there exists some $\epsilon >0$ such that $\Vert h\Vert _{L^{1}}+\epsilon <1$. In addition, note that H(M) is decreasing in M to 0 as $M\rightarrow \infty$. Hence, for sufficiently large M we have $H(M)\le \epsilon$. This implies

$$\begin{aligned} \Vert \phi \Vert _{L^{1}}&\le \frac{\Vert h\Vert _{L^{1}}}{1-\Vert h\Vert _{L^{1}}} +H(0)\int _{0}^{M}\phi (v)\mathrm{d}v +\epsilon \Vert \phi \Vert _{L^{1}} +\Vert h\Vert _{L^{1}}\Vert \phi \Vert _{L^{1}}. \end{aligned}$$

Note that $\text {Var}(N^{1}(t))<\infty$ (see, for example, Lemma 2 in [58]) and hence

$$\begin{aligned} K(t)=\frac{t}{1-\Vert h\Vert _{L^{1}}}+2\int _{0}^{t}\int _{0}^{t_{2}}\phi (t_{2}-t_{1})\mathrm{d}t_{1}\mathrm{d}t_{2}<\infty . \end{aligned}$$

(A.18)

Moreover, $\int _{0}^{t_{2}}\phi (t_{2}-t_{1})\mathrm{d}t_{1}=\int _{0}^{t_{2}}\phi (t_{1})\mathrm{d}t_{1}$ is non-decreasing in $t_{2}$, and $\int _{0}^{t}\int _{0}^{t_{2}}\phi (t_{1})\mathrm{d}t_{1}\mathrm{d}t_{2}<\infty$ for every t. This implies that $\int _{0}^{t_{2}}\phi (t_{1})\mathrm{d}t_{1}<\infty$ for every $t_{2}$. Thus, $\int _{0}^{M}\phi (v)\mathrm{d}v<\infty$. Hence, it follows that

$$\begin{aligned} \Vert \phi \Vert _{L^{1}} \le \frac{\Vert h\Vert _{L^{1}}}{(1-\Vert h\Vert _{L^{1}})(1-\Vert h\Vert _{L^{1}}-\epsilon )} +\frac{H(0)\int _{0}^{M}\phi (v)\mathrm{d}v}{1-\Vert h\Vert _{L^{1}}-\epsilon }<\infty . \end{aligned}$$

(A.19)

To establish the second part of Part (b), we first note from the definition of K(t) in (3.1) that

$$\begin{aligned} K'(t)=\frac{1}{1-\Vert h\Vert _{L^{1}}}+2\int _{0}^{t}\phi (t-t_{1})\mathrm{d}t_{1}. \end{aligned}$$

The differentiability of $K(\cdot )$ is due to the integrability of $\phi$ as given in (A.19). Since $\phi$ is nonnegative, it follows that $K'(\cdot )$ is non-decreasing. Hence, $K(\cdot )$ is convex. In addition, we note that $\int _{0}^{t}\phi (t-t_{1})\mathrm{d}t_{1}=\int _{0}^{t}\phi (t_{1})\mathrm{d}t_{1} \le \Vert \phi \Vert _{L^{1}}< \infty$ for all t. Thus, $K(\cdot )$ is Lipschitz continuous.

We next provide a proof of Part (a), which will be useful in proving Part (c). Write $N^1$ for a stationary Hawkes process with baseline intensity $\mu =1$. From the Bartlett spectrum for the stationary Hawkes process (see [27] or [14]), we know that

$$\begin{aligned} \text{ Var }\left( \int _{{\mathbb {R}}}\psi (s)N^1(\mathrm{d}s)\right) =\int _{{\mathbb {R}}}|{\hat{\psi }}(\omega )|^{2}\frac{1}{2\pi (1-\Vert h\Vert _{L^{1}})}\frac{1}{|1-{\hat{h}}(\omega )|^{2}}\mathrm{d}\omega , \end{aligned}$$

where ${\hat{\psi }}(\omega ) = \int _{{\mathbb {R}}} e^{i\omega t} \psi (t) \mathrm{d}t$. Note that

$$\begin{aligned} K(t)=\text{ Var }(N^1(t)) =\text{ Var }\left( \int _{{\mathbb {R}}}1_{[0,t]}(s)N^1(\mathrm{d}s)\right) , \end{aligned}$$

and

$$\begin{aligned} \int _{{\mathbb {R}}}e^{i\omega s}1_{[0,t]}(s)\mathrm{d}s =\int _{0}^{t}e^{i\omega s}\mathrm{d}s =\frac{e^{i\omega t}-1}{i\omega }. \end{aligned}$$

We can compute that

$$\begin{aligned} \left| \frac{e^{i\omega t}-1}{i\omega }\right| ^{2} =\frac{(\cos (\omega t)-1)^{2}+\sin ^{2}(\omega t)}{\omega ^{2}} =\frac{2-2\cos (\omega t)}{\omega ^{2}} =\frac{\sin ^{2}(\frac{1}{2}\omega t)}{\left( \frac{1}{2}\omega \right) ^{2}}. \end{aligned}$$

Therefore,

$$\begin{aligned} K(t)&=\frac{1}{2\pi (1-\Vert h\Vert _{L^{1}})} \int _{{\mathbb {R}}}\frac{\sin ^{2}(\frac{1}{2}\omega t)}{\left( \frac{1}{2}\omega \right) ^{2}} \frac{1}{|1-{\hat{h}}(\omega )|^{2}}\mathrm{d}\omega \nonumber \\&=\frac{t}{2\pi (1-\Vert h\Vert _{L^{1}})} \int _{{\mathbb {R}}}\frac{\sin ^{2}\left( \frac{1}{2}\omega \right) }{\left( \frac{1}{2}\omega \right) ^{2}} \frac{1}{|1-{\hat{h}}(\frac{\omega }{t})|^{2}}\mathrm{d}\omega . \end{aligned}$$

(A.20)

Notice that, for any t,

$$\begin{aligned} \left| 1-{\hat{h}}\left( \frac{\omega }{t}\right) \right| \ge 1-\left| {\hat{h}}\left( \frac{\omega }{t}\right) \right| \ge 1-\Vert h\Vert _{L^{1}}>0, \end{aligned}$$

and $\int _{{\mathbb {R}}}\frac{\sin ^{2}\left( \frac{1}{2}\omega \right) }{\left( \frac{1}{2}\omega \right) ^{2}}\mathrm{d}\omega =2\pi$. Thus, $\frac{\sin ^{2}\left( \frac{1}{2}\omega \right) }{\left( \frac{1}{2}\omega \right) ^{2}} \frac{1}{|1-{\hat{h}}(\frac{\omega }{t})|^{2}} \le \frac{\sin ^{2}\left( \frac{1}{2}\omega \right) }{\left( \frac{1}{2}\omega \right) ^{2}} \frac{1}{(1-\Vert h\Vert _{L^{1}})^{2}}$, which is integrable. On the other hand, for every $\omega$, $\lim _{t\rightarrow \infty }{\hat{h}}(\frac{\omega }{t}) ={\hat{h}}(0)=\Vert h\Vert _{L^{1}}$. Therefore, by the dominated convergence theorem, we obtain

$$\begin{aligned} \lim _{t\rightarrow \infty }\frac{K(t)}{t}=\frac{1}{(1-\Vert h\Vert _{L^{1}})^{3}}. \end{aligned}$$

We now prove Part (c). This requires a more delicate analysis of the Bartlett spectrum. Write ${\overline{z}}$ for the complex conjugate of a complex number z. We can obtain from (A.20) and ${\hat{h}}(0)=\Vert h\Vert _{L^{1}}$ that

$$\begin{aligned}&K(t)-\frac{t}{(1-\Vert h\Vert _{L^{1}})^{3}} \nonumber \\&\quad =\frac{t}{2\pi (1-\Vert h\Vert _{L^{1}})}\int _{{\mathbb {R}}}\frac{\sin ^{2}\left( \frac{1}{2}\omega \right) }{\left( \frac{1}{2}\omega \right) ^{2}} \left[ \frac{1}{|1-{\hat{h}}(\frac{\omega }{t})|^{2}}-\frac{1}{|1-{\hat{h}}(0)|^{2}}\right] \mathrm{d}\omega \nonumber \\&\quad =\frac{t}{2\pi (1-\Vert h\Vert _{L^{1}})^{3}}\int _{{\mathbb {R}}}\frac{\sin ^{2}\left( \frac{1}{2}\omega \right) }{\left( \frac{1}{2}\omega \right) ^{2}} \cdot \frac{-2\Vert h\Vert _{L^{1}}+\Vert h\Vert _{L^{1}}^{2} +{\hat{h}}(\frac{\omega }{t})+\overline{{\hat{h}}(\frac{\omega }{t})} -|{\hat{h}}(\frac{\omega }{t})|^{2}}{|1-{\hat{h}}(\frac{\omega }{t})|^{2}}\mathrm{d}\omega \nonumber \\&\quad =\frac{1}{2\pi (1-\Vert h\Vert _{L^{1}})^{3}}\int _{{\mathbb {R}}}\frac{\sin ^{2}(\frac{1}{2}\omega t)}{\left( \frac{1}{2}\omega \right) ^{2}} \cdot \frac{-2\Vert h\Vert _{L^{1}}+\Vert h\Vert _{L^{1}}^{2} +{\hat{h}}(\omega )+\overline{{\hat{h}}(\omega )} -|{\hat{h}}(\omega )|^{2}}{|1-{\hat{h}}(\omega )|^{2}}\mathrm{d}\omega \nonumber \\&\quad =\frac{1}{2\pi (1-\Vert h\Vert _{L^{1}})^{3}}\int _{{\mathbb {R}}}\sin ^{2}\left( \frac{1}{2}\omega t\right) f(\omega )\mathrm{d}\omega \nonumber \\&\quad =\frac{1}{4\pi (1-\Vert h\Vert _{L^{1}})^{3}}\int _{{\mathbb {R}}}f(\omega )\mathrm{d}\omega -\frac{1}{4\pi (1-\Vert h\Vert _{L^{1}})^{3}}\int _{{\mathbb {R}}}\cos (\omega t)f(\omega )\mathrm{d}\omega , \end{aligned}$$

(A.21)

where

$$\begin{aligned} f(\omega )=\frac{1}{\left( \frac{1}{2}\omega \right) ^{2}} \cdot \frac{-2\Vert h\Vert _{L^{1}}+\Vert h\Vert _{L^{1}}^{2} +{\hat{h}}(\omega )+\overline{{\hat{h}}(\omega )} -|{\hat{h}}(\omega )|^{2}}{|1-{\hat{h}}(\omega )|^{2}}. \end{aligned}$$

We claim that $f(\omega )\in L^{1}({\mathbb {R}})$, i.e., f is integrable on the real line. To see this, notice first that

$$\begin{aligned} |f(\omega )|\le \frac{1}{\left( \frac{1}{2}\omega \right) ^{2}} \frac{4\Vert h\Vert _{L^{1}}+2\Vert h\Vert _{L^{1}}^{2}}{(1-\Vert h\Vert _{L^{1}})^{2}}, \end{aligned}$$

and thus $\int _{|\omega |\ge \epsilon }|f(\omega )|\mathrm{d}\omega <\infty$ for any $\epsilon >0$. Moreover, by L’Hôpital’s rule, we can check that

$$\begin{aligned} \lim _{\omega \rightarrow 0}f(\omega )&=\lim _{\omega \rightarrow 0}\frac{1}{\left( \frac{1}{2}\omega \right) ^{2}} \frac{-2\Vert h\Vert _{L^{1}}+\Vert h\Vert _{L^{1}}^{2} +{\hat{h}}(\omega )+\overline{{\hat{h}}(\omega )} -|{\hat{h}}(\omega )|^{2}}{|1-{\hat{h}}(\omega )|^{2}}\\&=\frac{4}{(1-\Vert h\Vert _{L^{1}})^{2}} \lim _{\omega \rightarrow 0}\frac{-2\Vert h\Vert _{L^{1}}+\Vert h\Vert _{L^{1}}^{2} +{\hat{h}}(\omega )+\overline{{\hat{h}}(\omega )} -|{\hat{h}}(\omega )|^{2}}{\omega ^{2}} \\&=\frac{4}{(1-\Vert h\Vert _{L^{1}})^{2}} \lim _{\omega \rightarrow 0}\frac{{\hat{h}}'(\omega )+\overline{{\hat{h}}'(\omega )} -{\hat{h}}'(\omega )\overline{{\hat{h}}(\omega )} -{\hat{h}}(\omega )\overline{{\hat{h}}'(\omega )}}{2\omega }. \end{aligned}$$

Since

$$\begin{aligned} {\hat{h}}'(0)+\overline{{\hat{h}}'(0)} -{\hat{h}}'(0)\overline{{\hat{h}}(0)} -{\hat{h}}(0)\overline{{\hat{h}}'(0)} =(1-\Vert h\Vert _{L^{1}})({\hat{h}}'(0)+\overline{{\hat{h}}'(0)}) =0, \end{aligned}$$

we apply the L’Hôpital’s rule again and get

$$\begin{aligned} \lim _{\omega \rightarrow 0}f(\omega )&=\frac{4}{(1-\Vert h\Vert _{L^{1}})^{2}} \lim _{\omega \rightarrow 0}\frac{{\hat{h}}'(\omega )+\overline{{\hat{h}}'(\omega )} -{\hat{h}}'(\omega )\overline{{\hat{h}}(\omega )} -{\hat{h}}(\omega )\overline{{\hat{h}}'(\omega )}}{2\omega }\\&=\frac{2}{(1-\Vert h\Vert _{L^{1}})^{2}} \left[ {\hat{h}}''(0)+\overline{{\hat{h}}''(0)} -{\hat{h}}''(0)\Vert h\Vert _{L^{1}}\right. \\&\qquad \left. -2{\hat{h}}'(0)\overline{{\hat{h}}'(0)} -\overline{{\hat{h}}''(0)}\Vert h\Vert _{L^{1}}\right] \\&=\frac{4}{(1-\Vert h\Vert _{L^{1}})^{2}} \left[ (1-\Vert h\Vert _{L^{1}})\int _{0}^{\infty }t^{2}h(t)\mathrm{d}t -\left( \int _{0}^{\infty }th(t)\mathrm{d}t\right) ^{2}\right] . \end{aligned}$$

Note that the first and second derivatives of ${\hat{h}}(\omega )$ and $\overline{{\hat{h}}(\omega )}$ are well-defined due to the assumption $\int _{0}^{\infty }t^{2}h(t)\mathrm{d}t<\infty$, and are given by

$$\begin{aligned}&{\hat{h}}(\omega )=\int _{0}^{\infty }(\cos \omega t+i\sin \omega t)h(t)\mathrm{d}t,\\&\overline{{\hat{h}}(\omega )}=\int _{0}^{\infty }(\cos \omega t-i\sin \omega t)h(t)\mathrm{d}t,\\&{\hat{h}}'(\omega )=\int _{0}^{\infty }(-\sin \omega t+i\cos \omega t)th(t)\mathrm{d}t,\\&\overline{{\hat{h}}'(\omega )}=\int _{0}^{\infty }(-\sin \omega t-i\cos \omega t)th(t)\mathrm{d}t,\\&{\hat{h}}''(\omega )=\int _{0}^{\infty }(-\cos \omega t-i\sin \omega t)t^{2}h(t)\mathrm{d}t,\\&\overline{{\hat{h}}''(\omega )}=\int _{0}^{\infty }(-\cos \omega t+i\sin \omega t)t^{2}h(t)\mathrm{d}t. \end{aligned}$$

We deduce from the above that $f(\omega )\in L^{1}({\mathbb {R}})$, and then the Riemann–Lebesgue theorem gives

$$\begin{aligned} \lim _{t\rightarrow \infty }\int _{{\mathbb {R}}}e^{i\omega t}f(\omega )\mathrm{d}\omega =0, \end{aligned}$$

which implies that the limit of the real part is also zero:

$$\begin{aligned} \lim _{t\rightarrow \infty }\int _{{\mathbb {R}}}\cos (\omega t)f(\omega )\mathrm{d}\omega =0. \end{aligned}$$

Hence, we conclude from (A.21) that

$$\begin{aligned}&\lim _{t\rightarrow \infty }\left[ K(t)-\frac{t}{(1-\Vert h\Vert _{L^{1}})^{3}}\right] \\&\quad =\frac{1}{\pi (1-\Vert h\Vert _{L^{1}})^{3}}\int _{{\mathbb {R}}}\frac{1}{\omega ^{2}} \frac{-2\Vert h\Vert _{L^{1}}+\Vert h\Vert _{L^{1}}^{2} +{\hat{h}}(\omega )+\overline{{\hat{h}}(\omega )} -|{\hat{h}}(\omega )|^{2}}{|1-{\hat{h}}(\omega )|^{2}}\mathrm{d}\omega \\&\quad =\frac{1}{\pi (1-\Vert h\Vert _{L^{1}})^{3}}\int _{{\mathbb {R}}}\frac{1}{\omega ^{2}} \frac{(1-\Vert h\Vert _{L^{1}})^{2}-|1-{\hat{h}}(\omega )|^{2}}{|1-{\hat{h}}(\omega )|^{2}}\mathrm{d}\omega <0, \end{aligned}$$

where the fact that this constant is negative follows from the observation that, for each $\omega ,$

$$\begin{aligned} |1-{\hat{h}}(\omega )|\ge 1-|{\hat{h}}(\omega )| \ge 1-\Vert h\Vert _{L^{1}}>0. \end{aligned}$$

Hence we complete the proof of Part (c). $\square$

1.3 Proof of Proposition 4

Proof of Proposition 4

We first prove that G has stationary increments. One can directly verify this fact by noting that, for $\tau >0, s \ge 0,$ $G(s+ \tau )- G(s)$ is a mean-zero Gaussian random variable, with variance given by

$$\begin{aligned} \text{ Var }(G(s +\tau )-G(s)) = K(s+\tau ) +K(s) - 2 \text{ Cov } (G(s+ \tau ), G(s)). \end{aligned}$$

Using (3.1) and (3.3), it is easily checked that $\text{ Var }(G(s +\tau )-G(s))= K(\tau )= \text{ Var }(G(\tau ))$, which is independent of s.

We next show that the Gaussian process G is not Markovian unless $h\equiv 0$. To see this, recall (see, for example, Revuz and Yor [51, p.86]) that a centered Gaussian process $\Upsilon$ with covariance function $\Gamma (s,t):={\mathbb {E}}[\Upsilon _s \Upsilon _t]$ is Markovian if and only if

$$\begin{aligned} \Gamma (s,u) \Gamma (t, t) = \Gamma (s,t) \Gamma (t,u), \end{aligned}$$

(A.22)

for every $0 \le s< t < u$. Given the covariance function of G in (3.3), one can directly check that (A.22) does not hold for any nonzero exciting function h.

Finally, we prove that the paths of G are Hölder continuous of order $\gamma$ for every $\gamma <\frac{1}{2}$. To see this, note that $G(s+ \tau )- G(s)$ is a mean-zero Gaussian random variable with variance $K(\tau )$, which implies that, for $p>0$,

$$\begin{aligned} {\mathbb {E}}[|G(s+\tau )-G(s)|^p] = C \cdot {K(\tau )}^{p/2}, \end{aligned}$$

where $C={\mathbb {E}}|Z|^p$ and Z follows a standard normal distribution. By the Lipschitz property of K in Proposition 3, we infer from the Kolmogorov–Chentsov theorem that the sample paths of G are Hölder continuous with order less than $\frac{1}{2}$. $\square$

Proof of results in Section 5

1.1 Proof of Theorem 12

Proof of Theorem 12

For notational simplicity, we write, for each t,

$$\begin{aligned} {\mathbb {N}} (t)=(N^{1}(t),\ldots ,N^{k}(t)), \end{aligned}$$

to stand for a k-dimensional stationary Hawkes processes where the intensity is given in (5.1) with $\mu =1$. It should be self-evident that here the notation $N^{i}$ stands for the i-th process (in Appendix A we used this notation to represent a univariate Hawkes process with baseline intensity i). Let us define

$$\begin{aligned} \tilde{{\mathbb {N}}}(t):={\mathbb {N}}(t)-at, \end{aligned}$$

(B.1)

and let $\tilde{{\mathbb {N}}}_{j}, j =1, 2, \ldots ,$ be independent copies of $\tilde{{\mathbb {N}}}$, where ${\mathbb {E}} [\tilde{{\mathbb {N}}}(t)] =0$ for each t. We obtain from the immigration-birth representation of multivariate Hawkes processes that

$$\begin{aligned} \mathbb {{{\hat{N}}}}^{(\mu )} (t) := \frac{{\mathbb {N}}^{(\mu )}(t)-{{\bar{\lambda }}} t}{\sqrt{\mu }} =\frac{1}{\sqrt{\mu }}\sum _{j=1}^{\mu }\tilde{{\mathbb {N}}}_{j}(t), \end{aligned}$$

(B.2)

where, as in the proof of Theorem 2, it suffices to establish the weak convergence of the sequence $\left( \mathbb {{{\hat{N}}}}^{(\mu )} \right)$ for positive integer-valued $\mu$.

We first establish the tightness of the sequence of processes $\left( \mathbb {{{\hat{N}}}}^{(\mu )} \right)$. We use the tightness criteria in [33, Chapter VI. Theorem 4.1] and verify the three conditions there. Condition (i) trivially holds. To verify Condition (ii) and (iii), it suffices to check the following two conditions: for every $0<T<\infty$, there exist some positive constants $C_{1},C_{2}$ such that

$$\begin{aligned} {\mathbb {E}}\left[ \Vert \tilde{{\mathbb {N}}}(u)-\tilde{{\mathbb {N}}}(s)\Vert ^{2}\right] \le C_{1} \cdot (u-s), \end{aligned}$$

(B.3)

and

$$\begin{aligned} {\mathbb {E}}\left[ \Vert \tilde{{\mathbb {N}}}(u)-\tilde{{\mathbb {N}}}(t)\Vert ^{2}\Vert \tilde{{\mathbb {N}}}(t)-\tilde{{\mathbb {N}}}(s)\Vert ^{2}\right] \le C_{2} \cdot (u-s)^{2}, \end{aligned}$$

(B.4)

for all $0\le s\le t\le u\le T$ with $u-s<1$, where the notation $|| \cdot ||$ stands for the usual Euclidean norm of a vector in ${\mathbb {R}}^k.$ To see this, first notice that using the Markov inequality, it is straightforward to verify that (B.3) implies Condition (ii) in [33, Chapter VI. Theorem 4.1]. In addition, following the proof of Theorem 2 in [26] (the processes considered there are real-valued, but the argument in that proof also works for ${\mathbb {R}}^k$-valued processes), one can immediately deduce that (B.4) implies Condition (iii) in [33, Chapter VI. Theorem 4.1].

We now prove (B.3) and (B.4). As the dimension k of the multivariate Hawkes process ${\mathbb {N}}$ is finite, in order to prove (B.3) and (B.4), it suffices to check that for every $0<T<\infty$, there exist some positive constants $C_{1},C_{2}$ such that, for all $0\le s\le t\le u\le T$ with $u-s<1$ and every $1\le i,j\le k$,

$$\begin{aligned} {\mathbb {E}}\left[ (N^{i}(s,u])^{2}\right] \le C_{1}\cdot (u-s), \end{aligned}$$

(B.5)

and

$$\begin{aligned} {\mathbb {E}}\left[ (N^{i}(t,u])^{2}(N^{j}(s,t])^{2}\right] \le C_{2} \cdot (u-s)^{2}. \end{aligned}$$

(B.6)

We next prove (B.5) and (B.6). Similarly to (A.12), we can compute that

$$\begin{aligned} {\mathbb {E}}\left[ (N^{i}(s,u])^{2}\right]&\le 2{\mathbb {E}}\left[ \left( N^{i}(s,u]-\int _{s}^{u}\lambda ^{i}(v)\mathrm{d}v\right) ^{2}\right] +2{\mathbb {E}}\left[ \left( \int _{s}^{u}\lambda ^{i}(v)\mathrm{d}v\right) ^{2}\right] \\&=2a_{i}(u-s) +2{\mathbb {E}}\left[ \left( \int _{s}^{u}\lambda ^{i}(v)\mathrm{d}v\right) ^{2}\right] \\&\le 2a_{i}(u-s) +2(u-s){\mathbb {E}}\left[ \left( \int _{s}^{u}(\lambda ^{i}(v))^{2}\mathrm{d}v\right) \right] \\&=2a_{i}(u-s) +2(u-s)^{2}{\mathbb {E}}[(\lambda ^{i}(0))^{2}]\le C_{1}(u-s), \end{aligned}$$

for some positive constant $C_{1}$, provided that ${\mathbb {E}}[(\lambda ^{i}(0))^{2}]<\infty$.

Moreover, similarly to the derivations in (A.13), (A.14) and (A.15), we have

$$\begin{aligned}&{\mathbb {E}}\left[ (N^{i}(t,u])^{2}(N^{j}(s,t])^{2}\right] \\&\quad \le 2{\mathbb {E}}\left[ \left( N^{i}(t,u]-\int _{t}^{u}\lambda ^{i}(v)\mathrm{d}v\right) ^{2}(N^{j}(s,t])^{2}\right] \\&\qquad +2{\mathbb {E}}\left[ \left( \int _{t}^{u}\lambda ^{i}(v)\mathrm{d}v\right) ^{2}(N^{j}(s,t])^{2}\right] \\&\quad = 2{\mathbb {E}}\left[ \int _{t}^{u}\lambda ^{i}(v)\mathrm{d}v\cdot (N^{j}(s,t])^{2}\right] +2{\mathbb {E}}\left[ \left( \int _{t}^{u}\lambda ^{i}(v)\mathrm{d}v\right) ^{2}(N^{j}(s,t])^{2}\right] \\&\quad \le 2\left( {\mathbb {E}}\left[ \left( \int _{t}^{u}\lambda ^{i}(v)\mathrm{d}v\right) ^{2}\right] \right) ^{1/2} \left( {\mathbb {E}}\left[ (N^{j}(s,t])^{4}\right] \right) ^{1/2}\\&\qquad +2(u-t){\mathbb {E}}\left[ \left( \int _{t}^{u}(\lambda ^{i}(v))^{2}\mathrm{d}v\right) (N^{j}(s,t])^{2}\right] \\&\quad \le 2(u-t)^{1/2}\left( \int _{t}^{u}{\mathbb {E}}(\lambda ^{i}(v))^{2}\mathrm{d}v\right) ^{1/2} {\mathbb {E}}\left[ (N^{j}(s,t])^{4}\right] ^{1/2}\\&\qquad +(u-t){\mathbb {E}}\left[ \int _{t}^{u}\left( (\lambda ^{i}(v))^{4}+(N^{j}(s,t])^{2}\right) \mathrm{d}v\right] \\&\quad = 2(u-t)\left( {\mathbb {E}}\left[ (\lambda ^{i}(0))^{2}\right] \right) ^{1/2} \left( {\mathbb {E}}\left[ (N^{j}(0,t-s])^{4}\right] \right) ^{1/2}\\&\qquad +(u-t)^{2}\left( {\mathbb {E}}\left[ (\lambda ^{i}(0))^{4}\right] +{\mathbb {E}}\left[ (N^{j}(0,t-s])^{2}\right] \right) \\&\quad \le 2(u-t)\left( {\mathbb {E}}\left[ (\lambda ^{i}(0))^{2}\right] \right) ^{1/2} \left( {\mathbb {E}}\left[ (N^{j}(0,t-s])^{4}\right] \right) ^{1/2}\\&\qquad +(u-t)^{2}\left( {\mathbb {E}}\left[ (\lambda ^{i}(0))^{4}\right] +C_{1}(t-s)\right) . \end{aligned}$$

Similarly to (A.16) in the proof of Theorem 2, we have

$$\begin{aligned} {\mathbb {E}}\left[ (N^{j}(0,t-s])^{4}\right] \le 8{\bar{C}}(t-s)^{2}{\mathbb {E}}[(\lambda ^{j}(0))^{2}] +8(t-s)^{4}{\mathbb {E}}[(\lambda ^{j}(0))^{4}]. \end{aligned}$$

Hence, we obtain

$$\begin{aligned} {\mathbb {E}}\left[ (N^{i}(t,u])^{2}(N^{j}(s,t])^{2}\right] \le C_{2}(u-s)^{2}, \end{aligned}$$

for some positive constant $C_{2}$, provided that ${\mathbb {E}}[(\lambda ^{j}(0))^{4}]<\infty$ for $1 \le j \le k$.

It remains to prove that ${\mathbb {E}}[(\lambda ^{j}(0))^{4}]<\infty$ for $1 \le j \le k$, as this implies ${\mathbb {E}}[(\lambda ^{j}(0))^{2}]<\infty$. Similarly to the derivations in (A.9)–(A.11) in the proof of Theorem 2, since $h_{ij}$ are locally bounded and Riemann integrable by Assumption 1, for sufficiently small $\delta >0$ we obtain

$$\begin{aligned} {\mathbb {E}}[(\lambda ^{j}(0))^{4}]&={\mathbb {E}}\left[ \left( p_{j}+\sum _{\ell =1}^{k}\int _{-\infty }^{0-}h_{j\ell }(-s)N^{\ell }(\mathrm{d}s)\right) ^{4}\right] \\&\le C \left[ (p_{j})^{4} +\sum _{\ell =1}^{k}{\mathbb {E}}\left[ \left( \int _{-\infty }^{0-}h_{j\ell }(-s)N^{\ell }(\mathrm{d}s)\right) ^{4}\right] \right] \\&\le C \left[ (p_{j})^{4} +\sum _{\ell =1}^{k}\left( \sum _{i=0}^{\infty }\max _{t\in [-(i+1)\delta ,-i\delta ]}h_{j\ell }(t)\right) ^{4} {\mathbb {E}}\left[ \left( N^{\ell }(0,\delta ]\right) ^{4}\right] \right] , \end{aligned}$$

for some positive constant C. Thus, it remains to show that for every $1\le \ell \le k$, ${\mathbb {E}}\left[ \left( N^{\ell }(0,1]\right) ^{4}\right] <\infty$. It suffices to show that there exists some constant $c_{\ell }>0$ so that ${\mathbb {E}}[e^{c_{\ell }N^{\ell }(0,1]}]<\infty$. Let us define ${\mathbb {E}}^{\emptyset }$ as the expectation under which the process ${\mathbb {N}}=(N^{1},\ldots ,N^{k})$ (with slight abuse of notation) is a multivariate Hawkes process starting from an empty history, that is, $\lambda ^{i}(t)=p_{i}+\sum _{j=1}^{k}\int _{0}^{t-}h_{ij}(t-s)N^{j}(\mathrm{d}s)$. For any $\psi _{i}>0$, $e^{\sum _{i=1}^{k}\psi _{i}N^{i}(t)-\sum _{i=1}^{k}(e^{\psi _{i}}-1)\int _{0}^{t}\lambda ^{i}(s)\mathrm{d}s}$ is a martingale (see, for example, [52]), and thus

$$\begin{aligned} 1&={\mathbb {E}}^{\emptyset }\left[ e^{\sum _{i=1}^{k}\psi _{i}N^{i}(t)-\sum _{i=1}^{k}(e^{\psi _{i}}-1)\int _{0}^{t}\lambda ^{i}(s)\mathrm{d}s}\right] \\&={\mathbb {E}}^{\emptyset }\left[ e^{\sum _{i=1}^{k}\psi _{i}N^{i}(t) -\sum _{i=1}^{k} (e^{\psi _{i}}-1) \int _{0}^{t}(p_{i}+\int _{0}^{s}\sum _{j=1}^{k}h_{ij}(s-u)N^{j}(\mathrm{d}u))\mathrm{d}s} \right] \\&\ge {\mathbb {E}}^{\emptyset }\left[ e^{\sum _{i=1}^{k}\psi _{i}N^{i}(t)-\sum _{i=1}^{k}(e^{\psi _{i}}-1)(p_{i}t+\sum _{j=1}^{k}\Vert h_{ij}\Vert _{L^{1}}N^{j}(t))}\right] , \end{aligned}$$

which implies that

$$\begin{aligned} {\mathbb {E}}^{\emptyset }\left[ e^{\sum _{i=1}^{k}(\psi _{i}-\sum _{j=1}^{k}(e^{\psi _{j}}-1)\Vert h_{ji}\Vert _{L^{1}})N^{i}(t)}\right] \le e^{\sum _{i=1}^{k}(e^{\psi _{i}}-1)p_{i}t}. \end{aligned}$$

(B.7)

Since the spectral radius of the matrix ${\mathbb {H}}= (\Vert h_{ij}\Vert _{L^{1}})_{1 \le i,j \le k}$ is strictly less than 1, we know that $({\mathbb {I}}-{\mathbb {H}})^{-1}$ exists and $({\mathbb {I}}-{\mathbb {H}})^{-1}=\sum _{n=0}^{\infty }{\mathbb {H}}^{n}$. Thus, for any fixed positive column vector $m=(m_{1},\ldots ,m_{k})^{t}\in {\mathbb {R}}_{>0}^{k}$, we have $(({\mathbb {I}}-{\mathbb {H}})^{-1}m)_{i}>0$ for every i, where $(({\mathbb {I}}-{\mathbb {H}})^{-1}m)_{i}$ is the i-th component of the vector $({\mathbb {I}}-{\mathbb {H}})^{-1}m$. Let $\psi _{i}=\epsilon (({\mathbb {I}}-{\mathbb {H}})^{-1}m)_{i}$, where $\epsilon >0$ is sufficiently small so that we can find some constant C(m) that depends on m such that

$$\begin{aligned} \psi _{i}-\sum _{j=1}^{k}(e^{\psi _{j}}-1){\mathbb {H}}_{ji} \ge \psi _{i}-\sum _{j=1}^{k}\psi _{j}{\mathbb {H}}_{ji}-C(m)\epsilon ^{2} =\epsilon m_{i}-C(m)\epsilon ^{2}>0. \end{aligned}$$

Hence, we deduce from (B.7) that there exists $c_{i}=\epsilon m_{i}-C(m)\epsilon ^{2}>0$ such that

$$\begin{aligned} {\mathbb {E}}^{\emptyset }\left[ e^{\sum _{i=1}^{k}c_{i}N^{i}(t)}\right] \le e^{t\sum _{i=1}^{k}p_{i}(e^{\epsilon (({\mathbb {I}}-{\mathbb {H}})^{-1}m)_{i}}-1)}. \end{aligned}$$

In particular, for every $1\le \ell \le k$, we have

$$\begin{aligned} {\mathbb {E}}^{\emptyset }\left[ e^{c_{\ell }N^{\ell }(t)}\right] \le e^{t\sum _{i=1}^{k}p_{i}(e^{\epsilon (({\mathbb {I}}-{\mathbb {H}})^{-1}m)_{i}}-1)}. \end{aligned}$$

Since the linear Hawkes process (either with an empty history or the stationary version) is associated, we then deduce that, for positive integer t,

$$\begin{aligned} \prod _{n=1}^{t}{\mathbb {E}}^{\emptyset }\left[ e^{c_{\ell }N^{\ell }(n-1,n]}\right] \le {\mathbb {E}}^{\emptyset }\left[ e^{c_{\ell }N^{\ell }(t)}\right] \le e^{t\sum _{i=1}^{k}p_{i}(e^{\epsilon (({\mathbb {I}}-{\mathbb {H}})^{-1}m)_{i}}-1)}. \end{aligned}$$

Hence, by the ergodicity of the Hawkes processes with an empty history where the exciting function satisfies Assumption 1, see, for example, [8], we obtain

$$\begin{aligned} \log {\mathbb {E}}\left[ e^{c_{\ell }N^{\ell }(0,1)}\right]= & {} \lim _{t\rightarrow \infty }\frac{1}{t}\sum _{n=1}^{t}\log {\mathbb {E}}^{\emptyset }\left[ e^{c_{\ell }N^{\ell }(n-1,n]}\right] \\\le & {} \sum _{i=1}^{k}p_{i}(e^{\epsilon ((I-{\mathbb {H}}^{t})^{-1}m)_{i}}-1)<\infty , \end{aligned}$$

where we recall that ${\mathbb {E}}$ is the expectation under which the Hawkes process is stationary. Hence, we have proved that there exists some constant $c_{\ell }>0$ such that ${\mathbb {E}}[e^{c_{\ell }N^{\ell }(0,1]}]<\infty$ for each $1 \le l \le k$.

Now we have proved the tightness of the sequence $\left( \mathbb {{{\hat{N}}}}^{(\mu )} \right)$, we next show that the finite-dimensional distributions of the sequence of processes $\left( \mathbb {{{\hat{N}}}}^{(\mu )} \right)$ converge in distribution to those of the limiting process ${\mathbb {G}}$ as $\mu \rightarrow \infty$. To this end, we note that one can readily compute the covariance function of $\tilde{{\mathbb {N}}}$ as in the univariate case [see Eqs. (A.7)–(A.8)], and find that, for $t\ge s$,

$$\begin{aligned} \text {Cov}(\tilde{{\mathbb {N}}} (t),\tilde{{\mathbb {N}}}(s)) =\int _{s}^{t}\int _{0}^{s}\Phi (u-v)\mathrm{d}u\mathrm{d}v+{\mathbb {K}}(s) = \text {Cov}({\mathbb {G}}(t),{\mathbb {G}}(s)). \end{aligned}$$

Hence, in view of (B.2) and (B.3), the weak convergence of the finite-dimensional distributions of this sequence $\left( \mathbb {{{\hat{N}}}}^{(\mu )} \right)$ immediately follows from the central limit theorem for sums of i.i.d. random vectors and the Cramér–Wold device (see, for example, Section 4.3.2 in [53]).

Finally, as the covariance of ${\mathbb {G}}$ is the same as that of $\tilde{{\mathbb {N}}}$, we immediately infer from (B.3) that, for $s<u,$

$$\begin{aligned} {\mathbb {E}}\left[ \Vert {\mathbb {G}}(u)-{\mathbb {G}}(s)\Vert ^{2}\right] \le C_{1} \cdot (u-s), \end{aligned}$$

which implies that the limiting Gaussian process ${\mathbb {G}}$ has continuous sample paths [51, p.37]. The proof is therefore complete. $\square$

1.2 Proof of Proposition 13

Proof of Proposition 13

To prove Proposition 13, we rely on [40]. For notational simplicity, we prove the joint weak convergence of $({\mathbb {X}}_1^{\mu }, \ldots , {\mathbb {X}}_k^{\mu })$ as $\mu \rightarrow \infty$ for the case $k=2.$ The general case $k \ge 2$ follows similarly.

First, as in the proof of Theorem 3 in [40], we obtain that, for $i=1, 2$,

$$\begin{aligned} {\mathbb {X}}_i^{\mu }(t)= & {} {\sqrt{\mu }} \left( \frac{{\mathbb {Q}}_i^{\mu }(t)}{\mu } - q_{i0} (1-F_{i0}(t)) - a_i \cdot \int _{0}^t (1-F_i(t-u)) \mathrm{d}u \right) , \\= & {} \mu ^{-1/2} \sum _{j=1}^{{\mathbb {Q}}^{\mu }_i(0)} (1_{{{\bar{\eta }}}_{i,j} >t} - (1- F_{i0}(t))) + (1 - F_{i0}(t)) \mu ^{1/2}(\mu ^{-1}{\mathbb {Q}}^{\mu }_i(0) - q_{i0}) \\&+[M^{\mu }_{i1}(t) -M^{\mu }_{i2}(t)], \end{aligned}$$

where, for $t \ge 0,$

$$\begin{aligned} M^{\mu }_{i1}(t)&= \int _{0}^t (1 - F_i(t-s)) \mathrm{d} \left[ \frac{{\mathbb {N}}_i^{(\mu )}(s)-{{\bar{\lambda }}}_i s}{\sqrt{\mu }} \right] ,\\ M^{\mu }_{i2}(t)&= \int _{0}^t \int _{0}^t 1_{s+x \le t} \mathrm{d} U_i^\mu \left( \frac{{\mathbb {N}}_i^{(\mu )}(s)}{{\mu }} , F_i(x)\right) ,\\ U_i^\mu (t, x)&=\mu ^{-1/2} \sum _{j=1}^{\lfloor \mu t \rfloor } (1_{\zeta _{ij} \le x} -x). \end{aligned}$$

Here $\zeta _{ij}$ are all independent and uniformly distributed random variables on [0, 1] and the service times $\eta _{ij}=F_i^{-1}(\zeta _{ij})$, where $F_i^{-1}(x):= \inf \{y: F_i(y) \ge x \}$. For each fixed i, it was proved in [40] (see (6.1)–(6.3) there) that the following weak convergence of processes hold:

$$\begin{aligned}&\mu ^{-1/2} \sum _{j=1}^{{\mathbb {Q}}^{\mu }_i(0)} (1_{{{\bar{\eta }}}_{i,j} >t} - (1- F_{i0}(t))) \Rightarrow \sqrt{q_{i0}} W^{i0}(F_{i0}(t)), \end{aligned}$$

(B.8)

$$\begin{aligned}&(1 - F_{i0}(t)) \mu ^{1/2}(\mu ^{-1}{\mathbb {Q}}^{\mu }_i(0) - q_{i0}) \Rightarrow (1 - F_{i0}(t)) \xi _i, \end{aligned}$$

(B.9)

$$\begin{aligned}&\left( M^{\mu }_{i1}(t), M^{\mu }_{i2}(t) \right) \Rightarrow \left( \int _{0}^t (1- F_i(t-u)) \mathrm{d}{\mathbb {G}}_i(u), \int _0^t \int _0^t 1_{s +x \le t} \mathrm{d}U_i \left( a_i s, F_i(x) \right) \right) . \nonumber \\ \end{aligned}$$

(B.10)

In addition, there is clearly a joint weak convergence of the left-hand sides of (B.8)–(B.10) to the right-hand sides [40]. Now, by the hypothesis, the number of customers in the system at time zero, ${\mathbb {Q}}^{\mu }_i(0)$ for $i=1, 2$, as well as their respective service requirements ${{\bar{\eta }}}_{ij}$ for $i=1, 2,$ are mutually independent. Moreover, the arrivals of new customers and the service requirements of those new customers are independent of the initial number of customers ${\mathbb {Q}}^{\mu }(0)$ and their service times. Hence, in order to prove the joint weak convergence of $({\mathbb {X}}_1^{\mu }, {\mathbb {X}}_2^{\mu })$, it suffices to prove the weak convergence of $(M^{\mu }_{11}, M^{\mu }_{12}, M^{\mu }_{21}, M^{\mu }_{22})$ as $\mu \rightarrow \infty .$

To this end, let us define

$$\begin{aligned} {\tilde{M}}^{\mu }_{i2} = \int _{0}^t \int _{0}^t 1_{s+x \le t} \mathrm{d} U_i^\mu \left( a_i s, F_i(x)\right) , \quad \text {for }i=1, 2. \end{aligned}$$

Note that by Theorem 12, we have that the sequence of processes $\left( \frac{{\mathbb {N}}^{(\mu )}}{\mu }\right)$ converges in distribution to a deterministic limit process $\omega$, where $\omega (t): =at$ for each $t \ge 0.$ As $\omega$ has continuous paths and the Skorohod $J_1$ topology relativized to the space of continuous functions coincides with the uniform topology there [5, p.124], we obtain that, for each $T>0$, as $\mu \rightarrow \infty$,

$$\begin{aligned} \sup _{0 \le t \le T} ||{{\mathbb {N}}^{(\mu )}(t)}/{\mu } - at || \rightarrow 0 \quad \text {in probability}. \end{aligned}$$

Then, using a similar argument as in the proof of Lemma 5.3 in [40], we can establish that, for each $T>0$ and $\epsilon >0$,

$$\begin{aligned} \lim _{\mu \rightarrow \infty } P\left( \sup _{t \le T} |{\tilde{M}}^{\mu }_{i2}(t) - {M}^{\mu }_{i2} (t)|> \epsilon \right) =0, \quad \text {for }i=1, 2. \end{aligned}$$

(B.11)

In addition, using integration by parts, we can write $(M^{\mu }_{11}, M^{\mu }_{21}) = \left( g_1( \mathbb {{{\hat{N}}}}_1^{(\mu )} ), g_2 (\mathbb {{{\hat{N}}}}_2^{(\mu )} ) \right)$, where $\mathbb {{{\hat{N}}}}_i^{(\mu )} (s) := \frac{{\mathbb {N}}_i^{(\mu )}(s)-{{\bar{\lambda }}}_i s}{\sqrt{\mu }}$ for each $s \ge 0$, $g_i: D([0,\infty ),{\mathbb {R}}) \rightarrow D([0,\infty ),{\mathbb {R}})$ is defined by

$$\begin{aligned} g_i(x(\cdot )) (t) = x(t) - \int _{0}^t x(t-s) \mathrm{d}F_i(s), \quad \hbox { for}\ i=1, 2, \end{aligned}$$

and $g_i$ is continuous at points $x(\cdot ) \in C([0,\infty ),{\mathbb {R}})$. See the proof of Lemma 3.3 in [40]. Now in Theorem 12 we have established that $\left( \mathbb {{{\hat{N}}}}_1^{(\mu )}, \mathbb {{{\hat{N}}}}_2^{(\mu )} \right)$ converges in distribution to the Gaussian process $({\mathbb {G}}_1, {\mathbb {G}}_2)$ under the Skorohod $J_1$ topology, where the limiting Gaussian process has continuous paths; it then immediately follows that

$$\begin{aligned} (M^{\mu }_{11}, M^{\mu }_{21}) \Rightarrow (M_{11}, M_{21}), \quad \text {as }\mu \rightarrow \infty , \end{aligned}$$

(B.12)

where

$$\begin{aligned} M_{i1} (t) = \int _{0}^t (1 - F_i(t-s)) \mathrm{d} {\mathbb {G}}_i(s), \quad \text {for }i=1, 2. \end{aligned}$$

Furthermore, as the service processes of each class i customer are independent, we deduce that the two processes $U_1^\mu$ and $U_2^\mu$ are independent for each $\mu$, which further implies that ${{\tilde{M}}}^{\mu }_{12}$ and ${{\tilde{M}}}^{\mu }_{22}$ are two independent processes. By Lemma 3.1 of [40], we have $U_i^\mu \Rightarrow U_i$ in $D([0,\infty ),D[0,1])$ for each i as $\mu \rightarrow \infty$. Hence, we deduce from Lemma 5.3 of [40] that

$$\begin{aligned} ({{\tilde{M}}}^{\mu }_{12}, {{\tilde{M}}}^{\mu }_{22}) \Rightarrow (M_{12}, M_{22}), \quad \text {as }\mu \rightarrow \infty , \end{aligned}$$

(B.13)

where

$$\begin{aligned} M_{i2}(t) = \int _{0}^t \int _{0}^t 1_{s+x \le t} \mathrm{d} U_i \left( a_i t, F_i(x)\right) , \quad \text {for } i=1, 2. \end{aligned}$$

Then, we can obtain from (B.12), (B.13) and the independence of the service processes and the arrival processes of each class of customers that

$$\begin{aligned} (M^{\mu }_{11}, {{\tilde{M}}}^{\mu }_{12}, M^{\mu }_{21}, {{\tilde{M}}}^{\mu }_{22}) \Rightarrow (M_{11}, M_{12}, M_{21}, M_{22}), \quad \text {as } \mu \rightarrow \infty . \end{aligned}$$

Together with (B.11), which implies that $\left( {M}^{\mu }_{12} - {\tilde{M}}^{\mu }_{12}, {M}^{\mu }_{22} - {\tilde{M}}^{\mu }_{22} \right) \Rightarrow (0, 0)$, we infer that the process $(M^{\mu }_{11}, M^{\mu }_{12}, M^{\mu }_{21}, M^{\mu }_{22})$ converges in distribution to $(M_{11}, M_{12}, M_{21}, M_{22})$ as $\mu \rightarrow \infty$. Therefore, we obtain the weak convergence of $({\mathbb {X}}_1^{\mu }, {\mathbb {X}}_2^{\mu })$ to the desired limit process $({\mathbb {X}}_1, {\mathbb {X}}_2)$. The proof is completed. $\square$