# On a new method for the stochastic perturbation of the disease transmission
coefficient in SIS Models

Alberto Lanconelli Dipartimento di Scienze Statistiche Paolo Fortunati,
Università di Bologna, Bologna, Italy. e-mail<EMAIL_ADDRESS>Berk Tan Perçin Dipartimento di Fisica e Astronomia Augusto Righi, Università
di Bologna, Bologna, Italy. e-mail<EMAIL_ADDRESS>
###### Abstract

In this study we investigate a novel approach to stochastically perturb the
disease transmission coefficient, which is a key parameter in susceptible-
infected-susceptible (SIS) models. Motivated by the papers [5] and [2], we
perturb the disease transmission coefficient with a Gaussian white noise,
formally modelled as the time derivative of a mean reverting Ornstein-
Uhlenbeck process. We remark that, thanks to a suitable representation of the
solution to the deterministic SIS model, this perturbation is rigorous and
supported by a Wong-Zakai approximation argument that consists in smoothing
the singular Gaussian white noise and then taking limit of the solution from
the approximated model. We prove that the stochastic version of the classic
SIS model obtained this way preserves a crucial feature of the deterministic
equation: the reproduction number dictating the two possible asymptotic
regimes for the infection, i.e. extinction and persistence, remains unchanged.
We then identify the class of perturbing noises for which this property holds
and propose simple sufficient conditions for that. All the theoretical
discoveries are illustrated and discussed with the help of several numerical
simulations.

Key words and phrases: SIS epidemic model, Itô and Stratonovich stochastic
differential equations, Wong-Zakai approximation, extinction, persistence.

AMS 2000 classification: 60H10, 60H30, 92D30.

##  1 Introduction

There are key biological parameters to analyse in epidemiology. The pathogen’s
type, infection target, pathogen’s vulnerability to immune system,
contagiousness and its life time outside of the host are certainly some of
these key parameters. However, they are not enough to fully describe how an
epidemic evolves. A key element, which makes epidemiology a very complex area,
is the interaction inside the population and its unpredictability.

To address this feature several mathematical models have been proposed and
analysed in the literature: see e.g. [6, 5, 15, 7] and the references quoted
there. In the classic susceptible-infected-susceptible (SIS) model there is no
long term immunity for the infection. One example of such case is Gonorrhea,
as described in [6], where individuals who are recovered from the infection
become infected again and again. This means that long term immunity is not
effective in prevention from re-infection.

From a mathematical point of view, the SIS model (see e.g. [3]) is a very
handy tool that describes the average evolution of an infection with no
immunity. It consists of the system of ordinary differential equations

$\displaystyle\begin{cases}\frac{dS(t)}{dt}=\mu N-\beta S(t)I(t)+\gamma
I(t)-\mu S(t),&S(0)=s_{0}\in]0,N[;\\\ \frac{dI(t)}{dt}=\beta
S(t)I(t)-(\mu+\gamma)I(t),&I(0)=i_{0}\in]0,N[,\end{cases}$ (1.1)

where $S(t)$ and $I(t)$ denote the number of _susceptibles_ and _infecteds_ at
time $t$, respectively. Here, $N:=s_{0}+i_{0}$ is the initial size of the
population amongst whom the disease is spreading, $\mu$ denotes the per capita
death rate, $\gamma$ is the rate at which infected individuals become cured
and $\beta$ stands for the disease transmission coefficient. If we sum the
equations in (1.1), we get that

$\displaystyle\frac{d}{dt}(S(t)+I(t))=\mu(N-(S(t)+I(t))),\quad S(0)+I(0)=N,$

which yields

$\displaystyle S(t)+I(t)=S(0)+I(0)=N,\quad\mbox{ for all $t\geq 0$}.$

Therefore, system (1.1) can be reduced to the differential equation

$\displaystyle\frac{dI(t)}{dt}=\beta I(t)(N-I(t))-(\mu+\gamma)I(t),\quad
I(0)=i_{0}\in]0,N[,$ (1.2)

with $S(t):=N-I(t)$. Furthermore, equation (1.2) can be solved explicitly as

$\displaystyle I(t)=\frac{i_{0}e^{[\beta
N-(\mu+\gamma)]t}}{1+\beta\int_{0}^{t}i_{0}e^{[\beta
N-(\mu+\gamma)]s}ds},\quad t\geq 0.$ (1.3)

This explicit representation easily identifies two different asymptotic
regimes for the solution $I(t)$, namely

$\displaystyle\lim_{t\to+\infty}I(t)=\begin{cases}0,&\mbox{ if $R^{D}_{0}\leq
1$};\\\ N\left(1-\frac{1}{R_{0}^{D}}\right),&\mbox{ if
$R^{D}_{0}>1$},\end{cases}$ (1.4)

where

$\displaystyle R^{D}_{0}:=\frac{\beta N}{\mu+\gamma}.$

This ratio is known as _basic reproduction number_ of the infection and
determines whether the disease will become extinct, i.e. $I(t)$ will tend to
zero as $t$ goes to infinity, or will be persistent, i.e. $I(t)$ will tend to
a positive limit as $t$ increases.

As seen from (1.3) and identity $S(t)=N-I(t)$, the infected and susceptible
populations will have _smooth_ flows to each other with rates specified by the
model. This description is very good at yielding overall results for a
population but it fails to capture its heterogeneity, hence realism. In
reality, individuals have different recovering rates or get infected at
different rates, thus creating a distortion in the smooth flow between the two
populations. To describe this randomness, several approaches have been
proposed in the literature.
One example is the approach of Allen [1]: here one starts with a discrete
Markov chain whose transition probabilities reflect the dynamical behaviour of
the deterministic model; then, via a suitable scaling on the one-step
transition probability, one obtains a forward Fokker-Planck equation which is
canonically associated with a stochastic differential equation.
Another common method for introducing stochasticity is the so-called
_parameter perturbation_ approach [5, 14, 4]: it amounts at perturbing one of
the parameters of the model equation with a suitable source of randomness,
usually a Gaussian white noise. One of the most representative papers in this
direction is [5]: here the authors formally perturb equation (1.2), rewritten
in the form

$\displaystyle dI(t)=\beta I(t)(N-I(t))dt-(\mu+\gamma)I(t)dt,\quad
I(0)=i_{0}\in]0,N[,$

through the replacement

$\displaystyle\beta dt\mapsto\beta dt+\sigma dB_{t}$ (1.5)

with $\\{B_{t}\\}_{t\geq 0}$ being a standard one dimensional Brownian motion
and $\sigma$ an additional parameter of the model. Since $\beta$ is a
parameter for disease transmission rate, the term $\beta dt$ can be
interpreted as number of transmissions in time interval $[t,t+dt]$, as stated
in [5]. This way, the authors propose the model

$\displaystyle dI(t)=[\beta I(t)(N-I(t))-(\mu+\gamma)I(t)]dt+\sigma
I(t)(N-I(t))dB(t),$ (1.6)

interpreted as an Itô-type stochastic differential equation, which will
encapsulate the randomness in the disease transmission; moreover, they
identify a _stochastic reproduction number_

$\displaystyle R_{0}^{S}:=R^{D}_{0}-\frac{\sigma^{2}N^{2}}{2(\mu+\gamma)},$

which, in contrast to (1.4), characterizes the following asymptotic
behaviours:

  * •

if $R_{0}^{S}<1$ and $\sigma^{2}<\frac{\beta}{N}$ or if
$\sigma^{2}>\max\\{\frac{\beta}{N},\frac{\beta^{2}}{2(\mu+\gamma)}\\}$, then
the _disease will become extinct_ , i.e.

$\displaystyle\lim_{t\to+\infty}I(t)=0;$

  * •

if $R_{0}^{S}>1$, then _the disease will be persistent_ , i.e.

$\displaystyle\liminf_{t\to+\infty}I(t)\leq\xi\leq\limsup_{t\to+\infty}I(t),$

where
$\xi:=\frac{1}{\sigma^{2}}\left(\sqrt{\beta^{2}-2\sigma^{2}(\mu+\gamma)}-(\beta-\sigma^{2}N)\right)$.

(see also [14]). It is worth mentioning that going from (1.2) to (1.6), as
described in [5], one has to accept some reasonable but heuristic
manipulations of the infinitesimal quantities $dt$ and $dB_{t}$.

The aim of this paper is to propose a different method for perturbing the
disease transmission rate in the SIS model (1.2). Our idea stems from the
following simple observation: if we let $\beta$ in (1.2) to be a function of
time, then the solution formula (1.3) takes the form

$\displaystyle I(t)$
$\displaystyle=\frac{i_{0}e^{N\int_{0}^{t}\beta(s)ds-(\mu+\gamma)t}}{1+i_{0}\int_{0}^{t}\beta(s)e^{N\int_{0}^{s}\beta(r)dr-(\mu+\gamma)s}ds}$
$\displaystyle=\frac{i_{0}e^{N\int_{0}^{t}\beta(s)ds-(\mu+\gamma)t}}{1+\frac{i_{0}}{N}\left(e^{N\int_{0}^{t}\beta(r)dr-(\mu+\gamma)t}-1+\int_{0}^{t}e^{N\int\limits_{0}^{s}\beta(r)dr-(\gamma+\mu)s}(\gamma+\mu)ds\right)},$
(1.7)

where in the second equality we performed an integration by parts in the
denominator. Equation (1) now depends on the function $\beta(t)$ only through
its integral $\int_{0}^{t}\beta(s)ds$. In [2] the authors utilized this
approach to mimic the perturbation proposed in [5]; in this case, the
_singular_ perturbation

$\displaystyle\beta(t)\mapsto\beta(t)+\sigma\frac{dB_{t}}{dt},$

formally employed on the differential equation (1.2) by the authors in [5],
becomes the well defined transformation

$\displaystyle\int_{0}^{t}\beta(s)ds\mapsto\int_{0}^{t}\beta(s)ds+\sigma
B_{t},$

if directly applied on the explicit solution (1). As shown in [2], this
different procedure of parameter perturbation results in an alternative
stochastic SIS model which surprisingly exhibits the same asymptotic regimes
of its deterministic counterpart (1.2). It is important to remark that this
new parameter perturbation approach, which directly acts on the explicit
solution (1), is also justified via Wong-Zakai theorem at the level of
differential equations, thus ruling out the necessity of having a closed form
expression for the solution.
In the current paper we employ the just mentioned approach to the case where
the perturbation is modelled as a mean reverting Ornstein -Uhlenbeck process.
This choice is suggested, but not investigated, both in [1] and [5]. From a
modelling point of view it is motivated by the fact that the variance of a
mean reverting Ornstein -Uhlenbeck process is bounded in time, while the one
of a Brownian motion, utilized in [5], is not. This feature seems to be more
realistic and hence desirable. However, from a mathematical point of view, the
perturbation with a mean reverting Ornstein -Uhlenbeck process makes the
analysis of the model more demanding since in this case equation (1.6) becomes
a stochastic differential equation with random coefficients.
This problem is discussed in [4] with an approach that follows [5]. Here, we
introduce the model working directly on the explicit representation (1) and
cross-validate the proposal from a differential equations’ perspective passing
through the Wong-Zakai theorem. We prove that our model fulfils some basic
biological constraints, i.e. the solution is global and lives in the interval
$]0,N[$ with probability one. Then, we analyse the asymptotic behaviour and
discover that the threshold for the different regimes coincides with the one
for the deterministic SIS model; in other words, the parameters describing the
mean reverting Ornstein -Uhlenbeck process do not play any role in the
limiting behaviour of the solution. We also identify a class of perturbations
for which this invariance is preserved thus offering a complete analysis of
our approach.

The paper is organized as follows: in Section 2 we introduce the model and its
cross-validation via the Wong-Zakai theorem; Section 3 is devoted to the
analysis of our model: support of the solution, extinction, persistence and
discussion of several numerical simulations. Lastly, in Section 4 we address
the problem of finding a general class of perturbations for which the results
from Section 3 remain the same; numerical simulations are also presented for
this enlarged framework.

##  2 Stochastic parameter perturbation with a mean reverting Ornstein-
Uhlenbeck process

Let $\\{Y_{t}\\}_{t\geq 0}$ be a mean reverting Ornstein-Uhlenbeck process
driven by a standard one dimensional Brownian motion $\\{B_{t}\\}_{t\geq 0}$;
this means that $\\{Y_{t}\\}_{t\geq 0}$ is the unique strong solution of the
stochastic differential equation

$\displaystyle dY_{t}=-\alpha Y_{t}dt+\sigma dB_{t},\quad Y_{0}=0$ (2.1)

where the parameters $\alpha$ and $\sigma$ are positive real number. the
process $\\{Y_{t}\\}_{t\geq 0}$ can be explicitly represented as

$\displaystyle Y_{t}=\sigma\int\limits_{0}^{t}e^{-\alpha(t-s)}dB_{s},\quad
t\geq 0,$ (2.2)

entailing that $Y_{t}$ is a Gaussian random variable with mean zero and
variance $\frac{1-e^{-2\alpha t}}{2\alpha}\sigma^{2}$. We also recall the
ergodic property of $\\{Y_{t}\\}_{t\geq 0}$:

$\displaystyle\lim_{t\to+\infty}\frac{1}{t}\int_{0}^{t}Y_{s}ds=0\quad\mbox{almost
surely}.$ (2.3)

We now perturb (1) via the substitution

$\displaystyle\beta(t)\mapsto\beta+\frac{dY_{t}}{dt}$ (2.4)

or more rigorously

$\displaystyle\int_{0}^{t}\beta(s)ds\mapsto\int_{0}^{t}\left(\beta+\frac{dY_{s}}{ds}\right)ds=\beta
t+Y_{t}.$ (2.5)

This gives

$\displaystyle\mathtt{I}_{t}:=\frac{i_{0}e^{\nu
t+NY_{t}}}{1+\frac{i_{0}}{N}\left(e^{\nu t+NY_{t}}-1+\int_{0}^{t}e^{\nu
s+NY_{s}}(\gamma+\mu)ds\right)},\quad t\geq 0,$ (2.6)

where to ease the notation we set $\nu:=N\beta-(\gamma+\mu)$ and
$\mathtt{I}_{t}$ instead of $I(t,Y_{t})$; note that $R_{0}^{D}\leq 1$ is
equivalent to $\nu\leq 0$. The stochastic process (2.6) is the object of our
investigation. Observe that an application of the Itô formula gives

$\displaystyle d\mathtt{I}_{t}=$
$\displaystyle\left[\mathtt{I}_{t}(N-\mathtt{I}_{t})\left(\frac{\nu}{N}-\alpha
Y_{t}+\frac{N-2\mathtt{I}_{t}}{2}\sigma^{2}\right)-\frac{\gamma+\mu}{N}\mathtt{I}_{t}^{2}\right]dt+\sigma\mathtt{I}_{t}(N-\mathtt{I}_{t})dB_{t}.$
(2.7)

This equation can be considered either as a one dimensional stochastic
differential equation with random coefficients (for the presence of
$\\{Y_{t}\\}_{t\geq 0}$) or, if coupled with (2.1), as a two dimensional
system of stochastic differential equations. The local Lipschitz continuity of
the coefficients of such system entails path-wise uniqueness and hence that
the couple $\\{(\mathtt{I}_{t},Y_{t})\\}_{t\geq 0}$, with $\mathtt{I}_{t}$
defined in (2.6) and $Y_{t}$ defined in (2.2), is its unique solution (see
e.g. Theorem 2.5, Chapter 5 in [8]).

###  2.1 Cross-validation of the model via Wong-Zakai theorem

We obtained the stochastic process (2.6) perturbing the explicit solution (1)
with the transformation (2.5). One can however derive the stochastic
differential equation (2.7), which is uniquely solved by (2.6), through a
parameter perturbation procedure acting on the deterministic equation (1.2),
which resembles the approach employed in [5].
Let $\\{B^{\pi}_{t}\\}_{t\in[0,T]}$ be the polygonal approximation of the
Brownian motion $\\{B_{t}\\}_{t\in[0,T]}$, relative to the partition $\pi$.
This means that $\\{B^{\pi}_{t}\\}_{t\in[0,T]}$ is a continuous piecewise
linear random function converging to $\\{B_{t}\\}_{t\in[0,T]}$ almost surely
and uniformly on $[0,T]$, as the mesh of the partition tends to zero. Now,
replace (2.1) with

$\displaystyle\frac{dY^{\pi}_{t}}{dt}=-\alpha
Y^{\pi}_{t}+\sigma\frac{dB^{\pi}_{t}}{dt},\quad Y_{0}^{\pi}=0,$ (2.8)

which gives a smooth approximation of $\\{Y_{t}\\}_{t\geq 0}$. Using
$\\{Y^{\pi}_{t}\\}_{t\geq 0}$ instead of $\\{Y_{t}\\}_{t\geq 0}$ allows for a
rigorous implementation of the transformation (2.4) in (1.2), that means

$\displaystyle\frac{d\mathtt{I}_{t}^{\pi}}{dt}=$
$\displaystyle\beta\mathtt{I}_{t}^{\pi}(N-\mathtt{I}_{t}^{\pi})-(\mu+\gamma)\mathtt{I}_{t}^{\pi}+\mathtt{I}_{t}^{\pi}(N-\mathtt{I}_{t}^{\pi})\frac{dY^{\pi}_{t}}{dt}$
$\displaystyle=$
$\displaystyle\beta\mathtt{I}_{t}^{\pi}(N-\mathtt{I}_{t}^{\pi})-(\mu+\gamma)\mathtt{I}_{t}^{\pi}-\alpha\mathtt{I}_{t}^{\pi}(N-\mathtt{I}_{t}^{\pi})Y^{\pi}_{t}$
$\displaystyle+\sigma\mathtt{I}_{t}^{\pi}(N-\mathtt{I}_{t}^{\pi})\frac{dB^{\pi}_{t}}{dt}.$
(2.9)

According to the Wong-Zakai Theorem [13],[11] the unique solution
$\\{\mathtt{I}_{t}^{\pi}\\}_{t\in[0,T]}$ of the random ordinary differential
equation (2.1) converges, as the mesh of the partition $\pi$ tends to zero, to
the solution of the Stratonovich stochastic differential equation

$\displaystyle d\mathtt{I}_{t}=$
$\displaystyle\left[\beta\mathtt{I}_{t}(N-\mathtt{I}_{t})-(\mu+\gamma)\mathtt{I}_{t}-\alpha\mathtt{I}_{t}(N-\mathtt{I}_{t})Y_{t}\right]dt+\sigma\mathtt{I}_{t}(N-\mathtt{I}_{t})\circ
dB_{t},$

which in turn is equivalent to the Itô SDE

$\displaystyle d\mathtt{I}_{t}=$
$\displaystyle\left[\beta\mathtt{I}_{t}(N-\mathtt{I}_{t})-(\mu+\gamma)\mathtt{I}_{t}-\alpha\mathtt{I}_{t}(N-\mathtt{I}_{t})Y_{t}+\frac{\sigma^{2}}{2}\mathtt{I}(t)(N-\mathtt{I}(t))(N-2\mathtt{I}(t))\right]dt$
$\displaystyle+\sigma\mathtt{I}_{t}(N-\mathtt{I}_{t})dB_{t}.$

The stochastic differential equation above coincides with (2.7) thus
validating our parameter perturbation approach also from a model equation pont
of view.

##  3 Analysis of the stochastically perturbed SIS Model

In this section we analyse the stochastic process (2.6) which we recall to be
the unique strong solution of the SDE (2.7). We will in particular show that
such process lives in the interval $]0,N[$, for all $t\geq 0$, almost surely
and we will provide sufficient conditions for extinction and persistence.

###  3.1 Support of the solution

We start with the following.

######  Proposition 3.1.

For the stochastic process $\\{\mathtt{I}_{t}\\}_{t\geq 0}$ defined in (2.6)
we have

$\displaystyle\mathbb{P}(\mathtt{I}_{t}\in]0,N[)=1,\quad\mbox{for all $t\geq
0$}.$

###### Proof.

First of all, we observe that $\mathtt{I}_{t}$ can be rewritten as

$\displaystyle\mathtt{I}_{t}=\frac{Ni_{0}e^{\nu t+NY_{t}}}{N-i_{0}+i_{0}e^{\nu
t+NY_{t}}+i_{0}\int_{0}^{t}e^{\nu s+NY_{s}}(\gamma+\mu)ds}.$

Since by assumption $0<i_{0}<N$, we see that $\mathtt{I}_{t}$ is a ratio of
almost sure positive quantities; this yields $\mathtt{I}_{t}>0$ for all $t\geq
0$ almost surely. On the other hand, the last identity also gives

$\displaystyle\mathtt{I}_{t}$ $\displaystyle<\frac{Ni_{0}e^{\nu
t+NY_{t}}}{i_{0}e^{\nu t+NY_{t}}+i_{0}\int_{0}^{t}e^{\nu
s+NY_{s}}(\gamma+\mu)ds}$
$\displaystyle=\frac{N}{1+(\gamma+\mu)\int\limits_{0}^{t}e^{\nu(s-t)+N(Y_{s}-Y_{t})}ds}$
$\displaystyle<N.$

The proof is complete. ∎

###  3.2 Extinction of the infection

We now provide a sufficient condition for extinction; remarkably, the
parameters describing the stochastic perturbation, i.e. $\alpha$ and $\sigma$,
do not play role in that.

######  Theorem 3.2.

If $R^{D}_{0}=\frac{\beta N}{\gamma+\mu}\leq 1$, or equivalently $\nu=\beta
N-(\gamma+\nu)\leq 0$, then

$\displaystyle\lim_{t\to\infty}\mathtt{I}_{t}=0\quad\mbox{ almost surely}.$

###### Proof.

We take $G(x):=\ln\left(\frac{x}{N-x}\right)$ for $x\in]0,N[$ and observe that
$G$ is a strictly increasing function that maps the interval $]0,N[$ into
$]-\infty,+\infty[$. An application of the Itô formula gives

$\displaystyle dG(\mathtt{I}_{t})=\left[\nu-N\alpha
Y_{t}-(\gamma+\mu)\frac{\mathtt{I}_{t}}{N-\mathtt{I}_{t}}\right]dt+\sigma
NdB_{t},$

which corresponds to the integral equation

$\displaystyle G(\mathtt{I}_{t})=\ln\left(\frac{i_{0}}{N-i_{0}}\right)-\alpha
N\int\limits_{0}^{t}Y_{s}ds+\int\limits_{0}^{t}f(\mathtt{I}_{s})ds+\sigma
NB_{t},$ (3.1)

where

$\displaystyle f(x):=\nu-(\gamma+\mu)\frac{x}{N-x},\quad\mbox{ for
$x\in]0,N[$}.$ (3.2)

It is useful to note that $f$ is monotone decreasing on the interval $]0,N[$
and that $f(x)<\nu$, for $x\in]0,N[$. Therefore, from equation (3.1) we get

$\displaystyle\limsup_{t\to\infty}\frac{1}{t}G(\mathtt{I}_{t})\leq$
$\displaystyle\limsup_{t\to\infty}\frac{1}{t}\ln\left(\frac{i_{0}}{N-i_{0}}\right)-\lim_{t\to\infty}\frac{\alpha
N}{t}\int\limits_{0}^{t}Y_{s}ds$
$\displaystyle+\limsup_{t\to\infty}\frac{1}{t}\int\limits_{0}^{t}f(\mathtt{I}_{s})ds+\sigma
N\limsup_{t\to\infty}\frac{B_{t}}{t}.$ (3.3)

By recalling (2.3), it is easy to notice that the first two terms in right
hand side above are equal to zero. Moreover, by the strong law of large
numbers for martingales (see for instance [10]) we also have

$\displaystyle\lim_{t\to\infty}\frac{B_{t}}{t}=0,\quad\mbox{almost surely}.$

Therefore, inequality (3.2) now reads

$\displaystyle\limsup_{t\to\infty}\frac{1}{t}G(\mathtt{I}_{t})\leq\limsup_{t\to\infty}\frac{1}{t}\int\limits_{0}^{t}f(\mathtt{I}_{s})ds<\nu,$

that means

$\displaystyle\limsup_{t\to\infty}\frac{1}{t}\ln\left(\frac{\mathtt{I}_{t}}{N-\mathtt{I}_{t}}\right)<0,\quad\mbox{almost
surely}.$

Since the last statement implies our thesis, the proof is complete. ∎

In the stochastic SIS model obtained by this parameter perturbation method
from the deterministic one, the limiting behaviour of both models are the same
for $R^{D}_{0}\leq 1$.

###  3.3 Persistence of the Infection

We now turn to the problem of finding sufficient condition for persistence of
the disease. Again, the parameters describing the stochastic perturbation,
i.e. $\alpha$ and $\sigma$, do not influence the threshold.

######  Theorem 3.3.

If $R^{D}_{0}>1$, or equivalently $\nu>0$, then we have with probability one

$\displaystyle\limsup_{t\to\infty}\mathtt{I}_{t}\geq x^{*}\quad\mbox{ and
}\quad\liminf_{t\to\infty}\mathtt{I}_{t}\leq x^{*},$ (3.4)

where $x^{*}=N\left(1-\frac{1}{R^{D}_{0}}\right)$.

###### Proof.

The assumption of the theorem implies that the function in (3.2) has a unique
root in $]0,N[$ given by $x^{*}=N(1-1/R^{D}_{0})$, as shown in figure 1.

Figure 1: The plot of $f(x)$ from (3.2) with $\nu=20$, $\gamma+\mu=0.2$ and
$N=10$

We follow the proof of Theorem 5.1 in [5]. Assume the first inequality in
(3.4) to be false. Then, there exists $\varepsilon>0$ such that

$\displaystyle\mathbb{P}(\Omega_{1})>\varepsilon\quad\mbox{ where
}\quad\Omega_{1}:=\left\\{\limsup_{t\to\infty}\mathtt{I}_{t}\leq
x^{*}-\varepsilon\right\\}.$ (3.5)

Therefore, for all $\omega\in\Omega_{1}$ there exists $T(\omega)\geq 0$ such
that

$\displaystyle\mathtt{I}_{t}\leq x^{*}-\varepsilon,\quad\mbox{ for all $t\geq
T(\omega)$},$

and the monotonicity of $f$ yields

$\displaystyle f(\mathtt{I}_{t})\geq f(x^{*}-\varepsilon)>0,\quad\mbox{ for
all $t\geq T(\omega)$}.$ (3.6)

Therefore, using identities (3.1), (2.3) and the strong law of large numbers
for martingales, we can write for all $\omega\in\Omega_{1}$ that

$\displaystyle\liminf_{t\to\infty}\frac{1}{t}G(\mathtt{I}_{t})\geq$
$\displaystyle\liminf_{t\to\infty}\frac{1}{t}\ln\left(\frac{i_{0}}{N-i_{0}}\right)-\lim_{t\to\infty}\frac{\alpha
N}{t}\int\limits_{0}^{t}Y_{s}ds$
$\displaystyle+\liminf_{t\to\infty}\frac{1}{t}\int\limits_{0}^{t}f(\mathtt{I}_{s})ds+N\sigma\lim_{t\to\infty}\frac{B_{t}}{t}$
$\displaystyle=$
$\displaystyle\liminf_{t\to\infty}\frac{1}{t}\int\limits_{0}^{t}f(\mathtt{I}_{s})ds$
$\displaystyle\geq$
$\displaystyle\liminf_{t\to\infty}\frac{1}{t}\int\limits_{0}^{T}f(\mathtt{I}_{s})ds+f(x^{*}-\varepsilon)\liminf_{t\to\infty}\frac{t-T}{t}$

and hence

$\displaystyle\liminf_{t\to\infty}\frac{1}{t}\ln\left(\frac{\mathtt{I}_{t}}{N-\mathtt{I}_{t}}\right)\geq
f(x^{*}-\varepsilon)>0.$

This gives $\lim\limits_{t\to\infty}\mathtt{I}_{t}=N$ which contradicts (3.5).

To prove the second inequality in (3.4) we proceed as before and assume that
there exists $\varepsilon>0$ such that

$\displaystyle\mathbb{P}(\Omega_{2})>\varepsilon\quad\mbox{ where
}\quad\Omega_{2}=\left\\{\liminf_{t\to\infty}\mathtt{I}_{t}\geq
x^{*}+\varepsilon\right\\}.$ (3.7)

Therefore, for all $\omega\in\Omega_{2}$, there exists $T(\omega)\geq 0$ such
that

$\displaystyle\mathtt{I}_{t}\geq x^{*}-\varepsilon,\quad\mbox{ for all $t\geq
T(\omega)$}.$

The monotonicity of $f$ gives

$\displaystyle f(\mathtt{I}_{t})\leq f(x^{*}-\varepsilon)<0,\quad\mbox{ for
all $t\geq T(\omega)$},$ (3.8)

and

$\displaystyle\limsup_{t\to\infty}\frac{1}{t}\ln\left(\frac{\mathtt{I}_{t}}{N-\mathtt{I}_{t}}\right)\leq\limsup_{t\to\infty}\frac{1}{t}\int\limits_{0}^{t}f(\mathtt{I}_{s})ds\leq
f(x^{*}+\varepsilon)<0.$

This implies $\lim\limits_{t\to\infty}\mathtt{I}_{t}=0$, contradicting (3.7).
∎

###  3.4 Trajectory simulations

In this section we present various simulations; we consider two different
values of $\sigma$ to emphasize that, according to our theoretical results,
the limiting behaviour of the solution doesn’t depend on them.
For the first example, we set $N=200$, $i_{0}=100$, $\beta=0.06$,
$\gamma+\mu=14$, $\alpha=0.4$, so that $R^{D}_{0}=0.8$ and $\nu=-2$. According
to Theorem 3.2 the infection should extinct almost surely. See Figure 2.

(a) $\sigma=0.005$

(b) $\sigma=0.05$

Figure 2: The plot with parameters $N=200$, $i_{0}=100$, $\alpha=0.4$,
$R^{D}_{0}=0.857$ and hence $\nu=-2$. The label for y-axis $I(t,Y_{t})$ stands
for $\mathtt{I}_{t}$.

(a) $\sigma=0.005$

(b) $\sigma=0.05$

Figure 3: The plot with parameters $N=200$, $i_{0}=100$, $\alpha=0.4$,
$R^{D}_{0}=1$ and hence $\nu=0$. The label for y-axis $I(t,Y_{t})$ stands for
$\mathtt{I}_{t}$.

(a) $\sigma=0.005$

(b) $\sigma=0.05$

Figure 4: The plot with parameters $N=200$, $i_{0}=100$, $\alpha=0.4$ and
$R^{D}_{0}=1.2$ so $\nu=2$. The deterministic limit is 33. The label for
y-axis $I(t,Y_{t})$ stands for $\mathtt{I}_{t}$.

Next, the same simulation is performed with a different set of parameters. We
take $N=200$, $i_{0}=100$, $\beta=0.06$, $\gamma+\mu=12$, $\alpha=0.4$ which
gives $R^{D}_{0}=1$ and $\nu=0$. According to Theorem 3.2 the infection should
extinct almost surely. See Figure 3.

Lastly, we consider $N=200$, $i_{0}=100$, $\beta=0.06$, $\gamma+\mu=10$,
$\alpha=0.4$ so $R^{D}_{0}=1.2$ and $\nu=2$. According to Theorem 3.2 the
infection should persist a.s., or more precisely it should oscillate above and
below the deterministic limit. See Figure 4 where fluctuations above and below
the deterministic limit $x^{*}$ is visible as well.

##  4 Stochastic parameter perturbation with a general process

In this section we try to understand to which extent the results of the
previous sections are determined by the choice of the particular perturbation
$\\{Y_{t}\\}_{t\geq 0}$ in (2.1). To this aim we consider $\\{Z_{t}\\}_{t\geq
0}$, solution of the stochastic differential equation

$\displaystyle dZ_{t}=b(t,Z_{t})dt+\sigma dB_{t},\quad Z_{0}=0,$ (4.1)

where $\sigma$ is positive real number and the function
$b:[0,T]\times\mathbb{R}$ is assumed to be globally Lipschitz continuous in
$z$, uniformly in $t$; we also introduce the corresponding perturbation of the
parameter $\beta$, namely

$\displaystyle\int_{0}^{t}\beta(s)ds\mapsto\int_{0}^{t}\left(\beta+\frac{dZ_{s}}{ds}\right)ds=\beta
t+Z_{t}.$ (4.2)

If we employ such transformation in (1), which is the solution of the
deterministic SIS model with a time dependent transmission coefficient
$\beta$, we get the stochastic process

$\displaystyle\mathcal{I}_{t}:=\frac{i_{0}e^{\nu
t+NZ_{t}}}{1+\frac{i_{0}}{N}\left(e^{\nu t+NZ_{t}}-1+\int_{0}^{t}e^{\nu
s+NZ_{s}}(\gamma+\mu)ds\right)},\quad t\geq 0.$ (4.3)

where to ease the notation we set $\mathcal{I}_{t}$ instead of $I(t,Z_{t})$.
Moreover, an application of the Itô formula yields

$\displaystyle d\mathcal{I}_{t}=$
$\displaystyle\left[\mathcal{I}_{t}(N-\mathcal{I}_{t})\left(\frac{\nu}{N}+b(t,Z_{t})+\frac{\sigma^{2}N}{2}-\sigma^{2}\mathcal{I}_{t}\right)-\frac{\gamma+\mu}{N}\mathcal{I}_{t}^{2}\right]dt$
$\displaystyle+\sigma\mathcal{I}_{t}(N-\mathcal{I}_{t})dB_{t}$ (4.4)

######  Remark 4.1.

The restriction to constant diffusion coefficients imposed in (4.1) is due the
necessity of cross validating the model (4.3) also from a differential
equations’ point of view. In fact, if we smooth the process
$\\{Z_{t}\\}_{t\geq 0}$ as

$\displaystyle\frac{dZ^{\pi}_{t}}{dt}=b(t,Z^{\pi}_{t})+\sigma\frac{dB^{\pi}_{t}}{dt},\quad
Z^{\pi}_{0}=0,$

and perturb correspondingly the parameter $\beta$ in equation (1.2), we obtain

$\displaystyle\frac{d\mathcal{I}_{t}^{\pi}}{dt}=$
$\displaystyle\beta\mathcal{I}_{t}^{\pi}(N-\mathcal{I}_{t}^{\pi})-(\mu+\gamma)\mathcal{I}_{t}^{\pi}+b(t,Z^{\pi}_{t})\mathcal{I}_{t}^{\pi}(N-\mathcal{I}_{t}^{\pi})$
$\displaystyle+\sigma\mathcal{I}_{t}^{\pi}(N-\mathcal{I}_{t}^{\pi})\frac{dB^{\pi}_{t}}{dt}.$

According to the Wong-Zakai Theorem the unique solution
$\\{\mathcal{I}_{t}^{\pi}\\}_{t\in[0,T]}$ of the random ordinary differential
equation above converges, as the mesh of the partition $\pi$ tends to zero, to
the solution of the Stratonovich stochastic differential equation

$\displaystyle d\mathcal{I}_{t}=$
$\displaystyle\left[\beta\mathcal{I}_{t}(N-\mathcal{I}_{t})-(\mu+\gamma)\mathcal{I}_{t}+b(t,Z_{t})\mathcal{I}_{t}(N-\mathcal{I}_{t})\right]dt+\sigma\mathcal{I}_{t}(N-\mathcal{I}_{t})\circ
dB_{t},$

which in turn is equivalent to the Itô SDE

$\displaystyle d\mathcal{I}_{t}=$
$\displaystyle\left[\beta\mathcal{I}_{t}(N-\mathcal{I}_{t})-(\mu+\gamma)\mathcal{I}_{t}+b(t,Z_{t})\mathcal{I}_{t}(N-\mathcal{I}_{t})+\frac{\sigma^{2}}{2}\mathcal{I}(t)(N-\mathcal{I}(t))(N-2\mathcal{I}(t))\right]dt$
$\displaystyle+\sigma\mathcal{I}_{t}(N-\mathcal{I}_{t})dB_{t}.$

This SDE coincides with (4) thus validating that model. If we allow $\sigma$
in (4.1) to depend also on $Z$, then this match wouldn’t take place for the
presence of an additional drift term in the equation for $\\{Z_{t}\\}_{t\geq
0}$.

We now start to analyse the properties of $\\{\mathcal{I}_{t}\\}_{t\geq 0}$ by
stating the analogue of Proposition 3.1.

######  Proposition 4.2.

For the stochastic process $\\{\mathcal{I}_{t}\\}_{t\geq 0}$ defined in (4.3)
we have

$\displaystyle\mathbb{P}(\mathcal{I}_{t}\in]0,N[)=1,\quad\mbox{for all $t\geq
0$}.$

###### Proof.

Looking through the proof of Proposition 3.1 one easily see that the same
conclusion holds for $\\{\mathcal{I}_{t}\\}_{t\geq 0}$. ∎

The next theorem provides a sufficient condition on the stochastic process
$\\{Z_{t}\\}_{t\geq 0}$ which guarantees extinction for
$\\{\mathcal{I}_{t}\\}_{t\geq 0}$.

######  Theorem 4.3.

Assume that $R^{D}_{0}=\frac{\beta N}{\gamma+\mu}\leq 1$, or equivalently
$\nu=\beta N-(\gamma+\nu)\leq 0$. If

$\displaystyle\limsup_{t\to\infty}\frac{Z_{t}}{t}\leq 0,\quad\mbox{ almost
surely},$ (4.5)

then

$\displaystyle\lim_{t\to\infty}\mathcal{I}_{t}=0\quad\mbox{ almost surely}.$

###### Proof.

We take $G(y):=\ln\left(\frac{y}{N-y}\right)$, for $y\in]0,N[$; an application
of the Itô formula gives

$\displaystyle
dG(\mathcal{I}_{t})=\left[\nu-(\gamma+\mu)\frac{\mathcal{I}_{t}}{N-\mathcal{I}_{t}}\right]dt+NdZ_{t}.$

this yields

$\displaystyle
G(\mathcal{I}_{t})=\ln\left(\frac{i_{0}}{N-i_{0}}\right)+\int\limits_{0}^{t}f(\mathcal{I}_{s})ds+NZ_{t},$

where $f$ is defined as (3.2) again. By utilizing the monotonicity of $f$, we
get

$\displaystyle\limsup_{t\to\infty}\frac{1}{t}G(\mathcal{I}_{t})<\nu+N\limsup_{t\to\infty}\frac{Z_{t}}{t}.$

Now, if

$\displaystyle\limsup_{t\to\infty}\frac{Z_{t}}{t}\leq 0,\quad\mbox{ almost
surely},$

then $\nu\leq 0$ implies immediately that

$\displaystyle\limsup_{t\to\infty}\frac{1}{t}\ln\left(\frac{\mathcal{I}_{t}}{N-\mathcal{I}_{t}}\right)<0.$

The last inequality entails the statement of our theorem and competes the
proof. ∎

We now state the analogue of Theorem 3.3.

######  Theorem 4.4.

Assume that $R^{D}_{0}>1$, or equivalently $\nu>0$. If

$\displaystyle
0\leq\liminf_{t\to\infty}\frac{Z_{t}}{t}\leq\limsup_{t\to\infty}\frac{Z_{t}}{t}<+\infty,$
(4.6)

then the infection is persistent. More precisely, for all $x\in]0,N[$ we have

$\displaystyle\liminf_{t\to\infty}\mathcal{I}_{t}\leq
x\leq\limsup_{t\to\infty}\mathcal{I}_{t},\quad\mbox{almost surely}.$ (4.7)

###### Proof.

We now proceed by contradiction as in the proof of Theorem 3.3. Similarly, we
take $G(y):=\ln\left(\frac{y}{N-y}\right)$, for $y\in]0,N[$; an application of
the Itô formula gives

$\displaystyle
dG(\mathcal{I}_{t})=\left[\nu-(\gamma+\mu)\frac{\mathcal{I}_{t}}{N-\mathcal{I}_{t}}\right]dt+NdZ_{t}.$

We now fix $0<x<N$ and first prove that

$\displaystyle\limsup_{t\to\infty}\mathcal{I}_{t}\geq x,\quad\mbox{ almost
surely}.$

Assume to the contrary that the event
$\left\\{\limsup_{t\to\infty}\mathcal{I}_{t}<x\right\\}$ has positive
probability. This means that there exists $T(\omega)\geq 0$ such that
$\mathcal{I}_{t}<x$, for all $t\geq T(\omega)$; this, together with the
monotonicity of (3.2), gives

$\displaystyle\liminf_{t\to\infty}\frac{1}{t}\ln\left(\frac{\mathcal{I}_{t}}{N-\mathcal{I}_{t}}\right)\geq$
$\displaystyle\liminf_{t\to\infty}\frac{1}{t}\int_{0}^{t}f(\mathcal{I}_{s})ds+N\liminf_{t\to\infty}\frac{Z_{t}}{t}$
$\displaystyle>$ $\displaystyle f(x)+N\liminf_{t\to\infty}\frac{Z_{t}}{t}.$

Now, the core of Theorem 3.3 is the left hand side above being positive and
this is guaranteed if

$\displaystyle N\liminf_{t\to\infty}\frac{Z_{t}}{t}\geq-f(x).$ (4.8)

In this case
$\liminf_{t\to\infty}\frac{1}{t}\ln\left(\frac{\mathcal{I}_{t}}{N-\mathcal{I}_{t}}\right)>0$
and hence $\lim_{t\to\infty}\mathcal{I}_{t}=N$, contradicting our initial
assumption. Since the range of $f$ is $]-\infty,\nu]$, condition (4.8) is
implied by the first inequality in (4.6) and positivity of $\nu$.
We now prove that

$\displaystyle\liminf_{t\to\infty}\mathcal{I}_{t}\leq x,\quad\mbox{almost
surely}.$

Assume to the contrary that the event
$\left\\{\liminf_{t\to\infty}\mathcal{I}_{t}>x\right\\}$ has positive
probability; this implies the existence of $T(\omega)>0$ such that
$\mathcal{I}_{t}>x$, for all $t\geq T(\omega)$. Then,

$\displaystyle\limsup_{t\to\infty}\frac{1}{t}\ln\left(\frac{\mathcal{I}_{t}}{N-\mathcal{I}_{t}}\right)$
$\displaystyle\leq\limsup_{t\to\infty}\int_{0}^{t}f(\mathcal{I}_{s})ds+N\limsup_{t\to\infty}\frac{Z_{t}}{t}$
$\displaystyle<f(x)+N\limsup_{t\to\infty}\frac{Z_{t}}{t}.$

Now, the core of Theorem 3.3 is the left hand side above being negative and
this is guaranteed if

$\displaystyle N\limsup_{t\to\infty}\frac{Z_{t}}{t}\leq-f(x)$ (4.9)

In this case
$\limsup_{t\to\infty}\frac{1}{t}\ln\left(\frac{\mathcal{I}_{t}}{N-\mathcal{I}_{t}}\right)<0$
and hence $\lim_{t\to\infty}\mathcal{I}_{t}=0$ which contradicts our initial
assumption. Since the range of $f$ is $]-\infty,\nu]$, condition (4.9) is
implied by the last inequality in (4.6). ∎

######  Remark 4.5.

If the process $\\{Z_{t}\\}_{t\geq 0}$ satisfies both (4.5) and (4.6), then
for $\nu\leq 0$ one has extinction and for $\nu>0$ one has persistence for the
associated model. It is useful to note this is the case only when

$\displaystyle\lim_{t\to\infty}\frac{Z_{t}}{t}=0.$ (4.10)

###  4.1 Trajectory simulations

(a) $R^{D}_{0}=0.857$

(b) $R^{D}_{0}=1.000$

(c) $R^{D}_{0}=1.200$

(d) $R^{D}_{0}=1.333$

Figure 5: The plot with parameters $N=200$, $i_{0}=100$, $\alpha=-0.011$ and
$\sigma=0.005$. The label for y-axis $I(t,Z_{t})$ stands for
$\mathcal{I}_{t}$.

In this section a simple example will be given to support the theoretical
results of Theorems 4.3 and 4.4. We choose

$\displaystyle Z_{t}:=\alpha t+\sigma B_{t}$ (4.11)

and note that

$\displaystyle\limsup\limits_{t\to\infty}\frac{Z_{t}}{t}=\liminf\limits_{t\to\infty}\frac{Z_{t}}{t}=\alpha,\quad\mbox{
almost surely}.$

If $\alpha<0$, then the assumption of Theorem 4.3 is fulfilled while those of
Theorem 4.4 are not. This means that for $R^{D}_{0}\leq 1$, or equivalently
$\nu\leq 0$, the extinction is guaranteed; however, for $R^{D}_{0}>1$, or
equivalently $\nu>0$, the persistence of infection is not guaranteed.

Figure 5 supports this claim. As one can see the panels 5(a) and 5(b) shows
the extinction of the infection since Theorem 4.3 is satisfied. However,
panels 5(c) and 5(d) shows both the examples of extinction and persistence
even though the $\nu>0\>(R^{D}_{0}>1)$. Thus one looses the properties of the
deterministic model.

A complementary analysis can be made for $\alpha>0$. In this case the
assumption of Theorem 4.4 will be satisfied while those of Theorem 4.3 are
not. This indicates that the infection will be persistent as long as
$R^{D}_{0}>1$; however, the extinction of infection is not guaranteed for
$R^{D}_{0}\leq 1$. The results of Figure 6 support this claim.

(a) $R^{D}_{0}=0.800$

(b) $R^{D}_{0}=0.857$

(c) $R^{D}_{0}=1.000$

(d) $R^{D}_{0}=1.200$

Figure 6: The plot with parameters $N=200$, $i_{0}=100$, $\alpha=-0.011$ and
$\sigma=0.005$. The label for y-axis $I(t,Z_{t})$ stands for
$\mathcal{I}_{t}$.

In the panel 6(d), the infection is persistent as expected. However, the
extinction is not seen in panels 6(c) and 6(b) since the sufficient condition
in Theorem 4.3 is not satisfied. Because the condition was sufficient and not
necessary, although it is not satisfied, the extinction can still happen for
$R^{D}_{0}\leq 1$, as shown in panel 6(a).

One last important comment is that, when $\alpha=0$, then one has the same
kind of perturbation of the $\beta$ parameter utilized in [5], so now
$Z_{t}=\sigma B_{t}$. If the $\beta$ parameter is perturbed in the same way as
this study and not in [5], since $\lim\limits_{t\to\infty}Z_{t}/t=0$, by the
remark 4.5 it can be said that, for $R^{D}_{0}\leq 1\>(\nu\leq 0)$ cases there
will be extinction and for cases $R^{D}_{0}>1\>(\nu>0)$ there will be
persistence of infection. We would like to stress that conditions are not the
same as published in [5] because the perturbation method is different in this
study.

##  5 Discussion

In this study we propose a new perturbation method for the disease
transmission coefficient in SIS model. Our approach consists in acting
directly on the explicit solution of the deterministic problem, thus avoiding
delicate manipulations of differential quantities such as $dB_{t}$.
Nevertheless, our model is cross validated at the level of differential
equations once we smooth the perturbation and use the Wong-Zakai theorem. We
first use this method with a perturbation of mean reverting type and
properties of the corresponding model are analysed. Then, generalization to
different sources of randomness are investigated.
When the deterministic SIS model 1.1 is perturbed with a mean reverting
Ornstein-Uhlenbeck process $\\{Y_{t}\\}_{t\geq 0}$, the solution
$\\{\mathtt{I}_{t}\\}_{t\geq 0}$ is shown to preserve the deterministic
model’s regimes for extinction and persistence. Namely:

  * •

if $R^{D}_{0}\leq 1$, then we have extinction of infection;

  * •

if $R^{D}_{0}>1$, then we have persistence of infection.

Then, we identified some simple sufficient conditions on the class of possible
perturbations which entail the same key feature for the corresponding models.
This study emphasizes that the methodology of perturbation of deterministic
models is crucial in generating different stochastic versions. For further
studies one can apply the same rationale to other parameters in the
deterministic SIS model, i.e. $\gamma+\mu$, or to other models, such as Lotka-
Volterra type of equations ([9],[12]).

## References

  * [1] Allen, E., Modelling with Itô Stochastic Differential Equations, Springer-Verlag, 2007.
  * [2] E. Bernardi and A.Lanconelli, A note about the invariance of the basic reproduction number for stochastically perturbed SIS models, _arXiv:2005.04973v2_.
  * [3] F. Brauer, L.J.S. Allen, P. Van den Driessche and J. Wu, Mathematical Epidemiology, _Lecture Notes in Mathematics_ , No. 1945, Mathematical Biosciences Subseries, 2008\.
  * [4] Cai, Y., Jiao, J., Gui, Z., Liu, Y., & Wang, W. (2018). Environmental variability in a stochastic epidemic model. Applied Mathematics and Computation, 329, 210-226. doi:10.1016/j.amc.2018.02.009
  * [5] A. Gray, D. Greenhalgh, L. Hu, X. Mao, J. Pan, A stochastic differential equation SIS epidemic model, _SIAM J. Appl. Math._ 71 (3) (2011) 876–902.
  * [6] Hethcote, H.W. and Yorke, J.A., Gonorrhea Transmission Dynamics and Control, Lecture Notes in Biomathematics 56, Springer-Verlag, 1994.
  * [7] Jin, Y., Wang, W., & Xiao, S. (2007). An sirs model with a nonlinear incidence rate. Chaos, SolItôns & Fractals, 34(5), 1482-1497. doi:10.1016/j.chaos.2006.04.022
  * [8] I. Karatzas and S. E. Shreve, _Brownian motion and stochastic calculus_ , Springer-Verlag, New York, 1991.
  * [9] A. J. Lotka, Contribution to quantitative parasitology, _J. Wash. Acad. Sci._ 13 (1923) 152-158.
  * [10] Mao, X. (2011). Stochastic differential equations and applications. In Stochastic differential equations and applications (Second ed., pp. 12-13). Oxford: Woodhead.
  * [11] D. W. Stroock and S. R. S. Varadhan, On the support of diffusion processes with applications to the strong maximum principle, _Proceedings $6$-th Berkeley Symposium Math. Statist. Probab._ 3 (1972) University of California Press, Berkeley, 333-359.
  * [12] V. Volterra, Variazioni e fluttuazioni del numero d’individui in specie animali conviventi. _Mem. Acad. Lincei_ 2 (1926) 31-113.
  * [13] E. Wong and M. Zakai: On the relation between ordinary and stochastic differential equations, Intern. J. Engr. Sci. 3 (1965) 213-229.
  * [14] Xu, C. (2017). Global threshold dynamics of a stochastic differential equation sis model. Journal of Mathematical Analysis and Applications, 447(2), 736-757. doi:10.1016/j.jmaa.2016.10.041
  * [15] Zhou, Y., & Liu, H. (2003). Stability of periodic solutions for an Sis model with Pulse vaccination. Mathematical and Computer Modelling, 38(3-4), 299-308. doi:10.1016/s0895-7177(03)90088-4


# Machine Learning for
Detecting Malware in PE Files

C.C.Collin Connors  D.S.Dillip Sarkar

###### Abstract

The increasing number of sophisticated malware poses a major cybersecurity
threat. Portable executable (PE) files are a common vector for such malware.
In this work we review and evaluate machine learning-based PE malware
detection techniques. Using a large benchmark dataset, we evaluate features of
PE files using the most common machine learning techniques to detect malware.

Malware Detection,

Machine Learning,

Portable Executable Files,

EMBER Dataset,

Neural Network,

###### keywords:

Research

##  1 Introduction

_Mal_ icious soft _ware_ (Malware) infects digital systems (computers,
smartphones, networks, cloud servers, etc.), harming the systems or their
users. The degree of harm to a system varies from reduced performance to
unavailability for use or stealing sensitive information; the degree of harm
to the users varies from damaging their documents to stealing users’
information while using an infected system. Based on behavior and damage to an
infected system, malware is classified as adware, backdoor, bot, downloader,
launcher, ransomware, rootkit, spyware, Trojan, virus, and worm, etc. (see [2,
6], and references therein).

Early forms of malware, commonly known as traditional malware, had
characteristic signature strings, and they are identified by detecting those
signatures. However, soon a race started between malware developers and anti-
malware software developers. Newer malware utilized sophisticated and advanced
techniques; thus, they have no simple signature for detection; their
replicates have different appearances yet similar harmful capabilities.

Several methods were developed for analysis and mining of advanced malware to
discover their underlying characteristic features [8].

In addition to malware signatures, these mining tools analyzed DLL function
calls, hexadecimal sequences, assembly instructions, PE file header, or some
combination of these features [6, 2, 14, 8]. The features obtained after the
initial pass may be further processed before feeding them in a tool. For
example, entropy and hashing trick are used in [1]. Note that the so called
Hashing Trick is not a trick but rather a well studied dimension reduction
technique with a solid theoretical foundation [4].

Advanced malware is adaptive and persistent; they evade detection after they
infect systems. Thus, before any software file enters a digital system, it
ought to be evaluated for malignancy. Software evaluation methods are
categorized as static, dynamic, and hybrid [14]; hybrid analysis utilizes both
static and dynamic analysis.

Dynamic analysis of software is performed by executing it in a ‘ _sandbox_ ’
(an isolated environment that mimics end-user operating environments),
observing and recording its behavior. The recorded behavior of the software is
evaluated manually or by machine learning techniques for final classification
[3]. Dynamic analysis is prolonged, expensive, and time-consuming and thus is
impractical for most cases.

Static analysis of software is performed by extracting features present in the
software by analysis or mining and then classifying it based on the features
obtained. The features needed for classification depend on the tool to be
used. Five categories of tools for static analysis are signature-based,
behavior-based, heuristic-based, model checking-based, and machine learning
based [2].

The malware developers have access to all detection tools, and they swiftly
modify their code to evade one or more of these tools. Thus, malware detection
tools must be ready to identify newer variants or the arrival of a new malware
category. Unsurprising, the race between malware developers and malware
detection tool developers is expected to continue forever. Fortunately,
_machine learning_ (ML) based tools often successfully identify newer variants
and _learn to detect new malware category, by _incrementally_ retraining the
model after adding samples of newer variants and categories in the training
dataset_.

###  1.1 Machine Learning for Malware Detection

Signature-based malware detection tools’ malware detection rate decreased
because advanced malware lacked signature strings. Search for alternative
detection techniques converged to ML-based tools. The first ML-based tools
were reported in 2001 [13]. Since then, numerous ML-based models have been
reported and they have been compared and contrasted in many survey papers. The
most recent and comprehensive of these are [6, 2, 15].

Until recently, the ML-based models were developed using small datasets, and
their validation did not use unseen _future_ malware. The systems developed
and evaluated using the EMBER dataset [1] overcame both of these problems
because the dataset is large, data is divided into subsets based on malware’s
observation time; the training dataset is from the _past_. Moreover, the
evaluation dataset is from the _future_ , that is, each PE file in the
evaluation dataset was detected after all PE files in the training dataset
were detected. Since the EMBER data set has a time-based partition, work
reported in [5] evaluated deterioration over time.

In this work, first, we briefly describe the structure and elements of
Portable Executable (PE) files, review a benchmark dataset, and the extracted
feature sets provided with the dataset (see Table 2 for details). Also, we
briefly describe ML algorithms that have been proposed and evaluated in the
past, and then report the performance of ML-based systems using EMBER dataset.
Finally, we experimentally evaluate the contributions of the nine feature sets
towards the detection of malware.

##  2 Portable Executable Files

The portable executable (PE) file format is a windows based file format for
executable files based on the UNIX Common Object File Format. The PE format
contains instructions for the windows dynamic linker on loading the code
stored in the file and necessary libraries to execute the code into memory. A
few common PE file extensions include .exe, .dll, .sys, .cpl, and .tsp. A
simplified layout of the Windows PE file is shown in Figure 1.

PE files consist of two parts, the PE header, and the PE sections. The PE
sections contain the raw data associated with the PE file. While sections can
take any name there are a few common section names found in most PE files.
Common sections include:

  * •

.text - which includes the file’s code,

  * •

.data - which includes initialized data used by the file,

  * •

.rdata - which includes read-only data to be used by the file

  * •

.bss -which contains uninitialized data.

The PE header contains instructions on where libraries are located and how to
load the file’s code into memory. Most PE headers will include section
headers, a DOS Header, and Optional Headers. The Optional Header is required
for all executable files and contains the Data Directories where necessary
libraries are stored, the Import table, which contains the functions imported
by the PE file, and the export table, which contains functions that other PE
files can import.

Because PE files can be run on most versions of windows, they are the most
popular vector for malware. The malware detection site VirusTotal reported in
2016 that a plurality of files scanned, 47.8%, were PE files [7]. Since the PE
file format is a popular vector for malware, researchers have assembled
datasets containing both malicious and benign PE file samples.

##  3 The EMBER Dataset

The Elastic Malware Benchmark for Empowering Researchers (EMBER) dataset [1]
is a dataset consisting of preprocessed malicious and benign PE files. The
dataset contains 900,000 training samples and 200,000 testing samples. The
training set consists of 300,000 samples of malicious, benign, and unlabeled
files, while the testing set consists of 100,000 samples of malicious and
benign PE files.

The EMBER dataset v1 was initially released in 2017, followed in 2018 with
EMBER dataset v2. The most significant change is that EMBER v2 added a new
feature set, data directory, that was not present in EMBER v1. The paper
released with the dataset used v1 while our work uses v2; thus, some of our
results may vary slightly from what was initially reported.

EMBER consists of preprocessed PE files and software for vectorizing the
feature sets found in each file. The length of each feature can be found in
Table 2, which sums to 2381. Below is a brief description of these nine
feature sets. For the convenience of referring to these features, we
abbreviate each with two letters as shown within parenthesis.

  1. 1.

Byte histogram (BH) - The histogram of Byte distribution in the file

  2. 2.

Byte-entropy histogram (BE) - The histogram of the byte entropy in the file

  3. 3.

String information (ST) - Information about the strings found in the file

  4. 4.

General file information (GE) - General information about the file such as the
size of the file

  5. 5.

Header information (HE) - Information found in the section header of the file
such as number of sections and size of sections

  6. 6.

Section information (SE) - Information from each section, this feature set is
generated by using the Hashing Trick [4]

  7. 7.

Imported functions (IM) - The functions that the file calls, this feature set
is generated by using the Hashing Trick on the library calls and then on the
function calls

  8. 8.

Exported functions (EX) - The functions that the file exports, this feature
set is generated using the Hashing Trick

  9. 9.

Data Directory (DD) - The directories where critical files needed to run the
portable executable can be found.

The Byte Histogram, Byte-Entropy Histogram, and String Information are not
unique to the PE file format and could be used to make a more general malware
detection model. The EMBER dataset includes code to vectorize the data. We
will use the vectorized data in our models. The process for vectorization is
described in detail in the EMBER paper [1]. Several feature sets use the so
called ‘Hashing Trick’ for feature extractions, the fact is that it is not _a
trick_ at all; it is a well-studied dimension reduction technique with a solid
theoretical foundation [4].

Table 2 shows the description of each of the feature sets after vectorization.
Some of the feature sets, such as sections, have a wide range, while others,
such as byte histogram, have a narrow range. To account for this, all of our
models use feature scaling to keep feature sets in similar ranges and prevent
one feature set from dominating the models. Likewise, the length of each
feature set varies. We will look at methods for reducing the number of feature
sets to see if a smaller dataset results in better performance.

Feature | Abbr. | Feature | Abbr.
---|---|---|---
Histogram | BH | Sections | SE
ByteEntropy | BE | Imports | IM
Strings | ST | Exports | EX
General | GE | DataDirectory | DD
Header | HE |  |

Table 1: Abbreviations used for each feature set

Feature | Data | Minimum | Maximum | Size
---|---|---|---|---
| Type | value | value |
BH | Float | 0 | 0.9999558 | 256
BE | Float | 0 | 0.9997657 | 256
ST | Int+ | 0 | 61745440 | 104
GE | Int+ | 0 | 8589935000 | 10
HE | Int | -2 | 4294967300 | 62
SE | Int | -4294967300 | 6167139300 | 255
IM | Int | -882858 | 83917 | 1280
EX | Int | -2288 | 308 | 128
DD | Int+ | 0 | 4294967300 | 30

Table 2: Data type, minimum and maximum values, and size of each feature set
extracted from a PE file.

##  4 ML-Based Malware Detection

After a brief introduction to categories of ML algorithms, we review ML-
algorithms that have been used for developing ML-based malware detection
systems.

###  4.1 Categories of ML-algorithms

A Machine Learning (ML) algorithm (Readers not familiar with ML algorithms can
read them from one of the many ML books, for example [9]) creates an
intelligent system from sample data, known as _training data_ , presented to
it; the created system is evaluated using _unseen_ sample data, known as
_evaluation data_. A sample datum marked with the identity of the class where
it belongs is known as _labeled datum_. For simplifying the discussion, we
assume that a software file is either malicious or benign and thus, a
software’s label is either malware or benign. Unless otherwise stated, we
assume that a trained system classifies an input file as either malware or
benign.

Based on the types of training data used, ML algorithms are divided into three
categories: _unsupervised_ , _supervised_ , and _semi-supervised_. An
unsupervised ML algorithm is trained with unlabeled data. During training, the
ML algorithm finds parameter values of the model for clustering the training
data into two clusters: one cluster is for malware, and another cluster is for
benign software; simply stated, an unsupervised ML algorithm uses a distance
measure to determine cluster membership (see [14] for different distance
measures that have been used for training malware detection systems). When an
unseen datum is presented to a trained unsupervised system, it classifies it
to the cluster whose distance is minimum.

Algorithms in the supervised category are presented with _labeled_ training
data; that is, each sample has a class label. We use supervised models in this
work. The semi-supervised category of algorithms is presented with some
labeled and unlabeled data; this category’s definition varies in the
literature. Most of the successful ML-based malware detection tools are
trained with _extensive_ training datasets.

###  4.2 Algorithms Used for ML-based Detection Tools

First ML-based malware detection tools, proposed in 2001, used _Rule-based
classifier_ and _Naive Bayes_ learning algorithms [13]; subsequently,
potentially all machine learning algorithms have been used for creating ML-
based malware detection tools [14, 15, 2, 8, 6]. Among the survey papers, [6]
is most comprehensive; it has tabulated ML algorithms in three tables — one
each for supervised, semi-supervised, and unsupervised categories. Also, the
tables have sources of the datasets, their sizes, and features of the datasets
used for training and evaluation of the malware detection systems.

From these tables, we find that very few datasets were in the public domain.
Many of them are no longer available today, and most of them are too small for
developing any real-world malware detection system. Moreover, with one
exception [11] we suspect that these datasets _did not_ include the discovery
date of the malware; thus, the developed models were not evaluated for
endurance with evolving malware variants and new varieties of malware.

###  4.3 Previously used ML-algorithms

Machine learning algorithms used for training malware detection systems
include Artificial Neural Networks, Bayes classifier, Bayesian Networks,
Belief Propagation, different Clustering, Decision trees (DT), Gradient
Boosting Decision Trees (GBT), $k$-Nearest Neighbors ($k$-NN / KNN), Linear
Discriminant Analysis (LDA), Logistic Model Tree, Logistic Regression,
Multiple Kernel Learning, Naive Bayes (NB), Prototype-based Classification,
Quadratic Discriminant Analysis (QDA), Random Forest (RF), Rule-based
classifiers, Self-Organizing Maps, and Support Vector Machine (SVM).
Interestingly, 75% of the systems were developed using four ML algorithms — DT
(31%), SVM (29%), NB (10%), and RF (5%) [14].

Most recently, ANN, including Deep Neural Networks (DNN), and a variant of GBT
known as Light Gradient Boosting Machine (LGBM), have shown good performance
when trained and evaluated using a large dataset [1], EMBER (described
earlier).

ML-based systems developed and evaluated earlier have provided significant
insights, and these insights are valuable. However, the reported results
cannot be compared with each other because they were developed and evaluated
using myriad and relatively small datasets. Availability of the large EMBER
dataset makes a compelling case for developing and evaluating all those ML-
based systems again.

To that end, using Version 0.23.2 of Scikit-learn: Machine Learning in Python
package [10] and the EMBER dataset, we developed and evaluated malware
detection systems.

##  5 Performance of ML-based Systems Developed using EMBER dataset

Our study includes LGBM (original and normalized data), RF, LG, LDA, KNN, QDA,
SVM, and ANN/DNN ML-algorithms; all 9 feature sets of the EMBER dataset were
used to develop these systems.

Their malware detection rates are shown in Table 3. (The program codes will be
made available from our GitHub once the paper is accepted for publication.)
Because the paper describing the EMBER dataset (benchmark) reported that LGBM
has the best detection rate, we consider the performance of LGBM as the
baseline performance. The table entries are sorted from the highest accuracy
to the lowest accuracy (except ANN/DNN). The accurate detection rate for the
RF is 1.67% below LGBM, which is not a surprise because LGBM is an optimized
version of RF. The detection accuracy of LG, LDA, and KNN systems is 87.37%,
85.82%, and 82.23%, respectively. The systems developed using SVM, and QDA
performed very poorly — 54.18% and 50.29%, respectively; these are almost
identical to a coin toss.

| ACC | AUC | Presision | Recall | F1
---|---|---|---|---|---
LGBM |  |  |  |  |
Baseline | 0.9363 | 0.9845 | 0.9244 | 0.9504 | 0.9372
LGBM |  |  |  |  |
Normalized | 0.9351 | 0.9844 | 0.9267 | 0.9450 | 0.9357
RandomForest | 0.9196 | 0.9196 | 0.9208 | 0.9183 | 0.9195
Logistic |  |  |  |  |
Classifier | 0.8737 | 0.8737 | 0.8394 | 0.9243 | 0.8798
LDA | 0.8582 | 0.8582 | 0.8175 | 0.9222 | 0.8667
KNN | 0.8223 | 0.8223 | 0.7979 | 0.8632 | 0.8293
SVM | 0.5418 | 0.5418 | 0.5222 | 0.9853 | 0.6826
QDA | 0.5029 | 0.5029 | 0.9770 | 0.0059 | 0.0118
ANN/DNN | 0.9522 | 0.9815 | 0.9592 | 0.9445 | 0.9518

Table 3: Performances of ML-based systems on the EMBER dataset, from the
highest to the lowest ACCs. The AUC is reported for a false-positive rate less
than 1%.

As can be seen from the last row of the Table 3, the performance of our
proposed ANN/DNN system is better than the performance of the baseline LGBM-
based system. We describe our ANN/DNN architecture and compare its
performances with the performances of the baseline system.

###  5.1 Neural Network Approach to Malware Detection

We extensively evaluated various ANN/DNN architectures. We found that large
deep networks do not increase detection rate above the baseline LGBM, but they
needed significantly higher training time and slightly higher detection time.
We found a simple ANN architecture (shown in Fig. 2) that outperformed all
systems. As can be seen from the figure, it consists of three dense layers of
512, 128, and 8 neurons with _tanh_ activation functions. The output layer
consists of a two-unit $softmax$ layer. A batch normalization layer is used
before feeding the inputs to the network; also, we use another batch
normalization layer between the 2nd and 3rd dense layers. Performance of this
ANN-based system is shown in the last row; we can see that it improves
accurate detection rate from 93.363% to 95.22% – a 1.27% increase; but if one
considers the reduction of malware detection rate, the ANN-based system
reduces inaccurate detection rate to 4.78% from LGBM’s 6.05% — a reduction of
20.99%.

##  6 Evaluation of the Nine Feature Sets of the EMBER Dataset

Multiple features from PE files have been used in [12] and [1]. We were
interested to see how many contributions each of these features make towards
the detection rates. Because ANN and LGBM based systems have two top detection
rates, we only report our observations for those systems. It is important to
note that the length of the input to a system is the sum of lengths of the
feature sets used to train it.

###  6.1 Single Feature Set

Performance of the ANNs and LGBMs with only one feature set is shown in Table
4. First, two columns of the table are for accuracies (ACCs) and AUCs (AUC
stands for area under the curve and is the area under the receiver operator
characteristic (ROC) curve.) at a false-positive rate less than 1%, and the
3rd and the 4th columns are for those of the LGBM. While most of the features
show similar detection rates, there are a couple of notable differences.

| ANN Performance | LGBM Performance
---|---|---
| ACC | AUC | ACC | AUC
BH | 0.8949 | 0.9567 | 0.8657 | 0.9465
BE | 0.8912 | 0.9566 | 0.8529 | 0.9368
ST | 0.8859 | 0.9511 | 0.8645 | 0.9462
IM | 0.8662 | 0.9338 | 0.8300 | 0.9172
SE | 0.8648 | 0.9359 | 0.8757 | 0.9522
GE | 0.7316 | 0.7924 | 0.8285 | 0.9243
HE | 0.7312 | 0.8088 | 0.8102 | 0.9110
EX | 0.5571 | 0.5852 | 0.5531 | 0.5610
DD | 0.5050 | 0.4852 | 0.8473 | 0.9356

Table 4: The results of our ANN and the LGBM systems with only a single
feature set at a time. The AUC is for a false-positive rate less than 1%.

The ANN-based system has the highest detection rate of 89.49% with BH (byte
histogram), and that for the LGBM based system is 87.57 with SE (sections).
For most _single_ feature sets, we observe that performances of ANN and LGBM
systems are close to each other, with the only exception: DD (data directory);
ANN-based system performed very poorly (at 50.50%), but LGBM based system
achieved an 84.73% accuracy. Another notable difference is HE (header
information), where LGBM had a detection rate of 81.02% vs. 73.12% for the
ANN.

Next we present our observations for combinations of two feature sets at a
time.

| ANN Preformance | LGBM Preformance
---|---|---
| Acc | AUC | Acc | AUC
BH_SE | 0.9342 | 0.9745 | 0.9129 | 0.9737
BH_IM | 0.9313 | 0.9751 | 0.9064 | 0.9670
BE_IM | 0.9292 | 0.9702 | 0.9044 | 0.9686
BE_SE | 0.9275 | 0.9724 | 0.8962 | 0.9683
ST_SE | 0.9255 | 0.9718 | 0.9161 | 0.9762
$\vdots$ | $\vdots$ | $\vdots$ | $\vdots$ | $\vdots$
GE_EX | 0.7372 | 0.8031 | 0.8315 | 0.9233
HE_EX | 0.7304 | 0.8087 | 0.8107 | 0.9098
HE_DD | 0.7274 | 0.8083 | 0.8621 | 0.9489
GE_DD | 0.7203 | 0.7841 | 0.8796 | 0.9539
EX_DD | 0.6157 | 0.6367 | 0.8470 | 0.9341

Table 5: The results of our ANN and the LGBM systems with combinations of two
feature sets. The AUC is reported for a false-positive rate less than 1%.

###  6.2 Combination of Two Feature Sets

Out of nine feature sets, one can choose two features simultaneously in 36
distinct possible ways. We evaluated all 36 possible combinations of two
feature sets, but we show the top five and bottom five pairwise combinations
to conserve space and reduce clutter. It is important to note that the set of
included two-feature set combinations has the best performances of both ANN
and LGBM based systems. Our ANN-based system has the highest malware detection
rate of 93.42% with the feature sets BH (byte histogram) and SE( Section
information); the LGBM based system has the highest malware detection rate of
91.61% with the features ST (String information) and SE — both share the
feature set SE.

While differences in malware detection rates in the top five observations vary
a little, the picture for the bottom five is different; our ANN-based systems’
detection rate is drastically lower than that of the LGBM based systems’
detection rate. Our ANN-based systems have the lowest detection rate of 61.57%
with feature sets EX (exports) and DD (Data directories), but the lowest
detection rate for the LGBM based system is 84.70% — a difference of more than
23%; fortunately, there is no good reason to choose the worst-performing
combinations.

The top five pairwise combinations include features BH, BE, SE, IM, and ST.
Looking back at performances of single feature set based systems, these five
features produced the top five systems (see Table 4). Comparing Tables 4 and
5, we can find that BH and BE produced two top-performing systems with single
feature sets, but the combination of these two feature sets did not produce
any of the top-performing five systems for ANNs.

To create and evaluate systems with three or more features, we will consider
only BH, BE, SE, IM, and ST. Next, we discuss the performances of systems
created with three feature sets.

| ANN Preformance | LGBM Preformance
---|---|---
| Acc | AUC | Acc | AUC
BH_SE_IM | 0.9471 | 0.9810 | 0.9262 | 0.9800
BE_SE_IM | 0.9458 | 0.9799 | 0.9228 | 0.9770
ST_SE_IM | 0.9447 | 0.9803 | 0.9276 | 0.9809
BH_BE_SE | 0.9370 | 0.9754 | 0.9181 | 0.9766
BH_ST_SE | 0.9352 | 0.9749 | 0.9201 | 0.9784
BE_ST_SE | 0.9345 | 0.9740 | 0.9194 | 0.9781
BE_ST_IM | 0.9340 | 0.9729 | 0.9154 | 0.9740
BH_BE_IM | 0.9334 | 0.9749 | 0.9134 | 0.9728
BH_ST_IM | 0.9317 | 0.9735 | 0.9118 | 0.9738
BH_BE_ST | 0.9164 | 0.9630 | 0.8983 | 0.9639

Table 6: The results of our ANN and the LGBM systems with combinations of
three feature sets. Only Byte Histogram, Byte Entropy, Strings, Sections, and
Imports were tested. The AUC is reported for a false-positive rate less than
1%.

###  6.3 Combination of Three Features

Performances of systems created with three out of five features are shown in
table Table 6; we have included results for all ten possible combinations.
Among all the LGBM based systems created and evaluated with combinations of
three features, the system created with ST, SE, and IM has the highest malware
detection rate of 92.75% (this is less than 1% lower than the baseline system
with 93.63% detection rate). The most interesting observation is the
performances of the ANN-based systems. The top four systems outperform the
baseline system’s performance, and the best detection rate of 94.71% is
obtained by our ANN-based system created with BH, SE, and IM. Another
interesting observation is that SE is the common feature among all the top six
ANN-based systems.

| ANN Preformance | LGBM Preformance
---|---|---
| Acc | AUC | Acc | AUC
ALL | 0.9522 | 0.9815 | 0.9363 | 0.9845
BE_ST_SE_IM | 0.9500 | 0.9808 | 0.9310 | 0.9827
BH_BE_ST_SE_IM | 0.9497 | 0.9818 | 0.9339 | 0.9833
BH_ST_SE_IM | 0.9470 | 0.9797 | 0.9332 | 0.9830
BH_BE_SE_IM | 0.9464 | 0.9817 | 0.9297 | 0.9808
BH_BE_ST_SE | 0.9381 | 0.9745 | 0.9234 | 0.9799
BH_BE_ST_IM | 0.9362 | 0.9749 | 0.9183 | 0.9761

Table 7: The results of our ANN architecture and the LGBM against quadruplets
of feature sets, quintuplets of feature sets, and all of the feature sets.
Only Byte Histogram, Byte Entropy, Strings, Sections, and Imports were tested
for the quadruplets and quintuplets.

###  6.4 Combinations of Four or More Feature Sets

Results of systems created with four or more features are shown in Table 7. As
discussed and analyzed earlier, our ANN-based systems with all features have a
malware detection accuracy of 95.22% and significantly outperform the baseline
systems created with LGBM. Another observation is that our ANN-based systems
have outperformed the baseline LGBM based system.

##  7 Conclusion and Future Work

Availability of the large EMBER dataset [1] has opened the door for creating
and evaluating malware detection systems. The dataset includes more than just
the PE files; it includes the program code necessary to extract nine features
from each PE file and the program code for LGBM; we have used these programs
to obtain the baseline results we have reported. The program code for
extracting features from each PE file produces a long vector, concatenating
all nine features. We separated each of those features for our study reported
here.

We have seen that ANN-based systems can outperform other ML-based systems; we
have shown that all features do not contribute equally. Even when two features
perform well individually, they may not have any synergy; together, they may
not produce a system that performs much better than the systems created with
each feature individually; our study revealed that Byte histogram and Byte
entropy are such a pair of features.

Another interesting observation is that two ML algorithms may create systems
with the same set of features, but their performances vary widely; we have
seen that for our ANN-based and LGBM based systems.

## Declarations

### Availability of data and materials

Our models are available upon request.

### Competing interest

The authors declare that they have no competing interests.

### Funding

None.

### Authors’ contributions

The first author conceived the idea of the study and wrote the code and the
paper; both authors discussed the results and revised the final manuscript.
Both authors read and approved the final manuscript.

### Acknowledgements

Not applicable.

## References

  * [1] Hyrum S. Anderson and Phil Roth.  EMBER: An open dataset for training static PE malware machine learning models, 2018.
  * [2] Omer Aslan Aslan and Refik Samet.  A comprehensive review on malware detection approaches.  IEEE Access, 8:6249–6271, 2020.
  * [3] Ivan Firdausi, Charles lim, Alva Erwin, and Anto Satriyo Nugroho.  Analysis of machine learning techniques used in behavior-based malware detection.  In 2010 Second International Conference on Advances in Computing, Control, and Telecommunication Technologies, pages 201–203, 2010\.
  * [4] Casper Freksen, Lior Kamma, and Kasper Green Larsen.  Fully understanding the hashing trick.  In Proceedings of the 32nd International Conference on Neural Information Processing Systems, NIPS’18, pages 5394–5404, Red Hook, NY, USA, 2018. Curran Associates Inc.
  * [5] Colin Galen and Robert Steele.  Evaluating performance maintenance and deterioration over time of machine learning-based malware detection models on the EMBER PE dataset.  In 2020 Seventh International Conference on Social Networks Analysis, Management and Security (SNAMS), pages 1–7, 2020.
  * [6] Daniel Gibert, Carles Mateu, and Jordi Planes.  The rise of machine learning for detection and classification of malware: Research developments, trends and challenges.  Journal of Network and Computer Applications, 153:102526, 2020.
  * [7] Ajit Kumar, K.S. Kuppusamy, and G. Aghila.  A learning model to detect maliciousness of portable executable using integrated feature set.  Journal of King Saud University - Computer and Information Sciences, 31(2):252–265, 2019.
  * [8] Hoda El Merabet and Abderrahmane Hajraoui.  A survey of malware detection techniques based on machine learning.  International Journal of Advanced Computer Science and Applications, 10(1), 2019.
  * [9] Tom M. Mitchell.  Machine Learning.  McGraw-Hill, New York, 1997.
  * [10] F. Pedregosa, G. Varoquaux, A. Gramfort, V. Michel, B. Thirion, O. Grisel, M. Blondel, P. Prettenhofer, R. Weiss, V. Dubourg, J. Vanderplas, A. Passos, D. Cournapeau, M. Brucher, M. Perrot, and E. Duchesnay.  Scikit-learn: Machine learning in Python.  Journal of Machine Learning Research, 12:2825–2830, 2011.
  * [11] Joshua Saxe and Konstantin Berlin.  Deep neural network based malware detection using two dimensional binary program features.  In 10th International Conference on Malicious and Unwanted Software, MALWARE 2015, Fajardo, PR, USA, October 20-22, 2015, pages 11–20. IEEE Computer Society, 2015.
  * [12] Joshua Saxe, Richard Harang, Cody Wild, and Hillary Sanders.  A deep learning approach to fast, format-agnostic detection of malicious web content, 2018.
  * [13] M.G. Schultz, E. Eskin, F. Zadok, and S.J. Stolfo.  Data mining methods for detection of new malicious executables.  In Proceedings 2001 IEEE Symposium on Security and Privacy. S P 2001, pages 38–49, 2001.
  * [14] Alireza Souri and Rahil Hosseini.  A state-of-the-art survey of malware detection approaches using data mining techniques.  Hum.-Centric Comput. Inf. Sci., 8(1), December 2018.
  * [15] Daniele Ucci, Leonardo Aniello, and Roberto Baldoni.  Survey of machine learning techniques for malware analysis.  Computers & Security, 81:123–147, 2019.

## Figures

Figure 1: A simplified view of the notepad.exe PE file to highlight the
anatomy of the PE file format. Figure 2: A simple but effective ANN
architecture for malware detection,


#  An AI-Based System Utilizing IoT-Enabled Ambient Sensors and LLMs for
Complex Activity Tracking

Yuan Sun  WINLAB, Rutgers University Piscataway, NJ,  USA
<EMAIL_ADDRESS>and  Jorge Ortiz  WINLAB, Rutgers University
Piscataway, NJ,  USA<EMAIL_ADDRESS>
###### Abstract.

Complex activity recognition plays an important role in elderly care
assistance. However, the reasoning ability of edge devices is constrained by
the classic machine learning model capacity. In this paper, we present a non-
invasive ambient sensing system that can detect multiple activities and apply
large language models (LLMs) to reason the activity sequences. This method
effectively combines edge devices and LLMs to help elderly people in their
daily activities, such as reminding them to take pills or handling emergencies
like falls. The LLM-based edge device can also serve as an interface to
interact with elderly people, especially with memory issue, assisting them in
their daily lives. By deploying such a system, we believe that the smart
sensing system can improve the quality of life for older people and provide
more efficient protection.

††copyright: acmlicensed††journalyear:  ††doi: XXXXXXX.XXXXXXX††conference:
Make sure to enter the correct conference title from your rights confirmation
emai; ; ††isbn:

##  1\. Introduction

Non-intrusive sensors are crucial for modern sensing applications,
particularly in fields requiring continuous and unobtrusive monitoring, such
as elderly care(lentzas2020non, 11). These sensors, which do not rely on
cameras, offer significant advantages, including enhanced comfort and privacy
for users. By seamlessly integrating into the environment without requiring
direct interaction or visible placement, they minimize aesthetic and social
intrusiveness(chidurala2022iot, 7). Additionally, non-intrusive sensors reduce
deployment and maintenance efforts, leveraging existing infrastructure and
eliminating the need for frequent battery replacements(yu2013performance, 17).
This ease of installation and low maintenance makes them practical for long-
term use. Furthermore, they provide flexibility and scalability in sensing
applications, capturing a broad range of environmental data indirectly, thus
enabling comprehensive monitoring across diverse contexts.

In this work, we build a non-intrusive smart sensing system that relies on the
reasoning ability of large language models (LLMs) to assist in elderly care.
The system can detect complex activities composed of more than two atomic
activities. An atomic activity refers to unit-level activities that can be
captured by sensors within a short time window and cannot be broken down
further. Besides detecting normal atomic activities, we also use LLMs for
high-level explanations and reasoning. To adapt to the memory constraints of
our edge devices, we employ a local inference model on the IoT devices before
sending the data to the LLM. This approach reduces the transmission burden via
the wireless network, ensuring efficient and real-time processing.

##  2\. LLM and Complex Reasoning

Our system first collects non-intrusive sensor data from ambient sensors.
These sensors detect atomic-level activities. Due to the limited memory
capacity of the sensors, a small model runs locally to perform initial
inferences. Deploying a small model locally allows for real-time processing,
reducing latency and ensuring immediate response to critical activities.
Additionally, local processing minimizes the amount of data transmitted to the
cloud, enhancing privacy and reducing bandwidth usage. Once an atomic activity
is detected, it is sent to a cloud server where an LLM is deployed. The LLM
applies its reasoning and interaction capabilities with the user. For
instance, if the sequence of activities indicates eating and drinking without
sanitizing, the system can remind the user to wash their hands. Another
example is if the user drinks water but forgets to take their pills; the LLM
will interact with the user to remind them to take their medication.
Additionally, if a user is detected as getting dressed but skipping a critical
item like a coat on a cold day, the LLM can prompt them to wear it, ensuring
their well-being. Furthermore, in the event of a medical emergency, the system
can recognize distress signals or calls for help, and promptly alert emergency
services, providing vital assistance when needed.

Figure 1. The system employs both local inference using ambient sensors and
reasoning via a cloud-based LLM. The sensors detect atomic activities, and the
cloud server receives these activity sequences as context to further detect
higher-level meanings, make decisions, and interact with the user.

##  3\. sensor board setup

Our goal is to construct a device that enhances the generalizability of sensor
data explanation while safeguarding user privacy. The device (figure 3) ,
installed in the environment, utilizes non-invasive sensors to detect human
activities. These sensors include PIR (motion), IMU (accelerometers), audio,
RGB, pressure, humidity, magnetometer, gas, and temperature sensors. They are
widely available and well-suited for generalizing smart space sensing
capabilities.

Smart space sensors produce two types of outputs: binary and continuous.
Binary sensors, like PIR motion sensors, provide ”on” and ”off” states.
Although useful for indicating environmental status, they may be inadequate
for analyzing complex activities (gochoo2017dcnn, 8, 2). Our device uses PIR
sensors for motion detection, transmitting ”ON” and ”OFF” signals. In
contrast, accelerometers offer continuous outputs, tracking movement along the
X, Y, and Z axes. We utilize them to detect environmental vibrations, such as
different vibration frequencies for pouring water versus running
(ahmed2017improving, 1, 16, 14).

The RGB sensor detects colors within a specific range, identifying the
presence of colors based on the activated channels. Pressure sensors measure
barometric pressure, providing insights into weather patterns and indoor air
quality, especially in HVAC operations. Magnetometers detect magnetic fields,
but readings can be influenced by nearby electronic devices and human activity
due to RF electromagnetic fields (martinez2004comparison, 12, 5). Gas sensors
detect substances like ethanol, alcohol, and carbon monoxide, useful for
monitoring air quality changes due to activities such as heavy running.

Figure 2. The non-intrusive sensor board we design for our system

The Raspberry Pi Model B+ (Figure 3) serves as the central processing unit,
interfacing with a custom sensor board via GPIO pins. The sensor board,
equipped with various ambient sensors, collects data such as temperature,
humidity, and motion. The Raspberry Pi, powered by a Broadcom BCM2837B0 quad-
core ARM Cortex-A53 64-bit SoC at 1.4GHz with 1GB LPDDR2 SDRAM, runs a Python
script to read sensor values at regular intervals, process the data, and store
it for further analysis. This setup, featuring dual-band 802.11ac wireless
LAN, Bluetooth 4.2, Gigabit Ethernet, and multiple USB ports, enables real-
time monitoring and data collection essential for elderly care assistance and
activity tracking. Although our current sensor setup is insufficient to run a
model for data embedding locally, we plan to implement local atomic activity
detection in the next stage of this work.

Figure 3. Raspberry Pi Model B+ used in our ambient sensor setup, facilitating
seamless integration for elderly care assistance and activity tracking

##  4\.  initial experiment

In our initial experiment, we will use the ambient sensor to collect atomic
activities and test the LLM model on the cloud. The target is to determine if
it possesses the reasoning capabilities that can assist in the daily lives of
elderly individuals.

###  4.1. Data Collection

The atomic level activities we collect include 20 activities: eat, paperdis,
write, chop, hand wash, pour water, clean floor, knock, run, curtain, light
switch, type, door pass, wipe desk, chat, basketball, saw, shave, wash dish,
and brusing teeth. The sampling frequency is 90Hz for each sensor channel.
There are discussions about whether the high-fidelity audio sensor channel
(laput2017synthetic, 10) is the main channel for detecting patterns. Audio
data can also be intrusive (bandodkar2014non, 3, 15). To avoid privacy issues
associated with the audio channel, its frequency is decreased from 16kHz to
90kHz. According to Shannon’s sampling theorem, information contained in
frequencies greater than half the sampling rate cannot be recovered
(niazi2017statistical, 13, 9, 4). Figure 4 visualizes the sensor data
collected during eating activities. Figure 5 visualizes the data collected
during chopping activities.

Figure 4. Eating activity collected by the ambient sensor Figure 5. Data
collected from chopping activities

###  4.2. Initial Results

We developed a data encoder to process the sensor data before sending it to
the LLM. Our model primarily uses a channel-wise MLP to extract features from
each channel, followed by a sensor fusion module also based on MLP.
Additionally, we incorporated a fast Fourier convolution module (chi2020fast,
6) to improve the accuracy rate.

Table 1 presents the initial results of F1 scores, precision, and recall for
detecting various activities using the model. Most activities, such as ’eat’,
’write’, ’chop’, ’hand wash’, ’clean floor’, ’knock’, ’run’, ’curtain’, ’light
switch’, ’type’, ’door pass’, ’wipe desk’, ’basketball’, ’saw’, ’shave’, ’wash
dish’, and ’teeth’, achieve perfect scores of 1.00 across all metrics.
However, ’paperdis’, ’pour water’, and ’chat’ show lower performance with F1
scores of 0.43, 0.48, and 0.32, respectively. These variations indicate that
while the model performs exceptionally well for most activities, there is room
for improvement in detecting certain activities.

We also demonstrate the robustness of our model when the noise level 7 and
frequency6 of the data vary.

Figure 6. Model robustness in detecting data at different frequencies Figure 7. Model robustness with varying levels of noise in the sensor data Activity Name | F1 | Precision | Recall
---|---|---|---
eat | 1.00 | 1.00 | 1.00
paperdis | 0.43 | 0.50 | 0.38
write | 1.00 | 1.00 | 1.00
chop | 1.00 | 1.00 | 1.00
hand wash | 1.00 | 1.00 | 1.00
pour water | 0.48 | 0.50 | 0.47
clean floor | 1.00 | 1.00 | 1.00
knock | 1.00 | 1.00 | 1.00
run | 1.00 | 1.00 | 1.00
curtain | 1.00 | 1.00 | 1.00
light Switch | 1.00 | 1.00 | 1.00
type | 1.00 | 1.00 | 1.00
door pass | 1.00 | 1.00 | 1.00
wipe desk | 1.00 | 1.00 | 1.00
chat | 0.32 | 0.33 | 0.31
basketball | 1.00 | 1.00 | 1.00
saw | 1.00 | 1.00 | 1.00
shave | 1.00 | 1.00 | 1.00
wash dish | 1.00 | 1.00 | 1.00
teeth | 1.00 | 1.00 | 1.00
Table 1. Initial results of F1 scores for each activity

###  4.3. Applying LLMs to Sequence Detection and Complex Activity Recognition

After detecting atomic activities, we send the sequence of detected activities
to the LLM. For example, to remind elderly people to take medication before
eating, we send a sequence of activities such as ”brushing teeth → hand wash →
pour water → eat”. The corrected sequence from the LLM is ”brushing teeth →
hand wash → take medication → pour water → eat”, recognizing the complex
activity as ”forgetting medication”. If the original order detected is ”eat →
basketball → brush teeth”, the LLM will suggest ”brush teeth → wash hands →
eat → basketball” to ensure proper hygiene before eating, recognizing the
complex activity as ”unhygienic behavior”. We use GPT-4 to implement sequence
verification and reminders. For instance, if the input sequence is ”door pass
→ using paper dispenser”, the LLM suggests ”door pass → turn the switch →
paper dispenser” to remind the user to turn on the light, identifying the
complex activity as ”preventing slipping”. The order rules can be configured
in the prompt to meet detailed requirements. Our current results demonstrate
an implementable framework that can help elderly people, especially those with
memory issues, to live a high-quality life.

##  5\. Conclusion and Future Work

In this work, we designed a framework leveraging IoT and LLMs to assist
elderly people and enhance their quality of life. We employed a board with
non-intrusive sensors to avoid discomfort for the elderly while monitoring
their behavior to provide assistance. We implemented a low-frequency sensor
sampling strategy, particularly for the microphone channel, using an
unrecoverable sampling rate that cannot reproduce the original sound
information. A smart IoT box was built to collect initial data for
experimentation. Our initial results demonstrate that our model can
successfully identify daily activities from the sensor data. We further
utilized the GPT-4 interface to test the reasoning and complex activity
detection capabilities, showcasing a promising framework that can assist
elderly individuals.

However, the current smart sensor setup uses an older version of the Raspberry
Pi. We plan to use a higher version of the Raspberry Pi to enable real-time
inference on the edge device. Additionally, we will conduct user studies with
more participants to interact with our system, helping us evaluate the design
and implementation of the system.

## References

  * (1) Hamad Ahmed and Muhammad Tahir  “Improving the accuracy of human body orientation estimation with wearable IMU sensors”  In _IEEE Transactions on instrumentation and measurement_ 66.3 IEEE, 2017, pp. 535–542
  * (2) Luca Arrotta, Gabriele Civitarese and Claudio Bettini  “DeXAR: Deep Explainable Sensor-Based Activity Recognition in Smart-Home Environments”  In _Proceedings of the ACM on Interactive, Mobile, Wearable and Ubiquitous Technologies_ 6.1 ACM New York, NY, USA, 2022, pp. 1–30
  * (3) Amay J Bandodkar and Joseph Wang  “Non-invasive wearable electrochemical sensors: a review”  In _Trends in biotechnology_ 32.7 Elsevier, 2014, pp. 363–371
  * (4) Oresti Banos et al.  “Window size impact in human activity recognition”  In _Sensors_ 14.4 MDPI, 2014, pp. 6474–6499
  * (5) Paolo Bernardi, Marta Cavagnaro, Stefano Pisa and Emanuele Piuzzi  “Specific absorption rate and temperature increases in the head of a cellular-phone user”  In _IEEE transactions on microwave theory and techniques_ 48.7 IEEE, 2000, pp. 1118–1126
  * (6) Lu Chi, Borui Jiang and Yadong Mu  “Fast fourier convolution”  In _Advances in Neural Information Processing Systems_ 33, 2020, pp. 4479–4488
  * (7) Veena Chidurala, Xiaodong Wang, Xinrong Li and Jesse Hamner  “IoT based Sensor System Design for Real-Time Non-Intrusive Occupancy Monitoring”  In _2022 IEEE International Conference on Internet of Things and Intelligence Systems (IoTaIS)_ , 2022, pp. 252–257  IEEE
  * (8) Munkhjargal Gochoo et al.  “DCNN-based elderly activity recognition using binary sensors”  In _2017 international conference on electrical and computing technologies and applications (ICECTA)_ , 2017, pp. 1–5  IEEE
  * (9) Aftab Khan, Nils Hammerla, Sebastian Mellor and Thomas Plötz  “Optimising sampling rates for accelerometer-based human activity recognition”  In _Pattern Recognition Letters_ 73 Elsevier, 2016, pp. 33–40
  * (10) Gierad Laput, Yang Zhang and Chris Harrison  “Synthetic sensors: Towards general-purpose sensing”  In _Proceedings of the 2017 CHI Conference on Human Factors in Computing Systems_ , 2017, pp. 3986–3999
  * (11) Athanasios Lentzas and Dimitris Vrakas  “Non-intrusive human activity recognition and abnormal behavior detection on elderly people: A review”  In _Artificial Intelligence Review_ 53.3 Springer, 2020, pp. 1975–2021
  * (12) M Martinez-Burdalo, A Martin, M Anguiano and R Villar  “Comparison of FDTD-calculated specific absorption rate in adults and children when using a mobile phone at 900 and 1800 MHz”  In _Physics in Medicine & Biology_ 49.2 IOP Publishing, 2004, pp. 345
  * (13) Anzah H Niazi et al.  “Statistical Analysis of Window Sizes and Sampling Rates in Human Activity Recognition.”  In _HEALTHINF_ , 2017, pp. 319–325
  * (14) Samir A Rawashdeh, Derek A Rafeldt and Timothy L Uhl  “Wearable IMU for shoulder injury prevention in overhead sports”  In _Sensors_ 16.11 MDPI, 2016, pp. 1847
  * (15) Ioan Susnea et al.  “Unobtrusive monitoring the daily activity routine of elderly people living alone, with low-cost binary sensors”  In _Sensors_ 19.10 MDPI, 2019, pp. 2264
  * (16) Xuzhong Yan, Heng Li, Angus R Li and Hong Zhang  “Wearable IMU-based real-time motion warning system for construction workers’ musculoskeletal disorders prevention”  In _Automation in construction_ 74 Elsevier, 2017, pp. 2–11
  * (17) Xin Yu and Panos D Prevedouros  “Performance and challenges in utilizing non-intrusive sensors for traffic data collection”  Scientific Research Publishing, 2013