<figure><img src="image_1.png" /></figure>

HAL
open science

# Dynamic estimation of stochastic descriptors in the SPIR model using particle filtering

Vasileios E Papageorgiou, Irene Votsi, George Tsaklidis

<figure><img src="image_2.png" /></figure>

## To cite this version:

Vasileios E Papageorgiou, Irene Votsi, George Tsaklidis. Dynamic estimation of stochastic descriptors in the SPIR model using particle filtering. 2025. hal-05104613

HAL Id: hal-05104613

https://hal.science/hal-05104613v1

Preprint submitted on 10 Jun 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.
Dynamic estimation of stochastic descriptors in the SPIR model
using particle filtering

Vasileios E. Papageorgiou¹, Irene Votsi², George Tsaklis¹

¹Department of Mathematics, Aristotle University of Thessaloniki, Thessaloniki, Greece

²Laboratoire LIEC, CNRS, Université de Lorraine, Metz, F-57000, France

Correspondence: vpapageor@math.auth.gr

**Abstract**

In epidemiology, accurately assessing epidemic events deserves particular attention and becomes more and more urgent. Numerous studies have introduced different stochastic models and associated epidemiological descriptors, providing explicit formulas and algorithms that enable their calculation. Most of the studies consider constant epidemic dynamics, which limit the accuracy and the possibility of providing online estimates. In this paper, we introduce the SPIR (Susceptible, Presymptomatic, Infectious, Removed) model based on a 3-dimensional Markov chain. In the framework of this model, stochastic descriptors are introduced including the total number of infections, the occurrence time of a specific number of deaths and the number of infections generated from an index presymptomatic or an infectious case. Analytical expressions are obtained along with corresponding algorithms based on recursive algorithmic schemes. A new methodology is described to provide more precise estimates, which is based on an augmented state space model with dynamically updated parameters. The effectiveness of the proposed methodology is compared to the standard approach (with fixed parameters) by evaluating two descriptors: the total number of infectious cases and the time until the first death, using real monkeypox (mpox) data from Ghana in 2022. Additionally, the impact of infections due to contact with presymptomatic or an infectious case is investigated. In summary, incorporating time-varying parameter estimates through the filtering process enhances the precision of descriptors, allowing for more realistic epidemic scenarios. The ability to provide data-driven, online predictions while managing the inherent uncertainty of physical phenomena is a crucial tool in addressing emerging public health challenges.

**Keywords:** Particle Filtering, Markov Processes, State Space Models, Probability Distributions, Hitting Times, Epidemics

**MSC:** 62M20, 60J22, 65C40, 62G30, 62P10

**1. Introduction**

Deterministic models are widely used in epidemiological studies because of their clear interpretation and simple structure [1]. They are described by systems of differential equations corresponding to compartmental models, with the SIS and SIR models being fundamental in the field [2–4]. Other widely used models include the SEIR which incorporates the number of exposed individuals, the SEIRS which allows for multiple infections [5], and the SEIRD where infected cases may either die or recover from a disease [6]. Apparently, the aforementioned models can be extended in many ways including control descriptors, demographics or externals factors [7–11]. Despite their simplicity, deterministic approaches often fail to capture the dynamics of phenomena characterized by increased uncertainty such as epidemics [12,13]. On
the other hand, stochastic modelling proves valuable in describing epidemic events that exhibit
randomness in demography, disease transmission and environment rendering their utilization
necessary [14,15].

Recent advancements in stochastic modeling have introduced more sophisticated
frameworks to account for various complexities in epidemic scenarios. For example, the SEIQS
model with an external source of infection, explored in [16], integrates external infection
sources and quarantine states to model real-world complexities such as cross-border disease
transmission and containment descriptors. In parallel, the authors in [17] examine the influence
of random environmental factors on the SIS model, revealing the substantial effects of such
variability on the epidemic trajectory. These studies extend traditional models by incorporating
more realistic assumptions about disease transmission, providing a better understanding of
epidemic behavior in dynamically changing conditions.

In this paper we introduce the SPIR compartmental model, which considers the existence
of susceptible, presymptomatic, symptomatic, asymptomatic, recovered and deceased cases.
The total of symptomatic and asymptomatic cases constitutes the class of infectious individuals,
while an interesting attribute is that a new case can be generated either from a presymptomatic
or an infectious case. This model is well-suited for describing most epidemic outbreaks, as it
incorporates the main elements typically observed in disease transmission. Moreover, it avoids
overcomplicating the structure by limiting itself to the essential states, helping to prevent
overfitting. Furthermore, it can be easily adapted to fit the specific characteristics of the
epidemic being studied.

Further, research has focused on the implications of specific epidemiological features and
intervention strategies. Studies like [18–21] highlight the challenges of imperfect vaccination
and the reintroduction of infections, while [22] addresses the critical role of asymptomatic
carriers in disease propagation. The complexity of these factors is further expanded in works
such as [9] which considers the constraints on quarantine resources and their impact on disease
control. Additionally, analyses like [23] provide a deeper understanding of the implications of
discrete-time approximations for continuous-time models, which is crucial for practical
applications of these models.

Another field which gathered attention, especially during the last decade, is the modeling
of infectious diseases through Markov processes and the introduction of novel stochastic
descriptors. These methods adapt the Markovian perspective to epidemiological schemes suited
for analyzing diseases affecting smaller population sizes, where the consideration of constant
disease dynamics does not significantly diverge from reality. The foundational work by Artalejo
et al. in [24–27] provides an analytical framework for assessing the long-term behavior of
epidemics using stochastic performance descriptors, highlighting the significance of
uncertainty in model predictions. This study is complemented by research on recovery
dynamics, such as [27], which investigates the probabilistic properties of recovery processes in
fundamental epidemic models. These works underscore the necessity of incorporating
stochastic elements to better capture the real-world uncertainties in disease spread and recovery,
which deterministic models may overlook.

However, only a few of them incorporate real data to validate their modeling effectiveness.
The parameters in most of these models are assumed to be constant throughout the outbreak.
While this approach can adequately describe scenarios where disease transmission shows little
to no variation, Papageorgiou & Tsaklidis [12,28] suggest that a more comprehensive approach
would be to combine innovative stochastic descriptors with a dynamic parameter updating
mechanism.

For our purposes, in the context of SPIR, computational schemes and algorithms are
provided for a total of four stochastic descriptors, focusing on size-, time-, and individual-level
features. Specifically, emphasis is placed on calculating the total number of infectious cases,
the time until the occurrence of the first death, and the number of infections generated by a
presymptomatic or an infectious case. The selection of the first two descriptors facilitates the comparison between the introduced hybrid methodology with the standard Markovian approach in which parameters are considered to be constant. The other two descriptors are novel stochastic measures introduced in the literature, providing estimates about the future trajectory of the epidemic, whether it will expand or die out. Moreover, extensive comparisons between these two descriptors demonstrate their impact on the spread of the disease, enforcing public authorities with valuable information on where to focus their efforts.

To estimate the parameters of the stochastic models, Bayesian filtering approaches have been widely employed, especially for estimating the progression of the COVID-19 outbreak, with most articles focusing on the use of Kalman filtering and its non-linear extensions [29–40]. Fewer, but still significant, are the approaches that employ sampling techniques adapted to extensions of the SIRS model [41–44]. These data-driven methods account for the uncertainty inherent in the spread of epidemic phenomena by using prior distributions selected according to the characteristics of each outbreak. Moreover, this allows for the generation of predictions accompanied by confidence intervals at any specified confidence level.

In this article, efficient algorithmic schemes for the estimation of stochastic descriptors are enhanced by the integration of time-varying parameters deriving from particle filter. The incorporation of the particle filter has a dual positive impact, as it provides accurate estimates for the model's states and parameters. Moreover, through the sampling procedure we can capture more intricate data patterns, enabling the representation of density distributions of any form. Also, through the efficient selection of transition and measurement density functions, properties such as the non-negativity of states and parameters can be accounted for.

The present analysis seeks to address the existing gap in enhancing the precision of stochastic descriptor estimates by incorporating the fluctuating dynamics of epidemics into the estimation process, while explicitly accounting for the inherent randomness in the progression of infectious diseases. In contrast to previous studies that have primarily focused on predicting disease progression, our objective is to complement these approaches by emphasizing additional influential factors, which arise from the stochastic nature of the underlying Markovian framework and shape the future trajectory of an epidemic. Moreover, the proposed methodology is well-suited for generating online predictions, offering a distinct advantage over standard approaches that rely on fixed model dynamics and are better suited for offline estimation.

## 2. Methodology

This section is dedicated to the methodological framework which consists of three main components: (a) the introduction of the SPIR epidemiological model, (b) the presentation of robust, low-complexity algorithms for estimating four stochastic descriptors, and (c) the estimation process of epidemiological parameters using the particle filtering approach. The first part outlines a compartmental model that is particularly suitable for real epidemic scenarios, such as the mpox outbreak. Additionally, the introduction of a new Markovian model is detailed. Section 2.2 provides formulas and algorithms for estimating the total number of infectious cases, the time until the onset of a specific number of deaths and the number of infections generated by a presymptomatic or an infectious individual. Finally, Section 2.3 describes the integration of the particle filter for the online estimation of epidemiological parameters.

### 2.1. The SPIR epidemiological model

For the purposes of the present analysis, we introduce a fundamental epidemiological model which consists of four compartments (states) that describe the evolution of the number of
susceptible (S), presymptomatic (P), infectious (I) and removed (R) individuals with respect to time. It should be noted that the introduced model captures the primary characteristics of most epidemic outbreaks without imposing significant computational burden. More importantly, complex structures may lead to misleading conclusions because the parameter estimation procedure is typically based on datasets limited to records of cases and deaths. Consequently, as the number of estimated parameters increases, the model becomes more prone to overfitting.

Presymptomatic cases are individuals who have been infected with the virus and are capable of transmitting it, often at a lower rate, as they are still within the incubation period. Once the presymptomatic phase ends, these individuals transition to the infectious phase. Specifically, a susceptible individual moves to the compartment of presymptomatic cases (presymptomatic phase) after contact with an individual of this compartment (with rate $β_1$) or after contact with individuals of the infectious compartment (with rate $β_2$).

Parameters $β_1$ and $β_2$ represent the infection rates of the system. Based on the above, we have two sources of infection, as there are several epidemics where transmission involves not only infectious but also presymptomatic individuals of the system [45,46]. We use different transmission rates, as the values of these parameters are not necessarily identical. Then, the presymptomatic individuals, transition (on average) after a time interval of $1/a$ to the infectious compartment. The rate of transition from the presymptomatic to the infectious phase ($α$), is often referred to as the incubation rate of the disease [47]. After an interval of $1/γ$, the infectious individuals transition to the compartment of removed cases, which contains both recovered and deceased individuals.

The epidemiological model of the four compartments (states) is given by the following system of differential equations, as follows:

$$ 
\begin{aligned}
\frac{dS(t)}{dt} &= -\frac{\beta_1 P(t) + \beta_2 I(t)}{N} S(t), \\
\frac{dP(t)}{dt} &= \frac{\beta_1 P(t) + \beta_2 I(t)}{N} S(t) - aP(t), \\
\frac{dI(t)}{dt} &= aP(t) - \gamma I(t), \\
\frac{dR(t)}{dt} &= \gamma I(t),
\end{aligned}
$$

for $t \ge 0$. Figure 1 presents a diagrammatic representation of the introduced model.

<figure><img src="image_3.png" /><figcaption>Figure 1. Representation of the SPIR model.</figcaption></figure>

According to the structure of the SPIR, the number of deceased individuals is part of the compartment of removed cases (R). Knowing the number of cases in this compartment, we can obtain an estimate of the number of deceased and recovered individuals by multiplying with $q$ and $(1-q)$, respectively, where $q$ denotes the mortality ratio. Furthermore, both symptomatic and asymptomatic carriers of the disease can be characterized as infectious individuals. Thus, for any estimate corresponding to the number of infectious individuals, we can obtain an
estimate of the corresponding number of symptomatic and asymptomatic individuals by multiplying by $r$ and $(1-r)$, respectively, where $r$ denotes the symptomatic ratio.

<table><caption>Table 1. Notation of the epidemiological model</caption><thead><tr><th>Notation</th><th>Parameter / State</th></tr></thead><tbody><tr><td><i>N</i></td><td>Total Population</td></tr><tr><td><i>S</i></td><td>Susceptible</td></tr><tr><td><i>P</i></td><td>Presymptomatic</td></tr><tr><td><i>I</i></td><td>Infectious</td></tr><tr><td><i>R</i></td><td>Removed</td></tr><tr><td>β<sub>1</sub></td><td>Infection Rate of Presymptomatic Individuals</td></tr><tr><td>β<sub>2</sub></td><td>Infection Rate of Infectious Individuals</td></tr><tr><td>γ</td><td>Infection-Free Rate</td></tr><tr><td><i>a</i></td><td>Incubation Rate</td></tr><tr><td><i>r</i></td><td>Symptomatic Ratio</td></tr><tr><td><i>q</i></td><td>Mortality Ratio</td></tr></tbody></table>

The introduced model can satisfactorily describe the dynamics of transmission of various epidemic waves, such as that of mpox [48]. Moreover, it can be adapted to many different scenarios of epidemics with minor modifications. For instance, in cases where there is insufficient evidence that presymptomatic cases contribute to disease transmission, one may set $β_1 = 0$. Furthermore, in cases where the existence of asymptomatic infections is not supported by literature, it can be assumed that the infectious compartment consists exclusively of symptomatic carriers of the disease. Table 1 explains the notation of the proposed model.

## 2.2 Description of the stochastic SPIR

Let $\mathcal{X} = \{X(t) = (S(t), P(t), I(t)), t \ge 0\}$ be a continuous-time Markov chain which describes the number of susceptible, presymptomatic and infectious cases at time $t \ge 0$. The population size $N$ is considered to be constant throughout the examined period. Hence the number of removed cases can be described as a linear combination of the components of the chain. The assumption of a closed population is made in most studies that employ Markov chains in epidemiology (see [12,16–19,21,25,49–51]).

The state space of the process is given by $S = \{(s, p, i) : 0 \le s \le N-1, 0 \le p \le N-s, 0 \le i \le N-s-p\}$. It contains a subset of disease-free states $S_A = \{(s, 0, 0) : 0 \le s \le N-1\}$, which corresponds to the termination of the outbreak, while $|S| = N(N^2 + 6N + 11)/6$. Moreover, the state space $S$ can be arranged into levels. Namely, the first level is associated with the number of susceptible cases, where $S = \bigcup_{s=0}^{N-1} l(s)$, which contains states where the number of susceptible cases is $S(t) = s$. The second level is associated with the number of presymptomatic cases, where $l(s) = \bigcup_{p=0}^{N-s} l(s, p)$, while $l(s, p) = \{(s, p, i) : 0 \le i \le N-s-p\}$. Thus, the level $l(s, p)$ contains states with $(S(t), P(t)) = (s, p)$.

The non-zero elements of the generator matrix $Q = (q_{(s,p,i),(s',p',i')}:(s,p,i), (s', p', i') \in S)$ of the chain $\mathcal{X}$, are given by
$$ q_{(s,p,i),(s',p',i')} = \begin{cases} \frac{(\beta_1 \cdot p + \beta_2 \cdot i)s}{N}, & (s', p', i') = (s-1, p+1, i), \\ a \cdot p, & (s', p', i') = (s, p-1, i+1), \\ \gamma \cdot i, & (s', p', i') = (s, p, i-1), \\ -q_{s,p,i}, & (s', p', i') = (s, p, i), \end{cases} $$

where $q_{s,p,i} = \frac{(\beta_1 \cdot p + \beta_2 \cdot i)s}{N} + a \cdot p + \gamma \cdot i$, ensuring that the rows of the generator matrix sum to zero.

## 2.3 Stochastic descriptors

Completing the description of the proposed stochastic epidemiological model, we proceed to the introduction of the stochastic characteristics under study. Specifically, formulas and algorithms are developed for calculating the distribution and moments of characteristics such as the number of infectious individuals until the extinction of the disease, the time until the occurrence of a specific number of deaths, and the number of infections caused by a presymptomatic or an infectious individual.

### 2.3.1. Total number of infectious cases

Let $Z_{s,p,i}$ be the random variable that describes the number of infectious cases until the extinction of the epidemic, given that we initiate the examination of the phenomenon from the state $(s, p, i) \in S_T = S - S_A$. By calculating the probability mass function of this random variable, the moments of any order can be easily computed. Let $p_{s,p,i}(n) = P(Z_{s,p,i} = n)$ be the probability that a number of $n$ (additional) infectious cases will occur until the extinction of the disease, starting from the current state $(s, p, i) \in S_T$. For the absorbing states, we have

$$ p_{s,0,0}(n) = \delta_{n0}, \qquad (1) $$

for $0 \le n \le N$, since $(s, 0, 0) \in S_A$, while $\delta_{ij}$ denotes the Kronecker delta function. Conditioning on the event that determines the next transition from $(s, p, i) \in S_T$ to other states, we obtain

$$ p_{s,p,i}(n) = \frac{(\beta_1 p + \beta_2 i)s}{q_{s,p,i}N} (1 - \delta_{s0}) p_{s-1,p+1,i}(n) + \frac{ap}{q_{s,p,i}} (1 - \delta_{n0}) (1 - \delta_{p0}) p_{s,p-1,i+1}(n-1) + \frac{\gamma i}{q_{s,p,i}} (1 - \delta_{i0}) p_{s,p,i-1}(n). \quad (2) $$

Based on (1) and (2) we construct Algorithm 1, which allows the computation of the probabilities $p_{s,p,i}(n)$ for $0 \le n \le N$ and $(s, p, i) \in S$.
**Algorithm 1.** Distribution of total number of infectious cases

**Require:** $N, \beta_1, \beta_2, \alpha, \gamma$

//* Absorbing States*

**Step 1:** If $(s, p, i) \in S_A$, set

$$p_{s,p,i}(n) = \delta_{n0}.$$

//* Initialization*

**Step 2:** Set $n = 0$, (a) $s = 0$, (b) $p = 0$, (c) $i = \delta_{p0}$.

//* Recursive Procedure*

**Step 3a:** Compute $p_{s,p,i}(n)$ using the formula

$$p_{s,p,i}(n) = \frac{(\beta_1 p + \beta_2 i)s}{q_{s,p,i}N} (1 - \delta_{s0})p_{s-1,p+1,i}(n) + \frac{ap}{q_{s,p,i}} (1 - \delta_{p0})(1 - \delta_{n0})p_{s,p-1,i+1}(n-1) + \frac{\gamma i}{q_{s,p,i}} (1 - \delta_{i0})p_{s,p,i-1}(n).$$

**Step 3b:** Set $i = i + 1$. While $i \le N - s - p$, go to Step 3a.

**Step 3c:** Set $p = p + 1$. While $p \le N - s$, go to Step 2c.

**Step 3d:** Set $s = s + 1$. While $s \le N - 1$, go to Step 2b.

**Step 3e:** Set $n = n + 1$. While $n \le N$, go to Step 2a.

**Remark 1.** For each initial/current state $(s, p, i) \in S_T$, the $l$-th moments for the number of infectious individuals until the extinction of the epidemic can be computed according to

$$m_{s,p,i}^l = E[Z_{s,p,i}^l] = \sum_{n=0}^{s+p} n^l p_{s,p,i}(n),$$

for $l \in \mathbb{N}$. Moreover, for $l = 1$, by multiplying $m_{s,p,i}$ by the symptomatic ratio $(r)$, we can obtain an estimate of the expected number of symptomatic individuals. In the case where $s = p = 0$, we have $m_{0,0,i}^l = 0$ for $l \in \mathbb{N}$, while $m_{s,p,i}^0 = 1$ for $(s, p, i) \in S_T$.

### 2.3.2. Occurrence time of a specific number of deaths

Another key indicator for assessing the severity of an outbreak is the occurrence (hitting) time of the first death. This event serves as an initial warning signal for public health authorities, prompting close monitoring of disease spread and evaluation of the need for preventive descriptors. In this paper, we extend this concept by investigating the occurrence time of a specific number of deaths.

According to the structure SPIR model, the number of deaths is included in the compartment of removed cases (R). Therefore, for a given value of the mortality ratio $q$, the time of occurrence of $w$ deaths corresponds to the time of occurrence of $\omega = [w/q]$ recoveries-deaths. We note that $[x]$ denotes the ceiling function, which gives the smallest integer greater than or equal to $x$.
Let $T^\omega = \inf\{t \ge 0 : R(t) = \omega\}$ be the occurrence time of $\omega$ recoveries-deaths. Additionally, let $P(T_{s,p,i}^\omega < \infty)$ be the probability that this time is finite, meaning that the occurrence of $\omega$ recoveries-deaths will take place before the epidemic is eradicated. To ensure the finiteness of the respective $l$-th moments, we employ an indicator random variable, namely

$$M_{s,p,i}^l = E \left[ T_{s,p,i}^{\omega,l} 1_{\{T_{s,p,i}^{\omega} < \infty\}} \right], \quad l \in \mathbb{N} \cup \{0\}.$$

Quantity $1_A$ denotes the indicator random variable that takes the value 1 when the event $A$ occurs and 0 otherwise. Moreover, for $l = 0$ we obtain

$$M_{s,p,i}^0 = E \left[ T_{s,p,i}^{\omega,0} 1_{\{T_{s,p,i}^{\omega} < \infty\}} \right] = E \left[ 1_{\{T_{s,p,i}^{\omega} < \infty\}} \right] = P(T_{s,p,i}^{\omega} < \infty).$$

Thus, $M_{s,p,i}^0$ equals the probability that the occurrence time of $\omega$ recoveries-deaths is finite.

Next, we construct the set $S_\omega$ that contains states with $N - s - p - i \ge \omega$. Also, the set $S_{A'}$ contains states which are associated with the absence of presymptomatic or an infectious case while the number of removed cases is less than $\omega$. Next, we define the Laplace-Stieltjes transforms conditioning on the current state $(s, p, i) \in S_T - S_\omega$, that is

$$y_{s,p,i}(z) = E \left[ e^{-zT_{s,p,i}^{\omega}} 1_{\{T_{s,p,i}^{\omega} < \infty\}} \right], \quad \mathrm{Re}(z) \ge 0.$$

The following theorem describes the recursive computation of the moments of the random variable $T_{s,p,i}^\omega$.

**Theorem 1.** The $l$-th moments of the occurrence time of a specific number of removed cases $T_{s,p,i}^\omega$, can be computed recursively using the moments of one order less, namely

$$M^l = -lQ_\omega^{-1}M^{l-1}, \quad l \in \mathbb{N}. \qquad (3)$$

Matrix $Q_\omega$ originates from $Q$ after removing rows and columns which correspond to states of the sets $S_{A'}$ and $S_\omega$.

**Proof.** Conditioning on the next transition out of the current state $(s, p, i) \in S_T - S_\omega$, we obtain

$$y_{s,p,i}(z) = \frac{(\beta_1 p + \beta_2 i)s}{(z + q_{s,p,i})N} (1 - \delta_{s0}) y_{s-1,p+1,i}(z) \\ + \frac{ap}{(z + q_{s,p,i})} (1 - \delta_{p0}) y_{s,p-1,i+1}(z) + \frac{\gamma i}{(z + q_{s,p,i})} (1 - \delta_{i0}) y_{s,p,i-1}(z). \quad (4)$$

Multiplying both sides of (4) with $(z + q_{s,p,i})$, we obtain

$$(z + q_{s,p,i})y_{s,p,i}(z) = \frac{(\beta_1 p + \beta_2 i)s}{N} (1 - \delta_{s0}) y_{s-1,p+1,i}(z) \\ + ap(1 - \delta_{p0}) y_{s,p-1,i+1}(z) + \gamma i(1 - \delta_{i0}) y_{s,p,i-1}(z),$$

or
$$ (z + q_{s,p,i})y_{s,p,i}(z) - \frac{(\beta_1 p + \beta_2 i)s}{N} (1 - \delta_{s0})y_{s-1,p+1,i}(z) - ap(1 - \delta_{p0})y_{s,p-1,i+1}(z) - \gamma i(1 - \delta_{i0})y_{s,p,i-1}(z) = 0. \quad (5) $$

Equation (5) can be written in a matrix form as

$$ (z \mathbf{I}_{|S_T - S_\omega|} - \mathbf{Q}_\omega) \mathbf{y}(z) - \mathbf{c}_{|S_T - S_\omega|} = \mathbf{0}_{|S_T - S_\omega|}. \quad (6) $$

The vector $\mathbf{c}_{|S_T - S_\omega|}$ has dimension $|S_T - S_\omega|$, with the non-zero elements located in rows corresponding to states where $N - i - p - s = \omega - 1$ and these elements are equal to $\gamma i$.

By differentiating equation (6) $l$ times with respect to $z$, multiplying by $(-1)^l$ and evaluating at $z = 0$ yields the formula for the computation of the $l$-th order moments

$$ \mathbf{M}^l = -l \mathbf{Q}_\omega^{-1} \mathbf{M}^{l-1}, \quad l \in \mathbb{N}. $$

For $l = 0$ the moments can be calculated by setting $z = 0$ in equation (6), resulting in the formula $\mathbf{M}^0 = -\mathbf{Q}_\omega^{-1}\mathbf{c}_{|S_T - S_\omega|}$, since $M_{s,p,i}^0 = y_{s,p,i}(0)$. We note that $M_{s,p,i}^l = 0$, for $l \in \mathbb{N}$, in case where $N - i - p - s \ge \omega$, since the chain has already reached the specified level $\omega$.

### 2.3.3. Number of infections caused by a presymptomatic or an infectious case

In this section, we explore the impact of presymptomatic and infectious cases at an individual level. We select a specific presymptomatic or an infectious individual from the population under study and tag them. The random variables $N^P$ and $N^I$ represent the number of new infections caused by the tagged presymptomatic or an infectious individual, correspondingly. Their analysis is significant because they provide a direct estimate of the epidemic's progression toward either extinction or expansion [50]. Let

$$ p_{s,p,i}^P(m) = P(N_{s,p,i}^P = m), \quad 0 \le m \le s, \quad (s,p,i) \in S_T, $$

for $p > 0$, be the probability mass function of the number of infections caused by the tagged presymptomatic case. Similarly, we denote by

$$ p_{s,p,i}^I(m) = P(N_{s,p,i}^I = m), \quad 0 \le m \le s, \quad (s,p,i) \in S_T, $$

for $i > 0$, the probability function of the number of infections originating from the tagged infectious case, beginning the examination of the phenomenon from state $(s, p, i) \in S_T$.

**Theorem 2.** The distributions of the number of infections caused by a presymptomatic or an infectious individual of the population can be determined using recursive relations.

**Proof.** We select a presymptomatic individual to examine its impact on the spread of the epidemic. To separate the effect of the selected individual from the rest of the presymptomatic population, we express the infection rate of susceptible individuals as follows
$$ \frac{(\beta_1 p + \beta_2 i)s}{N} = \frac{(\beta_1(p-1) + \beta_2 i)s}{N} + \frac{\beta_1 s}{N}. $$

The first term of the right-hand side corresponds to the infection of susceptible individuals from any presymptomatic or an infectious case, while the second term represents the infection rate associated with the tagged individual. A similar separation is required for the incubation rate, where

$$ ap = a(p - 1) + a, $$

as the transition of the selected presymptomatic individual to the infectious compartment marks the end of its impact on disease transmission, at least from the presymptomatic phase.

According to first-step analysis, conditioning on the first transition out of the current state $(s, p, i) \in S_T$ with $p > 0$, we obtain

$$ \begin{aligned} p_{s,p,i}^{P}(m) = & \frac{(\beta_1(p-1) + \beta_2 i)s}{q_{s,p,i} N} (1 - \delta_{s0})(1 - \delta_{ms})p_{s-1,p+1,i}^{P}(m) \\ & + \frac{\beta_1 s}{q_{s,p,i} N} (1 - \delta_{m0})(1 - \delta_{s0})p_{s-1,p+1,i}^{P}(m-1) \\ & + \frac{a(p-1)}{q_{s,p,i}} p_{s,p-1,i+1}^{P}(m) + \frac{\gamma i}{q_{s,p,i}} (1 - \delta_{i0})p_{s,p,i-1}^{P}(m) + \delta_{m0} \frac{a}{q_{s,p,i}}. \end{aligned} \quad (7) $$

It is evident that $p_{s,0,i}^P(m) = \delta_{m0}$, for $0 \le m \le s$ and $0 < s \le N - 1$ due to the absence of presymptomatic cases. The term $\frac{a}{q_{s,p,i}}$ is multiplied by the Kronecker function $\delta_{m0}$, as it signifies the exit of the selected individual from the compartment of presymptomatic cases, marking the end of its impact on disease transmission and consequently on the production of new infections.

The respective recursive formula for the computation of the probability that $m$ infections are caused by the tagged infectious case is

$$ \begin{aligned} p_{s,p,i}^{I}(m) = & \frac{(\beta_1 p + \beta_2(i-1))s}{q_{s,p,i} N} (1 - \delta_{ms})(1 - \delta_{s0})p_{s-1,p+1,i}^{I}(m) \\ & + \frac{\beta_2 s}{q_{s,p,i} N} (1 - \delta_{m0})(1 - \delta_{s0})p_{s-1,p+1,i}^{I}(m-1) \\ & + \frac{ap}{q_{s,p,i}} (1 - \delta_{p0})p_{s,p-1,i+1}^{I}(m) + \frac{\gamma(i-1)}{q_{s,p,i}} p_{s,p,i-1}^{I}(m) + \delta_{m0} \frac{\gamma}{q_{s,p,i}}, \end{aligned} \quad (8) $$

for $i > 0$. Also, $p_{s,p,0}^I(m) = \delta_{m0}$, for $0 \le m \le s$ and $0 < s \le N - 1$. It should be noted that the construction of the above relationship requires the separation of the infection rates of susceptible individuals and the exit rate from the infectious compartment, namely, $\frac{(\beta_1 p + \beta_2 i)s}{N} = \frac{(\beta_1 p + \beta_2(i-1))s}{N} + \frac{\beta_2 s}{N}$ and $\gamma i = \gamma(i - 1) + \gamma$, respectively. Furthermore, in the absence of susceptible individuals, it holds that $p_{0,p,i}^P(m) = p_{0,p,i}^I(m) = \delta_{m0}$.

We observe that both equations (7) and (8) contain terms corresponding to a difference in the number of infections produced by one order.

<figure><img src="image_4.png" /></figure>
According to the above theorem, it becomes possible to develop an algorithm for calculating the distribution (and consequently the moments) of the random variables $N_{s,p,i}^P$ and $N_{s,p,i}^I$. Algorithm 2 describes the computation of the distribution of the random variable $N_{s,p,i}^P$, while the computational process for $N_{s,p,i}^I$ shows similarities.

**Algorithm 2.** Computation of the distribution of the number of infections caused by a presymptomatic case

**Require:** $N, \beta_1, \beta_2, \alpha, \gamma$

**Step 1:** If $(s, p, i) \in S_T$ with $p = 0$ or $s = 0$, set

$$p_{s,0,i}^P(m) = p_{0,p,i}^P(m) = \delta_{m0}.$$

//* Initialization

**Step 2:** Set $m = 0$, (a) $s = 1$, (b) $p = 1$, (c) $i = 0$.

//* Recursive Procedure

**Step 3a:** Compute $p_{s,p,i}^P(m)$ using the formula

$$p_{s,p,i}^P(m) = \frac{(\beta_1(p-1) + \beta_2i)s}{q_{s,p,i} N} (1 - \delta_{ms}) p_{s-1,p+1,i}^P(m) \\ + \frac{\beta_1 s}{q_{s,p,i} N} (1 - \delta_{m0}) p_{s-1,p+1,i}^P(m-1) \\ + \frac{a(p-1)}{q_{s,p,i}} p_{s,p-1,i+1}^P(m) + \frac{\gamma i}{q_{s,p,i}} (1 - \delta_{i0}) p_{s,p,i-1}^P(m) + \delta_{m0} \frac{a}{q_{s,p,i}}.$$

**Step 3b:** Set $i = i + 1$. While $i \le N - s - p$, go to Step 3a.

**Step 3c:** Set $p = p + 1$. While $p \le N - s$, go to Step 2c.

**Step 3d:** Set $s = s + 1$. While $s \le N - 1$, go to Step 2b.

**Step 3e:** Set $m = m + 1$. While $m \le N - 1$, go to Step 2a.

## 2.4 Estimation of epidemiological parameters using particle filter

Deterministic methods cannot account for the inherent randomness in epidemics. Standard compartmental models overlook the randomness associated with environmental factors, the age distribution of the population, or characteristics related to the immune system of infected individuals. Attempting to incorporate all this information into a system of differential equations would result in models with high complexity which are prone to overfitting. Additionally, measurements are subject to noise which aligns with the unpredictable nature of epidemics. Therefore, an online parameter estimation process is crucial to enhance the accuracy of stochastic descriptors, such as those discussed in the previous section.

From a statistical point of view, our objective is to estimate simultaneously the parameters of the model denoted by $\theta \in \Theta$ and the unobserved state vector $\mathbf{x}_{0:K} = (x_0, \dots, x_K)$ by means of the sequence of measurements $\mathbf{y}_{1:K} = (y_1, \dots, y_K)$. Concerning the parameter inference, in a frequentist approach, the maximum likelihood estimator of $\theta$ could be obtained by means of an expectation-maximization (EM) algorithm [52]. In a Bayesian approach, the target distribution is rather the posterior distribution of the parameters. Concerning the inference of the states and for a given parameter vector $\theta$, our purpose is to estimate the model's hidden
states. To do so, we employ a recursive Bayesian filter which is based on Monte-Carlo simulations [53].

The main idea of particle filtering (PF) is to use importance sampling to address the fact that the analytical representation of the density $p_{\theta}(\mathbf{x}_{0:k}|\mathbf{y}_{1:k})$ is unknown. To be more precise, since the likelihood depends on hidden states, it is intractable, which means that sampling from the posterior itself is infeasible. Hence, we sample from the so-called importance (or proposal) density function, which describes the transition law between the system's hidden states and is further denoted by $f_{\theta}(\mathbf{x}_{0:k}|\mathbf{y}_{1:k})$.

For our purposes, we have to define a state space model, where $p_{\theta}(\mathbf{x}_k|\mathbf{x}_{k-1})$ represents the state equation which dictates the progression of the vector of hidden states with respect to time, while $p_{\theta}(\mathbf{y}_k|\mathbf{x}_k)$ represents the measurement equation which connects the measurements with the hidden states $\mathbf{x}_k$.

For a given parameter vector $\theta$, the importance weights are defined and updated according to

$$w_k^{(i)} = \frac{p(\mathbf{x}_{0:k}^{(i)} | \mathbf{y}_{1:k})}{f(\mathbf{x}_{0:k}^{(i)} | \mathbf{y}_{1:k})} \propto \frac{p(\mathbf{y}_k | \mathbf{x}_k^{(i)}) p(\mathbf{x}_k^{(i)} | \mathbf{x}_{k-1}^{(i)}) p(\mathbf{x}_{0:k-1}^{(i)} | \mathbf{y}_{1:k-1})}{f(\mathbf{x}_k^{(i)} | \mathbf{x}_{0:k-1}^{(i)}, \mathbf{y}_{1:k}) f(\mathbf{x}_{0:k-1}^{(i)} | \mathbf{y}_{1:k-1})} \\ = w_{k-1}^{(i)} \frac{p(\mathbf{y}_k | \mathbf{x}_k^{(i)}) p(\mathbf{x}_k^{(i)} | \mathbf{x}_{k-1}^{(i)})}{f(\mathbf{x}_k^{(i)} | \mathbf{x}_{0:k-1}^{(i)}, \mathbf{y}_{1:k})},$$

for $k = 1, ..., K$, while index $(i)$ corresponds to the $i$-th particle. All the terms appearing in the above expression depend on $\theta$; however, we omit the dependence on $\theta$ in the notation for simplicity.

We further consider that the importance density function $f(\mathbf{x}_k^{(i)} | \mathbf{x}_{0:k-1}^{(i)}, \mathbf{y}_{1:k})$ is equivalent to the $f(\mathbf{x}_k^{(i)} | \mathbf{x}_{k-1}^{(i)}, \mathbf{y}_k)$, considering that the hidden process is a first-order Markov chain. Then the weights can be recursively computed based on

$$w_k^{(i)} \propto w_{k-1}^{(i)} \frac{p(\mathbf{y}_k | \mathbf{x}_k^{(i)}) p(\mathbf{x}_k^{(i)} | \mathbf{x}_{k-1}^{(i)})}{f(\mathbf{x}_k^{(i)} | \mathbf{x}_{k-1}^{(i)}, \mathbf{y}_k)} \quad (9)$$

The weights of equation (9) are then normalized. We notice that each weight is proportional to $p(\mathbf{y}_k | \mathbf{x}_k^{(i)})$. Consequently, the procedure retains the particles that are most likely to have generated the measurements [41].

Next, we specify the importance density function $f(\mathbf{x}_k^{(i)} | \mathbf{x}_{k-1}^{(i)}, \mathbf{y}_k)$, where binomial distributions are used since we refer to transitions of individuals between the compartments of the SPIR at each time step. Let $(d_{SP})_{k-1}^{(i)}$ denote the number of individuals that move from the compartment of susceptible to that of presymptomatic cases for the $i$-th particle at time step $k$. The interpretation of $(d_{PI})_{k-1}^{(i)}$ and $(d_{IR})_{k-1}^{(i)}$ is similar. Specifically,

$$(d_{SP})_{k-1}^{(i)} \sim B \left( S_{k-1}^{(i)}, 1 - e^{-\frac{\beta_1 P_{k-1}^{(i)} + \beta_2 I_{k-1}^{(i)}}{N}} \right), \quad (10)$$

$$(d_{PI})_{k-1}^{(i)} \sim B(P_{k-1}^{(i)}, 1 - e^{-\alpha}), \quad (11)$$
and

$$ (d_{IR})_{k-1}^{(i)} \sim B(I_{k-1}^{(i)}, 1 - e^{-\gamma}), \quad (12) $$

where the second parameter of (10)–(12) represents the probability of an individual transitioning from the susceptible to presymptomatic, presymptomatic to infectious and infectious to removed compartment. Thus, we introduce a state space representation, with state equation

$$ \mathbf{x}_{k}^{(i)} = \mathbf{x}_{k-1}^{(i)} + \mathbf{h}(\mathbf{x}_{k-1}^{(i)}; \boldsymbol{\theta}), \quad (13) $$

where $\mathbf{x}_k^{(i)} = (P_k^{(i)}, I_k^{(i)}, R_k^{(i)})^T$ and

$$ \mathbf{h}(\mathbf{x}_{k-1}^{(i)}; \boldsymbol{\theta}) = \begin{pmatrix} (d_{SP})_{k-1}^{(i)} - (d_{Pl})_{k-1}^{(i)} \\ (d_{Pl})_{k-1}^{(i)} - (d_{IR})_{k-1}^{(i)} \\ (d_{IR})_{k-1}^{(i)} \end{pmatrix}, $$

for $k = 1, ..., K$ and $i = 1, ..., M$. Equation (13) represents a discrete-time stochastic equivalent of the system of differential equations corresponding to the SPIR model. According to (13), the importance function for the $i$-th particle $f_{\theta}(\mathbf{x}_k^{(i)}|\mathbf{x}_{k-1}^{(i)}, \mathbf{y}_k)$, takes the form of the transition probability distribution $p(\mathbf{x}_k^{(i)}|\mathbf{x}_{k-1}^{(i)})$.

Usually, vector $\boldsymbol{\theta}$ is considered to be constant, that is $\boldsymbol{\theta}_k^{(i)} = \boldsymbol{\theta}^{(i)}$, for any $k = 1, ..., K$ and $i = 1, ..., M$. However, the constant dynamics generally lead to the divergence of particle filters. Focusing on the generation of time-varying parameter estimates, we augment the state space by adding the epidemic's infection rates and the mortality ratio [54,55]. We note that $(\beta_{1,k}, \beta_{2,k}, q_k) \in \mathbb{R}_+^2 \times (0,1)$. Thus, for the evolution of the infection rates we select truncated normal distributions depending on the previous value of the parameters, as

$$ \beta_{1,k}^{(i)} \sim g_1(\beta_{1,k-1}^{(i)}) = tN(\beta_{1,k-1}^{(i)}, \sigma_1, 0, \infty) \quad (14) $$

and

$$ \beta_{2,k}^{(i)} \sim g_2(\beta_{2,k-1}^{(i)}) = tN(\beta_{2,k-1}^{(i)}, \sigma_2, 0, \infty). \quad (15) $$

This approach allows $\beta_{1,k}^{(i)}$ and $\beta_{2,k}^{(i)}$ to evolve without enforcing any predetermined trends, as their paths are fully determined by the model's fit to the measurements. It should be noted that the introduction of noise into the evolution of the time-varying parameters mitigates the problem of particle impoverishment, which commonly arises when using particle filters [56]. For the mortality ratio we utilize a truncated normal distribution aiming to produce estimates in (0,1). Namely,

$$ q_k^{(i)} \sim g_3(q_{k-1}^{(i)}) = tN(q_{k-1}^{(i)}, \sigma_3, 0, 1). \quad (16) $$

In equations (14)–(16), the first term within the truncated normal distributions corresponds to the mean, the second to the standard deviation and the last two terms define the range of values of the random variable.
For brevity, the densities dictating the parameters' progression can be assembled in a vector $\mathbf{g}(\boldsymbol{\theta}_{k-1}^{(i)}) = (g_1(\beta_{1,k-1}^{(i)}), g_2(\beta_{2,k-1}^{(i)}), g_3(q_{k-1}^{(i)}))^T$. This vector can be extended with additional epidemiological rates if the integration of more time-varying parameters is required by the characteristics of the studied epidemic. The state and parameter estimates for the SPIR model can be obtained simultaneously by using an augmented state vector $\tilde{\mathbf{x}}_k^{(i)} = \text{vec}(\mathbf{x}_k^{(i)}, \boldsymbol{\theta}_k^{(i)})$ [55]. Considering equations (14)–(16), the importance function $f_{\theta}(\tilde{\mathbf{x}}_k^{(i)}|\tilde{\mathbf{x}}_{k-1}^{(i)}, \mathbf{y}_k)$ takes the form of the transition probability distribution $p(\tilde{\mathbf{x}}_k^{(i)}|\tilde{\mathbf{x}}_{k-1}^{(i)})$ [44], resulting in

$$w_k^{(i)} \propto w_{k-1}^{(i)} \frac{p(\mathbf{y}_k | \tilde{\mathbf{x}}_k^{(i)}) p(\tilde{\mathbf{x}}_k^{(i)} | \tilde{\mathbf{x}}_{k-1}^{(i)})}{f_{\theta}(\tilde{\mathbf{x}}_k^{(i)} | \tilde{\mathbf{x}}_{k-1}^{(i)}, \mathbf{y}_k)} = w_{k-1}^{(i)} p(\mathbf{y}_k | \tilde{\mathbf{x}}_k^{(i)}).$$

Next, we proceed to the determination of the measurement function, denoted by $p(\mathbf{y}_k|\tilde{\mathbf{x}}_k^{(i)})$. In the literature, measurements such as new infections have been considered as arrivals of a Poisson process [57–59]. Calvetti et al. [58] employed a univariate Poisson distribution, with its rate based on the estimated number of new cases at each time step as obtained from the particle filter. A similar measurement function is used by Childs et al. [59] connecting the estimated with the observed number of deaths.

Aiming to utilize the information provided by the recorded (symptomatic) infections and deaths, we choose as the measurement function the product of two Poisson distributions, given by the equation

$$p(\mathbf{y}_k|\tilde{\mathbf{x}}_k) = \prod_{j=1}^{2} Pois(y_{k,j}; \lambda_{k,j}), \qquad (17)$$

for $k = 1, ..., K$. The terms $\lambda_{k,j}$ represent the parameters of the Poisson distributions, whereas the terms $y_{k,j}$ denote the measurements at time step $k$. Furthermore, we choose as Poisson parameters, the estimated (via the particle filter) number of new symptomatic infections ($\lambda_{k,1}$) and deaths ($\lambda_{k,2}$) at each time step $k$. The first parameter can be calculated as the product of the estimated number of new infectious individuals times the ratio of symptomatic infections ($r$), while the second as the product of the estimated number of new recoveries-deaths times the mortality ratio ($q_k$). As a result, the parameters of the Poisson in (17) are obtained through the components of the vector $\tilde{\mathbf{x}}_k$.

Finally, we employ a resampling step to eliminate the weight degeneracy issue [43,60]. The decision on whether the existing particles should be resampled is based on the effective sample size measure of degeneracy

$$M_{\text{eff}}(k) = \frac{M}{1 + Var_{p(\cdot|\mathbf{y}_{1:k})}(w(\tilde{\mathbf{x}}_k))}.$$

As this quantity cannot be computed directly, we rely on the estimate [56,60]

$$\tilde{M}_{\text{eff}}(k) = \frac{1}{\sum_{i=1}^{M} (w_k^{(i)})^2}.$$

A threshold value must be selected for the resampling step, which is typically set as a percentage of the predetermined number of particles. We set $M_{th} = 0.75M$, as a relatively
frequent resampling is deemed appropriate to eliminate degeneracy issues [41,42,60]. Based on the sampling-resampling procedure, both states and parameters are dynamically estimated considering the evolution of recorded cases and deaths. Algorithm 3 shows the operation of the aforementioned estimation procedure.

<caption>**Algorithm 3.** The augmented particle filter with resampling</caption>

**Require:** $\mathbf{y}_{1:k}, g(\cdot), h(\cdot), p_{y_k|\tilde{x}_k}, K, M, M_{th}$

// * Initialization
**For** $i = 1, 2, ..., M$ **do**
  Sample: $\mathbf{x}_0^{(i)} \sim p_0(\mathbf{x}_0)$.
  Sample: $\boldsymbol{\theta}_0^{(i)} \sim p_0(\boldsymbol{\theta}_0)$.
  Set: $\tilde{\mathbf{x}}_0^{(i)} = \text{vec}(\mathbf{x}_0^{(i)}, \boldsymbol{\theta}_0^{(i)})$ and $\tilde{w}_0^{(i)} = \frac{1}{M}$.
**end**
**For** $k = 1, ..., K$ **do**
  // * Importance Sampling
  **For** $i = 1, 2, ..., M$ **do**
    Sample: $\mathbf{x}_k^{(i)} = \mathbf{x}_{k-1}^{(i)} + \mathbf{h}(\tilde{\mathbf{x}}_{k-1}^{(i)})$.
    Sample: $\boldsymbol{\theta}_k^{(i)} \sim \mathbf{g}(\boldsymbol{\theta}_{k-1}^{(i)})$.
    Set: $\tilde{\mathbf{x}}_k^{(i)} = \text{vec}(\mathbf{x}_k^{(i)}, \boldsymbol{\theta}_k^{(i)})$.
    Compute weights: $w_k^{(i)} = \tilde{w}_{k-1}^{(i)} p(\mathbf{y}_k | \tilde{\mathbf{x}}_k^{(i)})$.
**end**
**For** $i = 1, 2, ..., M$ **do**
  Normalize weights: $\tilde{w}_k^{(i)} = \frac{w_k^{(i)}}{\sum_{j=1}^{M} w_k^{(j)}}$.
**end**
// * Resampling (if necessary)
**If** $\hat{M}_{eff}(k) < M_{th}$
  **For** $i = 1, 2, ..., M$ **do**
    Resample with replacement new $M$ particles using: $p(j(i) = d) = \tilde{w}_k^{(d)}$, $d = 1, ..., M$.
    Set: $\tilde{\mathbf{x}}_k^{(i)} = \tilde{\mathbf{x}}_k^{(j(i))}$.
    Set: $\tilde{w}_k^{(i)} = \frac{1}{M}$.
**end**
**end**
**end**

**3. Sensitivity analysis and application to real data**

In this section, we first study the behavior of the number of infections caused by a presymptomatic or an infectious individual as a function of the model's parameters. Next, we examine the efficiency of the introduced estimation methodology, which considers the
existence of time-varying epidemiological parameters, based on the data of cases (symptomatic infections) and deaths from the mpox outbreak of 2022 in Ghana. The efficiency of the proposed time-varying method is compared with that of the “standard” method, which considers constant parameters for the entire duration of the epidemic.

**3.1 Sensitivity analysis for the number of infections caused by a presymptomatic or an infectious individual**

In this section, we examine the behavior of the number of infections caused by a presymptomatic or an infectious individual of the population, as a function of the model’s epidemiological parameters. In Figure 2, using the formulas from Section 2.3.3, the distribution of the number of infections caused by a presymptomatic individual is presented for $β_1 \in \{0.5, 1, 2\}$, $β_2 \in \{1, 2\}$ and $α \in \{0.1, 0.25\}$. The population size was set to $N = 100$, the exit rate from the infectious phase to $γ = 0.1$ and the initial state was defined as (90, 5, 5).

According to Figure 2, the distribution of the number of infections exhibits both decreasing and unimodal forms. The decreasing forms gradually transform into unimodal as the infection rate of presymptomatic cases ($β_1$) increases or the incubation rate ($α$) decreases. The expected value of the number of infections increases as a function of the infection rate ($β_1$), as the impact of presymptomatic individuals on disease spread intensifies.

Additionally, the expected value shows a decreasing trend as we move to higher values of $α$, because presymptomatic individuals transition more quickly to the infectious compartment, participating in disease transmission for a shorter period (at least from the presymptomatic phase). Furthermore, higher transmission rates of infectious individuals ($β_2$) lead to a decrease in the expected value, as the probability of susceptible individuals being infected by infectious cases increases. Finally, the highest expected value for the number of infections, based on the selected parameter combinations, corresponds to the triplet ($β_1, β_2, α$) = (2, 1, 0.1), while the lowest value is found for ($β_1, β_2, α$) = (0.5, 2, 0.25), where a presymptomatic individual is expected to cause 0.82 infections.

Figure 3 shows the distributions of the number of infections caused by an infectious individual of the population for $β_1 \in \{0.5, 1\}$, $β_2 \in \{0.5, 1, 2\}$ and $γ \in \{0.1, 0.25\}$. The population size was chosen to be $N = 100$, with an incubation rate $α = 0.1$. The initial state was set to (90, 5, 5). Distributions with decreasing forms are observed, which are replaced by unimodal ones, as the transmission rate of infectious individuals increases ($β_2$) or the infection-free ($γ$) decreases.

The expected value of the number of infections increases as a function of the transmission rate ($β_2$), since the impact of infectious individuals on the spread of the disease becomes more pronounced. Moreover, the expected value exhibits a decreasing trend for higher values of ($γ$), as infectious individuals exit this phase more quickly, affecting the spread of the disease for shorter periods. Higher transmission rates ($β_1$) lead to a decrease in the expected value because the probability of susceptible individuals being infected by presymptomatic individuals of the population increases. The highest expected value of the number of infections based on the selected parameter combinations corresponds to the triplet ($β_1, β_2, γ$) = (0.5, 2, 0.1), while the lowest corresponds to the triplet ($β_1, β_2, γ$) = (1, 0.5, 0.25) with a value of 0.9.
<figure><img src="image_5.png" /><figcaption>Figure 2. Number of infections caused by a presymptomatic individual of the population ($N^P$) for $\beta_1 \in \{0.5, 1, 2\}$, $\beta_2 \in \{1, 2\}$ and $\alpha \in \{0.1, 0.25\}$.</figcaption></figure>
<figure><img src="image_6.png" /><figcaption>Figure 3. Number of infections caused by an infectious individual of the population ($N^I$) for $\beta_1 \in \{0.5, 1\}$, $\beta_2 \in \{0.5, 1, 2\}$ and $\gamma \in \{0.1, 0.25\}$.</figcaption></figure>

## 3.2 Comparison of the estimation methods in terms of performance

In this section, we compare the results obtained by the standard approach, where the parameters are considered fixed, with those obtained by the proposed method, where the parameters are dynamically updated via the particle filter. The application of the particle filter requires the initialization of the states and parameters $\theta = (\beta_1, \beta_2, \alpha, \gamma, q, r)^T$, using prior distributions or point estimates. The weekly time scale was chosen to estimate the model
parameters, since the records of mpox cases and deaths in Ghana, as presented in the database, are updated on a weekly basis [61]. According to [62], the average incubation time of the disease is about 2 weeks. Thus, we choose $\hat{a} = 1/2$. Furthermore, we set $\hat{r} = 0.961$, as 1.3%–6.5% of individuals infected with mpox do not exhibit symptoms [63], leading to a symptomatic ratio of 93.5%–98.7%. As the mortality ratio from the disease ranges between 0% and 11%, we choose $\hat{q}_0 \sim \text{Unif}(0, 0.11)$ to initialize this parameter.

Based on the data, for the initial number of infectious individuals we set $\hat{l}_0 = 2$, while $\hat{R}_0 = 0$ as we refer to the beginning of the outbreak. In addition, there are no records available for the initial estimate of the number of presymptomatic individuals ($P_0$). Furthermore, the initial number of susceptible individuals ($S_0$) should be estimated, while it should be noted that mpox infected a small part of the country's population [64]. Furthermore, there is no clear range of values that can be used to initialize parameters, such as the transmission and infection-free rates. Therefore, their values have to be estimated. Using the least squares method we obtain $\hat{S}_0 = 118$, $\hat{P}_0 = 9$, $\hat{\gamma} = 0.24$, $\hat{\beta}_{1,0} = 0.37$ and $\hat{\beta}_{2,0} = 0.49$. These values are also utilized during the implementation of the standard approach. Poisson distributions are chosen for the initializations of the states. According to the above, for the prior distributions of the particle filter we set $\hat{P}_0 \sim \text{Pois}(9)$ and $\hat{l}_0 \sim \text{Pois}(2)$, while for the transmission rates we use uniform distributions, where $\hat{\beta}_{1,0} \sim \text{Unif}(0.24, 0.49)$ and $\hat{\beta}_{2,0} \sim \text{Unif}(0.32, 0.66)$.

Next, we estimate the standard deviations of the truncated normal distributions used for the epidemiological parameter trajectories (14)–(16). We use the estimate of the maximum log-likelihood according to the formula presented by Chen et al. (2003) [44]. The maximum value was given for the combination ($\sigma_1, \sigma_2, \sigma_3$) = (0.1, 0.25, 0.006). Figure 4 presents the estimated number of new cases and deaths obtained using the augmented particle filter with resampling, along with the 95% credible intervals, using $M = 2000$ particles. We note that increasing the number of particles produces similar estimates but at a higher computational cost. Therefore, using 2000 particles appears to be an appropriate choice for our case.

<figure><img src="image_7.png" /><figcaption>Figure 4. Estimates for the number of new cases (new symptomatic infectious cases) and deaths using the particle filter. The black dots correspond to the recorded cases and deaths from mpox in Ghana, while the red lines are the mean posterior estimates of the filter along with the 95% credible intervals.</figcaption></figure>
Next, the efficiency of the introduced time-varying method is examined, and comparisons
are made with the standard method which assumes constant dynamics. In the case of mpox (as
well as in many other cases of epidemics), the available records are often limited to the number
of cases and deaths. As a result, the evaluation of the estimation capability of the two
approaches was based on the following two stochastic characteristics: the number of infectious
cases until the extinction of the epidemic and the occurrence time of the first death.

In Figure 5 the number of infections (symptomatic infectious individuals) recorded until the extinction of the disease is shown in black. Additionally, the estimates of the proposed and the standard method, are presented in blue and orange, correspondingly. For the computation of the expected values, we employ the formula introduced in Section 2.3.1. Both methods demonstrate satisfactory accuracy up to the 5th week and from the 27th week onward, when the number of new infections decreases significantly. However, the expected values produced by the standard method show notable deviations between the 5th and 27th weeks, corresponding to underestimations of the number of infections leading up to the disease’s extinction. The expected values corresponding to the proposed method show smaller deviations from the actual infection counts, with these estimates appearing satisfactory up to the 20th week. Underestimations are also observed in the case of the proposed method from the 20th to the 26th week. However, the deviations of the estimates during this period are smaller compared to the standard approach.

Specifically, the expected values obtained by means of the proposed method correspond to a root mean squared error (RMSE) of 5.480, whereas the standard approach corresponds to an RMSE of 12.497. As a result, the estimates of the standard method lead to an increase of the RMSE by 128.047%. Based on above results, and according to the standard method, the epidemic wave will last for a shorter period, as the estimated number of remaining infections until the end of the outbreak decreases more rapidly compared to the proposed method.

The two methods are further compared in terms of the occurrence time of the first death. Using equation (3), estimates are examined up to the 9th week when the first death was recorded. The standard method underestimates this period, estimating that the first death will occur sooner (Figure 6). On the other hand, the estimates of the proposed method improve after the 2nd week. In particular, the proposed method leads to an RMSE of 1.588, whereas the standard method gives an RMSE of 2.082, which corresponds to an increase of 31.108%.

We note that for the estimation of the studied stochastic descriptors, precise estimates of the model's states and parameters at each time step are required. These values originate from the mean posterior estimates of the introduced time-varying method. Consequently, we state that the improved precision provided by the proposed method can be also attributed to the denoised state estimates.
<figure><img src="image_8.png" /><figcaption>Figure 5. Expected number of infections until the disease's extinction. The black diamonds represent the progression of the actual number of infections, while the orange circles and blue triangles indicate the expected values obtained using the standard and proposed methods, respectively.</figcaption></figure>

<figure><img src="image_9.png" /><figcaption>Figure 6. Expected time until the first death. The black diamonds represent the number of time steps (weeks) remaining until the occurrence of the first death. The orange circles and blue triangles indicate the expected values obtained using the standard and proposed methods, respectively.</figcaption></figure>

### 3.3 Impact of presymptomatic and infectious individuals on the spread of the disease

In this section, we examine the impact of presymptomatic and infectious cases on the spread of the disease through the stochastic descriptors presented in Section 2.3.3. Specifically, we explore the random variables corresponding to the number of infections caused by a presymptomatic or an infectious individual. It is noted that these characteristics play a central role in the field of epidemiology, as their values serve as indicators for estimating the possible course of the disease in terms of its spread or elimination [50]. Initially, we compare the
evolution of these indicators, taking into account the parameters estimated through the introduced method based on mpox data.

To obtain the expected values displayed in Figure 7, we choose as parameters the mean posterior estimates derived from the proposed method. In the case of Ghana, it is observed that the impact of infectious individuals on disease transmission is more pronounced compared to that of presymptomatic cases (Figure 7). Specifically, in the first week, an infectious individual is expected to cause approximately 1.25 infections, while the impact of a presymptomatic individual is significantly lower, as it is expected to cause about 0.5 infections. This behavior is attributed to the higher transmission rate of infectious ($β_2$) compared to that of presymptomatic individuals ($β_1$), combined with the fact that carriers of the disease remain longer in the infectious rather than in the presymptomatic phase, as indicated by the estimated values of the parameters $α = 0.5$ and $γ = 0.24$. However, as will be demonstrated later, the impact of presymptomatic individuals should not be considered negligible. Finally, similar behavior is observed in the following weeks, with the impact of infectious cases on disease spread being more significant compared to that of presymptomatic individuals.

In Figures 8 and 9, the trends of the expected number of infections caused by a presymptomatic or an infectious individual, respectively, are compared between the two methods. The estimates of the two indicators are shown in red for the standard and in blue for the proposed method. In the case of the expected number of infections caused by a presymptomatic individual, the two curves follow a similar trajectory (Figure 8). In the initial weeks, according to the standard approach, the impact of presymptomatic individuals on disease transmission is slightly greater. However, from the 7th week onward, the impact of presymptomatic individuals appears to be stronger under the proposed method.

<figure><img src="image_10.png" /><figcaption>Figure 7. Expected number of infections caused by a presymptomatic ($N^P$) or an infectious individual ($N^I$) according to the proposed method.</figcaption></figure>

Regarding the trends of the expected number of infections caused by an infectious individual (Figure 9), according to the standard method, the expected number declines more rapidly compared to the method with time-varying parameters. The increased values observed from the 3rd week onward in the case of the proposed method indicate that the impact of infectious individuals on the spread of the disease remains significant for longer periods, resulting in an extended epidemic duration. As a result, the findings of the proposed method
suggest that the spread of the disease remains non-negligible for a longer period. This information can be used to organize timely and targeted containment and control measures for the epidemic's spread, such as contact tracing of symptomatic carriers and the expansion of diagnostic testing.

<figure><img src="image_11.png" /><figcaption>Figure 8. Expected number of infections caused by a presymptomatic individual ($N^P$) derived from the standard and the proposed method.</figcaption></figure>

<figure><img src="image_12.png" /><figcaption>Figure 9. Expected number of infections caused by an infectious individual ($N^I$) derived from the standard and the proposed method.</figcaption></figure>

In the case of mpox in Ghana, the impact of presymptomatic individuals appears to be less significant compared to that of infectious individuals. This can be attributed to the lower transmission rate of individuals in the first group and the difference in the expected duration of the presymptomatic and infectious phases. However, in other epidemics disease carriers may exhibit higher transmission rates during the presymptomatic phase [65,66], making it important to monitor the evolution of $N^P$. In Figure 10, a theoretical example is presented to examine the
above claim using constant epidemiological parameters. Specifically, we choose $\beta_1 = \beta_2 = 0.25$ and $\gamma = \alpha = 0.10$. During the first week, a presymptomatic individual is expected to cause 1.4 infections, while the expected value of the indicator is close to 0.2 in the 17th week. Consequently, it is evident that presymptomatic individuals can play a critical role in initiating or amplifying epidemiological outbreaks. This measurement highlights the importance of limiting their impact, particularly due to the challenges arising from the difficulty of timely identification of such cases. Measures such as monitoring and isolating suspected cases or contacts of symptomatic individuals could help reduce the spread of the disease [67].

<figure><img src="image_13.png" /><figcaption>Figure 10. Expected number of infections caused by a presymptomatic individual ($N^P$) assuming constant epidemiological parameters.</figcaption></figure>

**4. Concluding remarks**

This paper introduces a novel approach that combines efficient algorithms for estimating stochastic descriptors, which are critical for assessing epidemic severity, with an online parameter estimation framework. A SPIR model is presented that captures the key dynamics of most epidemics, while incorporating several fundamental epidemic characteristics without significantly increasing the computational burden. The proposed methodology enables the online estimation of stochastic descriptors by incorporating time-varying parameters estimated through a particle filtering procedure.

Four stochastic descriptors are investigated, representing size-, time-, and individual-related characteristics of epidemics. These are the total number of infectious individuals until the extinction of the outbreak, the occurrence time of a specific number of deaths and the number of infections caused by a presymptomatic or an infectious case. The first two descriptors are used to compare the proposed hybrid method with the standard Markovian approach which considers constant dynamics, validated through mpox data from Ghana. The other two descriptors introduce novel stochastic indicators, providing insights into the spread of the disease and offering valuable information for public health strategies. Alongside the parameter update mechanism originating from particle filter, the introduced methodology is well-suited for generating accurate online predictions, offering a significant advantage over traditional models that rely on fixed dynamics.
The performance of the particle filter has been evaluated by incorporating additional epidemiological rates, such as the incubation ($\alpha$) and infection-free ($\gamma$) rates, too. Additional time-varying rates could have made the model more flexible, even though fluctuating incubation and infection-free rates are not supported by the literature in the case of mpox [41]. Although adding these parameters to the augmented state vector slightly improved the particle filter's fitting performance, the precision of the examined stochastic descriptors was inferior, sometimes even falling behind the approach which includes the estimation of $\beta_1, \beta_2$ and $q$. This suggests that more complex models are not necessarily more effective, reinforcing the idea that simpler approaches tend to be more robust and less prone to overfitting. Furthermore, we also experimented with using negative binomial distributions for the measurement model. This selection did not improve the estimation performance, while requiring the estimation of two additional size parameters (one for each negative binomial distribution). Therefore, we ultimately selected Poisson distributions.

As mentioned, estimating the number of infectious and recovered-deceased individuals enables the estimation of the number of symptomatic and deceased individuals by multiplying with the corresponding ratios of symptomatic infections and mortality, respectively. As a result, estimates for the number of individuals in these two categories can be obtained without the need to use compartmental Markov chains with additional components, which would otherwise increase the computational cost involved in determining the distributions and moments under study. Furthermore, the introduced model's structure allows for the development of recursive relations and algorithms that, in the majority of stochastic descriptors, do not require the use of matrix-based representations; an approach typically deemed necessary in cases where population members may undergo multiple cycles (transitions) between the compartments of the model [9, 12, 16, 49, 68-71].

In general, the effectiveness of the standard method decreases in cases of epidemics with significant time-varying dynamics. According to the results presented in Section 3.2, the standard method leads to increased RMSE values by 128.047% and 31.108% for the number of infectious individuals until the epidemic's extinction and the time of the first death, respectively. Nevertheless, we emphasize the efficiency of the standard approach, which has been widely used in the literature for epidemic modeling [16, 18, 21, 50, 70, 72, 73]. Its performance is quite satisfactory for epidemics characterized by low temporal variability. Moreover, it remains useful in epidemic scenarios where rapid estimation is preferred or where the assumption of an outbreak with constant parameters is considered adequate.

In future work, it would be interesting to explore the use of artificial intelligence for estimating stochastic descriptors and to compare its performance with the methodology presented in this paper. When epidemic dynamics follow distinct phases, parameter estimation can focus on these periods, reducing computational and storage demands. Change-point analysis can help identify such phases [74-76]. Moreover, likelihood-free methods such as the Approximate Bayesian methods [77] could be employed to estimate the approximate posterior distribution of the parameters in case where more flexibility in the structure of the model is requested. From a methodological point of view, other augmentation schemes could be used for statistical inference such as the one presented in [78].

## Declaration of Competing Interests

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.
## Data Availability Statement

The utilized dataset exists in the “Our World in Data” Monkeypox GitHub repository and can be accessed through the link https://raw.github usercontent.com/owid/monkeypox/2088d25b06d19276d07e7182870d771faf bea450/owid-monkeypox-data.csv

## Funding

The authors acknowledge the support of the French National Research Agency (ANR), under grant ANR-21-CE40-005 (project HSMM-INCA). The authors would also like to thank the A.G. Leventis Foundation for its funding support.

## References

1. Papageorgiou VE, Vasiliadis G, Tsaklidis G. Analyzing the Asymptotic Behavior of an Extended SEIR Model with Vaccination for COVID-19. Mathematics 2023;12:55.

2. Britton T. Stochastic epidemic models: A survey. Math Biosci 2010;225:24-35.

3. Brauer F, Castillo-Chavez C, Feng Z. Endemic Disease Models. 2019. page 63-116.

4. Vargas-De-León C. On the global stability of SIS, SIR and SIRS epidemic models with standard incidence. Chaos Solitons Fractals 2011;44:1106-10.

5. Malkov E. Simulation of coronavirus disease 2019 (COVID-19) scenarios with possibility of reinfection. Chaos Solitons Fractals 2020;139.

6. Korolev I. Identification and estimation of the SEIRD epidemic model for COVID-19. J Econom 2021;220:63-85.

7. Gamboa M, López-García M, López-Herrero MJ. A Stochastic SVIR Model with Imperfect Vaccine and External Source of Infection. In: Lecture Notes in Computer Science (including subseries Lecture Notes in Artificial Intelligence and Lecture Notes in Bioinformatics). Springer Science and Business Media Deutschland GmbH; 2021. page 197-209.

8. Tuite AR, Fisman DN, Odutayo A, Bobos P, Allen V, Bogoch II, et al. COVID-19 Hospitalizations, ICU Admissions and Deaths Associated with the New Variants of Concern. 2021.

9. Amador J, Gómez-Corral A. A stochastic epidemic model with two quarantine states and limited carrying capacity for quarantine. Physica A: Statistical Mechanics and its Applications 2020;544:121899.

10. Marshall P. The impact of quarantine on Covid-19 infections. Epidemiol Methods 2021;10.

11. Mohsen AA, AL-Husseiny HF, Zhou X, Hattaf K. Global stability of COVID-19 model involving the quarantine strategy and media coverage effects. AIMS Public Health 2020;7:587-605.
12. Papageorgiou VE. Novel stochastic descriptors of a Markovian SIRD model for the assessment of the severity behind epidemic outbreaks. J Franklin Inst 2024;107022.

13. Papageorgiou VE. Commentary on "Stochastic modeling of computer virus spreading with warning signals." J Franklin Inst 2024;361:106916.

14. Allen LJS. A primer on stochastic epidemic models: Formulation, numerical simulation, and analysis. Infect Dis Model 2017;2:128-42.

15. Artalejo JR, Economou A, Lopez-Herrero MJ. The stochastic SEIR model before extinction: Computational approaches. Appl Math Comput 2015;265:1026-43.

16. Amador J. The SEIQS stochastic epidemic model with external source of infection. Appl Math Model 2016;40:8352-65.

17. Economou A, Lopez-Herrero MJ. The deterministic SIS epidemic model in a Markovian random environment. J Math Biol 2016;73:91-121.

18. Gamboa M, Lopez-Herrero MJ. Measuring Infection Transmission in a Stochastic SIV Model with Infection Reintroduction and Imperfect Vaccine. Acta Biotheor 2020;68:395-420.

19. Gamboa M, Lopez-Herrero MJ. On the number of periodic inspections during outbreaks of discrete-time stochastic SIS epidemic models. Mathematics 2018;6.

20. Gamboa M, Lopez-Herrero MJ. The effect of setting a warning vaccination level on a stochastic SIVS model with imperfect vaccine. Mathematics 2020;8.

21. Gamboa M, Lopez-Herrero M. Measures to assess a warning vaccination level in a stochastic SIV model with imperfect vaccine. Studies in Applied Mathematics 2022;148:1411-38.

22. Li QL, Wang C, Yang F, Zhang C. Markov modeling and performance analysis of infectious diseases with asymptomatic patients. Mathematical Biosciences and Engineering 2023;20:17822-48.

23. Gómez-Corral A, López-García M, Rodríguez-Bernal MT. On time-discretized versions of the stochastic SIS epidemic model: a comparative analysis. J Math Biol 2021;82.

24. Artalejo JR, Economou A, Lopez-Herrero MJ. Stochastic epidemic models revisited: Analysis of some continuous performance measures. J Biol Dyn 2012;6:189-211.

25. Artalejo JR, Lopez-Herrero MJ. Stochastic epidemic models: New behavioral indicators of the disease spreading. Appl Math Model 2014;38:4371-87.

26. Artalejo JR, Lopez-Herrero MJ. Quasi-stationary and ratio of expectations distributions: A comparative study. J Theor Biol 2010;266:264-74.

27. Artalejo JR, Economou A, Lopez-Herrero MJ. On the number of recovered individuals in the SIS and SIR stochastic epidemic models. Math Biosci 2010;228:45-55.

28. Papageorgiou VE, Tsaklidis G. A stochastic SIRD model with imperfect immunity for the evaluation of epidemics. Appl Math Model 2023;124:768-90.

29. Alyami L, Das S. Extended Skew Kalman Filters for COVID-19 Pandemic State Estimation. In: 2023 Sixth International Conference of Women in Data Science at Prince Sultan University (WiDS PSU). IEEE; 2023. page 162-7.
30. Alyami L, Das S. State Estimation of the Spread of COVID-19 in Saudi Arabia using Extended Kalman Filter. In: 2022 Sensor Signal Processing for Defence Conference, SSPD 2022 - Proceedings. Institute of Electrical and Electronics Engineers Inc.; 2022.

31. Alyami L, Panda DK, Das S. Bayesian Noise Modelling for State Estimation of the Spread of COVID-19 in Saudi Arabia with Extended Kalman Filters. Sensors 2023;23.

32. Sebbagh A, Kechida S. EKF-SIRD model algorithm for predicting the coronavirus (COVID-19) spreading dynamics. Sci Rep 2022;12.

33. Zhu X, Gao B, Zhong Y, Gu C, Choi KS. Extended Kalman filter based on stochastic epidemiological model for COVID-19 modelling. Comput Biol Med 2021;137.

34. Papageorgiou VE, Tsaklidis G. An improved epidemiological-unscented Kalman filter (hybrid SEIHCRDV-UKF) model for the prediction of COVID-19. Application on real-time data. Chaos Solitons Fractals 2023;166.

35. Song J, Xie H, Gao B, Zhong Y, Gu C, Choi KS. Maximum likelihood-based extended Kalman filter for COVID-19 prediction. Chaos Solitons Fractals 2021;146.

36. Ghostine R, Gharamti M, Hassrouny S, Hoteit I. An extended SEIR model with vaccination for forecasting the COVID-19 pandemic in Saudi Arabia using an Ensemble Kalman Filter. Mathematics 2021;9.

37. Masum Bhuiyan M Al, Mahmud S, Islam MR, Tasnim N. Volatility estimation for COVID-19 daily rates using Kalman filtering technique. Results Phys 2021;26.

38. Singh KK, Kumar S, Dixit P, Bajpai MK. Kalman filter based short term prediction model for COVID-19 spread. Applied Intelligence 2021;51:2714-26.

39. Arroyo-Marioli F, Bullano F, Kucinskas S, Rondón-Moreno C. Tracking R of COVID-19: A new real-time estimation using the Kalman filter. PLoS One 2021;16.

40. Özbek L. A study on the estimation of COVID-19 case with the adaptive kalman filter (AKF) in USA, Germany, United Kingdom, Italy, France, Russia, Brazil, India, Türkiye, Spain, Peru, Colombia, South Africa, Argentina, Iran, Pakistan. Sigma Journal of Engineering and Natural Sciences 2022;40:323-43.

41. Papageorgiou VE, Kolias P. A novel epidemiologically informed particle filter for assessing epidemic phenomena. Application to the monkeypox outbreak of 2022. Inverse Probl 2024;

42. Papageorgiou VE, Tsaklidis G. A stochastic particle extended SEIRS model with repeated vaccination: Application to real data of COVID-19 in Italy. Math Methods Appl Sci 2024;47:6504-38.

43. Endo A, van Leeuwen E, Baguelin M. Introduction to particle Markov-chain Monte Carlo for disease dynamics modellers. Epidemics 2019;29.

44. CHEN SM, HSU YS, PEARN WL. Bayesian Filtering: From Kalman Filters to Particle Filters, and Beyond. Statistics: A Journal of Theoretical and Applied Statistics 2003;37:1-24.

45. Brosius I, Van Dijck C, Coppens J, Vandenhove L, Bangwen E, Vanroye F, et al. Presymptomatic viral shedding in high-risk mpox contacts: A prospective cohort study. J Med Virol 2023;95.
46. Wu P, Liu F, Chang Z, Lin Y, Ren M, Zheng C, et al. Assessing Asymptomatic, Presymptomatic, and Symptomatic Transmission Risk of Severe Acute Respiratory Syndrome Coronavirus 2. Clinical Infectious Diseases 2021;73:e1314-20.

47. Virlogeux V, Fang VJ, Wu JT, Ho LM, Peiris JSM, Leung GM, et al. Brief Report: Incubation Period Duration and Severity of Clinical Disease Following Severe Acute Respiratory Syndrome Coronavirus Infection. Epidemiology 2015;26:666-9.

48. World Health Organization (WHO). Monkeypox [Internet]. 2024 [cited 2024 Sep 30];Available from: https://www.who.int/news-room/fact-sheets/detail/mpox

49. Amador J, Artalejo JR. Stochastic modeling of computer virus spreading with warning signals. J Franklin Inst 2013;350:1112-38.

50. Artalejo JR, Lopez-Herrero MJ. On the Exact Measure of Disease Spread in Stochastic Epidemic Models. Bull Math Biol 2013;75:1031-50.

51. Amador J, Artalejo JR. Modeling computer virus with the BSDE approach. Computer Networks 2013;57:302-16.

52. Dempster AP, Laird NM, Rubin DB. Maximum Likelihood from Incomplete Data via the EM Algorithm. Journal of the Royal Statistical Society. Series B (Methodological) 1977;39:1-38.

53. Arulampalam MS, Maskell S, Gordon N, Clapp T. A tutorial on particle filters for online nonlinear/non-Gaussian Bayesian tracking. IEEE Transactions on Signal Processing 2002;50:174-88.

54. Higuchi T. Monte carlo filter using the genetic algorithm operators. J Stat Comput Simul 1997;59:1-23.

55. Kitagawa G. A Self-Organizing State-Space Model. J Am Stat Assoc 1998;93:1203.

56. Elfring J, Torta E, van de Molengraft R. Particle Filters: A Hands-On Tutorial. Sensors 2021;21:438.

57. Johndrow J, Ball P, Gargiulo M, Lum K. Estimating the Number of SARS-CoV-2 Infections and the Impact of Mitigation Policies in the United States. Harv Data Sci Rev 2020;

58. Calvetti D, Hoover A, Rose J, Somersalo E. Bayesian dynamical estimation of the parameters of an SE(A)IR COVID-19 spread model. 2020;Available from: http://arxiv.org/abs/2005.04365

59. Childs ML, Kain MP, Harris MJ, Kirk D, Couper L, Nova N, et al. The impact of long-term non-pharmaceutical interventions on COVID-19 epidemic dynamics and control: the value and limitations of early models. Proceedings of the Royal Society B: Biological Sciences 2021;288:20210811.

60. Lykou R, Tsaklis G. Particle filtering: A priori estimation of observational errors of a state-space model with linear observation equation. Mathematics 2021;9.

61. Edouard Mathieu, Fiona Spooner, Saloni Dattani, Hannah Ritchie, Max Roser. Mpox. Our World in Data2022;
62. Patel A, Bilinska J, Tam JCH, Da Silva Fontoura D, Mason CY, Daunt A, et al. Clinical features and novel presentations of human monkeypox in a central London centre during the 2022 outbreak: Descriptive case series. The BMJ 2022;

63. European Centre for Disease Prevention and Control (ECDC). Factsheet for health professionals on mpox [Internet]. 2024 [cited 2024 Sep 16];Available from: https://www.ecdc.europa.eu/en/all-topics-z/monkeypox/factsheet-health-professionals

64. Bettencourt LMA. An Ensemble Trajectory Method for Real-Time Modeling and Prediction of Unfolding Epidemics: Analysis of the 2005 Marburg Fever Outbreak in Angola. In: Mathematical and Statistical Estimation Approaches in Epidemiology. Dordrecht: Springer Netherlands; 2009. page 143–61.

65. Krieg SJ, Schnur JJ, Miranda ML, Pfrender ME, Chawla N V. Symptomatic, Presymptomatic, and Asymptomatic Transmission of SARS-CoV-2 in a University Student Population, August–November 2020. Public Health Reports 2022;137:1023–30.

66. He X, Lau EHY, Wu P, Deng X, Wang J, Hao X, et al. Temporal dynamics in viral shedding and transmissibility of COVID-19. Nat Med 2020;26:672–5.

67. Hellewell J, Abbott S, Gimma A, Bosse NI, Jarvis CI, Russell TW, et al. Feasibility of controlling COVID-19 outbreaks by isolation of cases and contacts. Lancet Glob Health 2020;8:e488–96.

68. Amador J. The stochastic SIRA model for computer viruses. Appl Math Comput 2014;232:1112–24.

69. Fernández-Peralta R, Gómez-Corral A. A structured Markov chain model to investigate the effects of pre-exposure vaccines in tuberculosis control. J Theor Biol 2021;509.

70. López-García M. Stochastic descriptors in an SIR epidemic model for heterogeneous individuals in small networks. Math Biosci 2016;271:42–61.

71. Papageorgiou VE. Estimating the prevalence of terrorism under control policies. A statistical modelling approach. Appl Math Model 2025;137:115642.

72. Artalejo JR, Lopez-Herrero MJ. Stochastic epidemic models: New behavioral indicators of the disease spreading. Appl Math Model 2014;38:4371–87.

73. Keeling MJ, Ross JV. On methods for studying stochastic disease dynamics. J R Soc Interface 2008;5:171–81.

74. Dass SC, Kwok WM, Gibson GJ, Gill BS, Sundram BM, Singh S. A data driven change-point epidemic model for assessing the impact of large gathering and subsequent movement control order on COVID-19 spread in Malaysia. PLoS One 2021;16:e0252136.

75. Texier G, Farouh M, Pellegrin L, Jackson ML, Meynard JB, Deparis X, et al. Outbreak definition by change point analysis: a tool for public health decision? BMC Med Inform Decis Mak 2016;16:33.

76. S.T. PK, Lahiri B, Alvarado R. Multiple change point estimation of trends in Covid-19 infections and deaths in India as compared with WHO regions. Spat Stat 2022;49:100538.
77. Koutroumpas K, Ballarini P, Votsi I, Cournède PH. Bayesian parameter estimation for the Wnt pathway: an infinite mixture models approach. Bioinformatics 2016;32:i781–9.

78. Votsi I, Cournède PH. A data augmentation scheme embedding a sequential Monte Carlo method for Bayesian parameter inference in state space models. 48èmes Journées de Statistique de la SFdS, May 2016, Montpellier, France. ffhal-01355334f. 48èmes Journées de Statistique de la SFdS 2016;