{'fraction_non_alphanumeric': 0.037381152966662654, 'fraction_numerical': 0.015284631122878806, 'mean_word_length': 4.812255176273084, 'pattern_counts': {'":': 0, '<': 4, '<?xml version=': 0, '>': 25, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 1, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Evolutionary dynamics in an SI epidemic model with phenotype-structured susceptible compartment

We present an SI epidemic model whereby a continuous variable captures variability in proliferative potential and resistance to infection among susceptibles. The occurrence of heritable, spontaneous changes in these phenotype and the presence of a fitness trade-off between resistance to infection and proliferative potential are incorporated into the model. The model comprises an ODE for the number of infected individuals that is coupled with a partial integrodifferential equation for the population density of susceptibles through an integral term. The expression for the basic reproduction number $\mathcal{R}_0$ is derived, the disease-free and endemic equilibrium of the model are characterised and a threshold theorem is proved. Analytical results are integrated with numerical simulations of a calibrated version of the model based on the results of artificial selection experiments in a host-parasite system. The results of our mathematical study disentangle the impact of different evolutionary parameters on the spread of infectious diseases and the consequent phenotypic adaption of susceptible individuals. In particular, these results provide a theoretical basis for the observation that infectious diseases exerting stronger selective pressures on susceptibles and being characterised by higher infection rates are more likely to spread. Moreover, our results indicate that spontaneous phenotypic changes in proliferative potential and resistance to infection can either promote or prevent the spread of diseases depending on the strength of selection acting on susceptible individuals prior to infection. Finally, we demonstrate that, when an endemic equilibrium is established, higher levels of resistance to infection and lower degrees of phenotypic heterogeneity are to be expected in the presence of infections which are characterised by lower rates of death and exert stronger selective pressures.


Introduction
Mathematical models have been widely used in epidemiology to gain a deeper understanding of infectious diseases, make predictions of disease dynamics and inform possible strategies for controlling disease outbreaks [9,26,34,55].
A well-established mathematical approach to describing the evolution of infectious diseases is to use compartment models formulated in terms of systems of ordinary differential equations (ODEs), whereby each differential equation governs the dynamic of the size of a compartment [28,29,30]. This approach traditionally relies on the assumption that phenotypic characteristics are homogeneously distributed among individuals in the same compartment (e.g. susceptible individuals are all characterised by the same disease susceptibility, and all infected individuals have the same death and recovery rates). However, such a simplifying assumption is only prima facie justified and may narrow down the application domain of these epidemic models considerably. For instance, a mounting body of empirical evidence shows that there are situations in which the evolution of infectious diseases may be importantly shaped by phenotypic heterogeneity among susceptible individuals [14,20,25,41,45,47,53].
A possible way of incorporating phenotypic heterogeneity in compartment models consists in structuring one or more compartments by a continuous variable that provides a mathematical representation of the phenotype of the individuals (i.e. individuals with different phenotypes are characterised by different values of the structuring variable). The phenotypic distribution of individuals in a structured compartment is described by a population density function, the dynamic of which is governed by an integrodifferential equation (IDE) or a partial integrodifferential equation (PIDE) that replaces the ODE for the size of the corresponding unstructured compartment in the original model. Moreover, some of the parameters of the unstructured counterpart of the model are replaced by suitable functions of the structuring variable and the model equations may be coupled through integral terms.
Epidemic models that comprise structured compartments have been proposed in many different contexts, from age-structured populations [13,31,33,54] to populations occupying spatially heterogeneous environments [32,40,46,57,58]. Closer to the topic of our study, we mention some papers in which susceptibility and/or infectivity depend on a continuous structuring variable [24,42,43,44,56], leading to an IDE for the population density of the corresponding compartment.
We focus here on heterogeneity in susceptibility to infection, which we assume to depend on a structuring variable that regulates also the proliferative potential of the individuals, giving rise to a trade-off between resistance to infection and reproduction rate. Such a fitness trade-off has been experimentally shown to play a key role in the adaptation of various species to the selective pressure exerted by different forms of infection. For instance, in [27] it was reported that collared flycatcher (Ficedula albicollis) exhibits a trade-off between reproduction rate and resistance to parasites. The same kind of trade-off was documented and studied in [8], using a more theoretical approach. In [22], the author found a positive correlation, in a population of bighorn ewes (Ouis Canadensis), between resources allocated for reproduction and the presence of a specific parasite. In [37] and [50], the authors report on various other examples of species in which this trade-off is present. In [59], the existence of a fitness trade-off between resistance to a parasite (Schistosoma mansoni ) and fertility was observed in a species of freshwater snail (Biomphalaria glabrata).
The variable controlling reproduction rate and resistance to infection is subjected to evolutionary pressures that vary depending on the number of infective individuals. We further assume the occurrence of heritable, spontaneous changes in these phenotypic characteristics. Phenotypic changes are incorporated into our model via a diffusion term [23].
Models coupling epidemic processes with evolutionary changes have been studied by Burie, Ducrot and co-workers in a series of papers [1,10,11,12,21]. For instance, in [21], the authors considered a system with three compartments, that is, a homogeneous compartment of susceptible hosts, a heterogeneous compartment of infective agents (pathogen spores) and an age-structured compartment of infected hosts. A model with a homogeneous susceptible compartment and a continuously-structured infective compartment was also considered in [19]; the authors then reduced the model to a system of ordinary differential equations, one dependent variable of which represented the mean virulence.
To the best of our knowledge, there are no published papers that consider a model that includes a continuous phenotypic structure in the susceptible compartment and allows for changes in this phenotypic structure to occur, under the evolutionary pressure associated with a fitness trade-off between reproductive potential and susceptibility. This enables a mathematical dissection of the roles played by different evolutionary parameters in the phenotypic adaptation of the host population to the selective pressure exerted by the infection.
We frame our analysis in the context of an SI epidemic model, as it is common in evolutionary epidemiology [2,48]. The model comprises an ODE for the number of infected individuals which is coupled with a PIDE for the population density function of susceptible individuals through an integral term. Exploiting the analytical tractability of the model equations, we obtain a detailed mathematical depiction of the interactions between susceptible and infected individuals. In summary, the expression for the basic reproduction number R 0 is derived, the disease-free equilibrium and endemic equilibrium of the model are characterised and a threshold theorem involving R 0 is proved. Analytical results are integrated with the results of numerical simulations of a calibrated version of the model based on the results of artificial selection experiments presented in [59]. The analytical and numerical results obtained clarify the role of different evolutionary parameters on the spread of infectious diseases and the consequent phenotypic adaption of susceptible individuals.
The paper is organised as follows. In Section 2 we describe the mathematical model. In Section 3 we present the main analytical results of our study. In Section 4 we integrate analytical results with numerical simulations, and we discuss the biological relevance of our theoretical findings. Section 5 concludes the paper and provides a brief overview of possible research perspectives.
2 Description of the model Model equations. We consider an SI epidemic model whereby the susceptible compartment is structured by a continuous variable x ∈ X, with X being an interval (possibly infinite) of R, which models the phenotypic state of every susceptible individual and takes into account interindividual variability in proliferative potential (i.e. the number of progeny produced per unit time) and resistance to infection. The phenotype distribution of susceptible individuals at time t ≥ 0 is described by the population density function s(x, t) ≥ 0, while the number of infected individuals at time t is modelled by the function I(t) ≥ 0. The size of the susceptible compartment (i.e. the number of susceptible individuals) and the total number of individuals are defined, respectively, as Moreover, the mean phenotypic state of susceptible individuals and the related variance are defined, respectively, as The function σ 2 (t) provides a measure of the degree of phenotypic heterogeneity in the susceptible compartment. The evolution of the functions s(x, t) and I(t) is governed by the following PIDE-ODE system subject to a biologically relevant initial condition of components and to zero Neumann (i.e. no-flux) boundary conditions at the endpoints of X in the case where the interval X is bounded. In the PIDE (3) 1 , the diffusion term models the effect of heritable, spontaneous phenotypic changes, which occur at rate β > 0. The function a(x) is the net per capita growth rate of susceptible individuals in the phenotypic state x in the absence of infected individuals (i.e. the intrinsic net per capita growth rate), while the parameter K > 0 is related to the carrying capacity of the system (i.e. the maximal number of individuals that can be accommodated within the system due to space limitations). Furthermore, the function b(x) is the rate of infection (i.e. the product between the contact rate and the relative susceptibility of individuals in the phenotypic state x). Finally, the parameter ν > 0 in the ODE (3) 2 is the rate of death caused by the infection. Notice that we neglect physiological death of infected individuals since we assume its rate, which might depend on N (t), to be much lower than the rate of death caused by the infection. Notice also that infected individuals do not reproduce since we consider the case where proliferation of infected individuals is impaired by the infection.
Assumptions on the functions a(x) and b(x). We focus on a biological scenario whereby there is a fitness trade-off between resistance to infection and proliferative potential [18,35,49,59,60]. Under this scenario, without loss of generality, we let X be such that [0, 1] ⊂ X and make the following biological assumptions.
Assumption 1. Susceptible individuals in the phenotypic state x = 0 are characterised by the highest proliferative potential, whilst susceptible individuals in the phenotypic state x = 1 have the highest resistance to infection (cf. the schematics displayed in Figure 1).

Assumption 2.
Natural selection is moderately intense in the absence of infection. As a result, individuals in the phenotypic state x = 1 have a positive intrinsic net per capita growth rate. Moreover, due to natural selection, phenotypic variants that are not proliferative enough nor sufficiently resistant to infection (i.e. individuals in phenotypic states x sufficiently smaller than 0 or sufficiently larger than 1) will have a negative intrinsic net per capita growth rate.
Considering a biological scenario corresponding to Assumptions 1 and 2, we let the functions a(x) and b(x) be such that where R + * denotes the set of positive real numbers.
Specific definitions of the functions a(x) and b(x). Under assumptions (5) and (6), building upon the modelling strategies presented in [4,17,38], we will focus on the case where In definitions (7), the parameter γ > 0 models the intrinsic net per capita growth rate of susceptible individuals in the highly proliferative phenotypic state x = 0 (i.e. γ is the maximum intrinsic net per capita growth rate), and the parameter η > 0 is a selection gradient that quantifies the selective pressure acting on susceptible individuals prior to infection (i.e. η is the intrinsic selection gradient). The parameter ξ > 0 is the rate at which the infection is transmitted from infected individuals to susceptible individuals in the highly resistant phenotypic state x = 1 (i.e. ξ is the minimal infection rate), and the parameter ζ > 0 is a selection gradient that quantifies the intensity of the selective pressure that the infection exerts on susceptible individuals (i.e. ζ is the selection gradient related to infection).

Analytical results
In this section, we characterise the disease-free equilibrium of the model and derive the expression for the basic reproduction number R 0 (see Section 3.1). Furthermore, we characterise the endemic equilibrium of the model (see Section 3.2). Finally, we study the long-time asymptotic behaviour of the model equations by proving a threshold theorem involving R 0 (see Section 3.3).

Disease-free equilibrium and basic reproduction number
Disease-free equilibrium. The population density function of susceptible individuals at disease-free equilibrium s F (x) satisfies the following steady-state problem subject to zero Neumann boundary conditions if X is bounded.
In the case where the intrinsic net per capita growth rate a(x) satisfies assumptions (5), a characterisation of s F (x) is provided by Theorem 1, the proof of which relies on the facts recalled in Lemma 1 about the principal eigenpair (λ, ψ) of the elliptic differential operator (i.e. L[ψ] = −λψ) acting on functions defined on X. These are well known facts in the case where X is a bounded interval and Neumann boundary conditions are considered. The case where X ≡ R is considered, for instance, in [7, Theorem 2.3.1], while possible generalisations can be found in [5,6]. Notice that the characterisation of the principal eigenvalue given by (12) follows immediately, once the existence of a complete set of discrete eigenfunctions of L has been shown.
Lemma 1. When a(x) satisfies assumptions (5), the principal eigenvalue λ of the elliptic differential operator L is simple and there is a unique normalised positive eigenfunction ψ associated with λ. The eigenfunction ψ is smooth and the principal eigenvalue λ is given by Theorem 1. Let assumptions (5) hold and assume a(x) and β to be such that Then, there exists a unique solution of problem (10), subject to zero Neumann boundary conditions if X is bounded, which is given by Proof. The differential equation (10) 1 can be rewritten as is an eigenvalue of the elliptic differential operator L. Hence, Lemma 1 allows us to conclude that, under assumptions (5) and (13), there exists a unique solution of problem (10), subject to zero Neumann boundary conditions if X is bounded, which is given by Furthermore, integrating both sides of the differential equation (10) 1 over X, dividing through by S F and substituting the expression (15) for s F (x) into the resulting equation gives This concludes the proof of Theorem 1.
If X ≡ R and the intrinsic net per capita growth rate a(x) is defined according to (7), the population density function s F (x) can be explicitly characterised, as shown by Proposition 1. The proof of Proposition 1 is similar to that of Lemma 3.2 in [17] and, therefore, it is omitted here.
then the unique solution of problem (10) is with Basic reproduction number R 0 . In the modelling framework of (3), the basic reproduction number R 0 is defined as Under assumptions (5) and (13), substituting the expression (14) for the population density function of susceptible individuals at disease-free equilibrium into (19) gives In particular, when the intrinsic net per capita growth rate a(x) and the infection rate b(x) are defined via (7) and assumption (16) holds, substituting the expression (17) for the population density function s F (x) and definition (7) of the infection rate b(x) into (19), with a little algebra one finds the following explicit expression for the basic reproduction number

Endemic equilibrium
The population density function of susceptible individuals at endemic equilibrium s E (x) and the corresponding number of infected individuals I E satisfy the following steady-state problem subject to zero Neumann boundary conditions if X is bounded.
In the case where the intrinsic net per capita growth rate a(x) and the infection rate b(x) satisfy assumptions (5) and (6) and the additional assumption (13) holds, a characterisation of s E (x) and I E is provided by Theorem 2. The proof of Theorem 2 relies on the well known facts recalled in Lemma 2 (corresponding to Lemma 1) about the principal eigenpair (λ y , ψ y ) of the elliptic differential operator (i.e. L y [ψ y ] = −λ y ψ y ) with zero Neumann boundary conditions (if X is bounded) and acting on functions defined on X.
Lemma 2. When a(x) and b(x) satisfy assumptions (5) and (6), for any y ≥ 0 the principal eigenvalue λ y of the elliptic differential operator L y is simple, and there is a unique normalised positive eigenfunction ψ y associated with λ y . The eigenfunction ψ y is smooth, and the principal eigenvalue λ y is given by with Q(L y , φ) defined via (12), and is monotonically increasing in y.
Theorem 2. Let assumptions (5), (6) and (13) hold and assume R 0 > 1. Then there exists at least a positive solution (s E , I E ) of problem (22), subject to zero Neumann boundary conditions if X is bounded, whereby I E is a positive solution of the algebraic equation and then the solution I E of the algebraic equation (24) is unique and the solution of problem (22) is unique as well.
Proof. The differential equation (22) 1 can be rewritten as is an eigenvalue of the elliptic differential operator L I E .
Lemma 2 allows us to conclude that, under assumptions (5), (6) and (13), for there exists a unique solution of (22) 1 complemented with (22) 3 , subject to (22) 4 , and to zero Neumann boundary conditions if X is bounded, which is given by Since S E = −K λ I E , using the fact that λ I E is monotonically increasing in I E , as established by Lemma 2, we conclude that S E is monotonically decreasing in I E .
Substituting expression (25) for s E (x) into (22) 2 gives with B I E defined via (24). Since S E = −K λ I E , we can introduce the function and rewrite (28) as Since L 0 is the same as L defined in (11), using (15) we find that Hence, under the assumption R 0 > 1, we have that F (0) > ν. On the other hand, there existsȳ > 0 such that λȳ = 0 and, therefore, F (ȳ) = 0. By continuity of the function F (y), there exists (26) is met, the solution I E of the algebraic equation (24) is unique.
If X ≡ R and the intrinsic net per capita growth rate a(x) and the infection rate b(x) are defined according to (7), the population density function s E (x) and the size of the infected compartment I E can be characterised more explicitly, as shown by Proposition 2.
Proposition 2. Let X ≡ R. If the functions a(x) and b(x) are defined via (7), assumption (16) holds and R 0 > 1, then the unique solution of problem (22) is given by with and by I E that is the unique positive solution of the algebraic equation with F (y) := ν ξ + ζ ζy η + ζy − 1 2 + ζ β η + ζy Proof. Throughout the proof we drop the subscript E for convenience of notation. Building upon the method of proof presented in [17,52], first we note that when a(x) and b(x) are defined according to (7) the differential equation (22) 1 can be rewritten as where the functions g, h and ϕ are defined via (9). Then we make the change of variables y = x − ϕ(I) and rewrite the differential equation (34) as Finally, we make the additional change of variables and in so doing we find that u(z) satisfies the differential equation with P := 1 2 K β h(I) 1 2 S − K g(I) + I .
The differential equation (37) From the expression (38) for S we see that if assumption (16) holds then for I small enough we have S > 0. The bounded and non-negative solutions of the differential equation (37) are of the form Combining (39)  and conclude that If I and s(x) satisfy conditions (22) 4 then equation (22) Substituting (40) into (41) gives Equating expressions (38) and (42) for S and rearranging terms we find and with F (I) defined via (42). The left-hand side of the algebraic equation (43) is a straight line with positive slope and vertical intercept 0, whereas the righthand side is a monotonically decreasing function of I. Hence, there is a unique I such that the algebraic equation (43) is satisfied and for the condition I > 0 to be met it suffices that the right-hand side of equation (43) is positive for I −→ 0 + . This gives the following condition Since the functions a(x) and b(x) are defined according to (7) and assumption (16) holds, the basic reproduction number R 0 is given by (21) and, therefore, condition (44) is equivalent to the condition R 0 > 1.

Long-time asymptotic behaviour
A characterisation of the asymptotic behaviour for t → ∞ of the solution to the PIDE-ODE system (3) under assumptions (5) and (6) is provided by the results of Theorem 3.
Theorem 3. Under assumptions (5), (6) and (13), the solution to the PIDE-ODE system (3) subject to (4), and to zero Neumann boundary conditions if X is bounded, satisfies the following: with s F (x) given by (14), and and either X is bounded or Proof. Preliminary results used in the proof of (i). Under assumptions (5), (6) and (13), we have that s(·, t) ≥ 0 and I(t) ≥ 0 for all t ∈ [0, ∞), and s(x, t) is a sub-solution of the Cauchy problem subject to zero Neumann boundary conditions if X is bounded. It is known (see, for instance, Theorem 3 in [39] for the case where X is bounded) that with (λ, ψ) being the principal eigenpair of the elliptic differential operator L (cf. Lemma 1). The asymptotic result (52) along with the results established by Theorem 1 allows us to conclude that with s F (x) given by (14). The asymptotic result (53) ensures that with R 0 defined according to (19).
Proof of (i). Assume by contradiction that I(t) is bounded away from zero for t → ∞. If so, s(x, t) is a strict sub-solution of the Cauchy problem (51) and the asymptotic result (54) along with assumption (45) allows one to conclude that there exists t * < ∞ sufficiently large and ε > 0 sufficiently small so that This along with the ODE (3) 2 gives which contradicts the original assumption and allows us to conclude that, under assumption (45), the asymptotic result (47) holds. Therefore, the PIDE (3) 1 can be rewritten as Due to the structure of the ODE (3) 2 , the asymptotic result (47) implies that Σ(x, t) −→ 0 exponentially fast in L ∞ (X) as t → ∞. As a result, starting from (56) and using arguments analogous to those used to obtain the asymptotic result (53) one can prove that the asymptotic result (46) holds.
Preliminary results used in the proof of (ii). Integrating both sides of the PIDE (3) 1 over X and adding the resulting ODE for S(t) and the ODE (3) 2 for I(t) we obtain the following ODE for N (t) := S(t) + I(t) Estimating from above using the nonnegativity of S and I, the positivity of ν and the fact that max x∈R a(x) = a(0) > 0 (cf. assumptions (5)) we obtain the following differential inequality Since N (0) < ∞ (cf. initial conditions (4)), a(0) < ∞ and 0 ≤ S ≤ N is nonnegative, the latter differential inequality implies that Proof of (ii). Assume by contradiction that I(t) −→ 0 as t → ∞. If this occurs, then the PIDE (3) 1 can be rewritten in the form of (56) and the same arguments used in the proof of (i) yield The latter asymptotic result along with assumption (48) allows one to conclude that there exists t * < ∞ sufficiently large and ε > 0 sufficiently small so that This along with the ODE (3) 2 gives which contradicts the assumption that I(t) −→ 0 as t → ∞, as well as the upper bound (57). Since I(t) ≥ 0, we conclude lim sup t→∞ I(t) > 0. The structure of the and assumptions (5) Hence, If X is bounded then and the conclusion follows immediately from (58). (50) is what is called weak persistence in [51]. It seems likely that more careful estimates could lead to uniform weak persistence so that the use of appropriate theorems from Chapter 4 of [51] could lead to uniform strong persistence. This is certainly possible in the case when definitions (7) are considered, because in this case the PIDE (3) 1 can be reduced to a finite-dimensional system, but it is not pursued here.

Numerical results and biological interpretation
In this section, we report on the results of numerical simulations that complement the analytical results established in Section 3. In Section 4.1, model calibration, set-up of numerical simulations and numerical methods are described. In Section 4.2, the results of base-case numerical simulations carried out assuming either R 0 ≤ 1 or R 0 > 1 are presented. In Section 4.3, these results are integrated with the results of additional numerical simulations that summarise the outcomes of a systematic sensitivity analysis with respect to the model parameters, and the biological implications of our theoretical findings are discussed.

Model calibration, set-up of numerical simulations and numerical methods
Model calibration. We define the functions a(x) and b(x) via (7). The parameter values used to carry out base-case numerical simulations are listed in Table 1. These parameters are estimated using the results of artificial selection experiments in a host-parasite system presented by Webster & Woolhouse [59], who demonstrated the existence of a fitness trade-off between resistance to parasitic Schistosoma mansoni and fertility in Biomphalaria glabrata. In summary, the values of the parameters γ and η are chosen such that the intrinsic net per capita growth rates of individuals in the phenotypic states x = 0 and x = 1 (i.e. the values of a(0) = γ and a(1) = γ − η) match with the average experimental net proliferation rates of host individuals with high fertility and high resistance to infection, respectively. Notice that the values of γ and η so obtained are such that, coherently with assumptions (5), we have γ > η. Given the values of γ and η, we identified the values of the parameters ξ, ζ and K through exploratory numerical simulations carried out assuming β = 0 (i.e. neglecting phenotypic changes), letting the value of the parameter ν match the experimental value of the death rate of infected individuals, and choosing ξ, ζ and K in such a way as to minimise the mean square error between the average experimental value of the number of infected individuals at different time instants and the corresponding values of I(t). Finally, the value of the parameter β is chosen such that condition (16) is met and the value of ν is chosen on a case-by-case basis so that R 0 given by (21) satisfies either (45) or (48).

Remark 3.
Since the functions a(x) and b(x) are defined via (7), throughout the rest of the paper the basic reproduction number R 0 is computed via (21). Set-up of numerical simulations and numerical methods. We select a uniform discretisation consisting of 1200 points on the interval (−L, L) with L = 10 as the computational domain of the independent variable x and impose zero Neumann boundary conditions in x = −L and x = L. Moreover, we assume t ∈ (0, T ], with T > 0 being the final time of simulations, and we discretise the interval (0, T ] with the uniform step ∆t = 10 −4 . The method for constructing numerical solutions to the PIDE (3) 1 is based on an explicit finite difference scheme in which a three-point stencil is used to approximate the diffusion term and an explicit finite difference scheme is used for the reaction term [36]. Numerical solutions to the ODE (3) 2 are constructed using the explicit Euler method. All numerical computations are performed in Matlab.
We consider an initial condition that satisfies (4), that is, and I(0) = 10 with c 0 ∈ R + * being a normalisation constant such that

Base-case numerical simulations
Base-case numerical simulations for R 0 ≤ 1. The numerical solutions displayed in Figure 2 demonstrate that, in agreement with the asymptotic results of claim (i) in Theorem 3, when condition (45) is met (i.e. if R 0 ≤ 1) I(t) decays to zero while s(x, t) converges to an equilibrium population density function. Coherently with Theorem 1 and Proposition 1, the equilibrium population density function is the function s F (x) given by (17) and, therefore, S(t) converges to the equilibrium value S F given by (18).  (45) is met (i.e. when R 0 ≤ 1). The dashed line in the central panel highlights the equilibrium value S F given by (18). The population density function s(x, t) at the final time t = 300 is displayed in the left panel (solid line), where the dashed line highlights the equilibrium population density function s F (x) given by (17).
The values of the model parameters are those reported in Table 1 with ν = 2.
Base-case numerical simulations for R 0 > 1. The numerical solutions displayed in Figure 3 demonstrate that, in agreement with the asymptotic results of claim (ii) in Theorem 3, when condition (48) is met (i.e. if R 0 ≤ 1) both I(t) and S(t) remain bounded away from zero. It can also be seen that I(t) converges to a stationary value, which is the unique equilibrium I E , the existence of which is shown in Theorem 2, which is computed by solving numerically the algebraic equation (32). Moreover, s(x, t) converges to the equilibrium population density function s E (x) given by (29) and, therefore, S(t) converges to the equilibrium value S E given by (31).

Sensitivity analysis with respect to the model parameters and main biological implications
Conditions underpinning spread of infection. The results of claim (i) in Theorem 3, along with the numerical solutions displayed in Figure 2, demonstrate that if R 0 ≤ 1 then an infection cannot spread and a disease-free equilibrium is ultimately attained. On the other hand, the results of claim (ii) in Theorem 3, along with the numerical solutions displayed in Figure 3, show that if R 0 > 1 then the infection will spread.  (29). The values of the model parameters are those reported in Table 1 with ν = 0.2.
From (21) one sees that the basic reproduction number R 0 is a decreasing function of the rate of death caused by infection ν and of the intrinsic selection gradient η, and an increasing function of the minimum infection rate ξ and of the selection gradient related to infection ζ. Moreover, R 0 is an increasing function of β for η sufficiently low and a decreasing function of β for η sufficiently high. This is illustrated by the heat maps in Figure 4.
From the disease-free equilibrium given in Proposition 1, one sees that the asymptotic value of the size of the susceptible compartment S at the diseasefree equilibrium and the basic reproduction number R 0 are both decreasing in η. Moreover, from the endemic equilibrium given in Proposition 2, we observe that the asymptotic value of the mean phenotypic state µ, bounded from below by 0, is decreasing in η as well, meaning that the higher is η the smaller the susceptible population will be, and it will be concentrated on a phenotypic range which makes it more vulnerable to the infection.
These results provide a mathematical formalisation of the idea that infections exerting stronger selective pressures on susceptible individuals, and being characterised by lower mortality rates and higher infection rates, are more likely to spread. Furthermore, when susceptible individuals are exposed to weaker intrinsic selective pressures, frequent heritable, spontaneous phenotypic changes may promote the spread of an infection. On the other hand, higher rates of heritable, spontaneous phenotypic changes may correlate with a lower risk of spreading infection in the case where susceptible individuals are subject to stronger intrinsic selective pressures.  Table 1 with ν = 2 for the plots in the central and right panels, while β = 0.1 for the plot in the central panel.
Dependence of the disease-free equilibrium on evolutionary parameters. Coherently with Theorem 1 and Proposition 1, the numerical solutions presented in Figure 2 indicate that, if a disease-free equilibrium is attained, the size of the susceptible compartment S converges to the equilibrium value S F given by (18), which increases with the maximum intrinsic net per capita growth rate γ and the rescaled carrying capacity K. Moreover, as shown also by the numerical solutions displayed in Figure 5 and Figure 6, the equilibrium value S F decreases with the intrinsic selection gradient η and the rate of heritable, spontaneous phenotypic variation β.
The phenotypic distribution of susceptible individuals at disease-free equilibrium is of the Gaussian type (17). The mean value of the equilibrium phenotypic distribution corresponds to the fittest phenotypic state ϕ(0) defined via (9), which coincides with the phenotypic state characterised by the highest proliferative potential, i.e.
The variance of the equilibrium phenotypic distribution is governed by the ratio between β and the effective selection gradient h(0) defined via (9), which coincides with the intrinsic selection gradient η, i.e.
This is confirmed by the numerical solutions displayed in Figure 5 and Figure 6.  (17). The values of the other model parameters are those reported in Table 1 with ν = 2.
These results demonstrate that, when a disease-free equilibrium is attained, the susceptible population will be mainly composed of individuals with a high proliferative potential and the degree of phenotypic heterogeneity will be an increasing function of the rate of heritable, spontaneous phenotypic variation and a decreasing function of the intrinsic selection gradient.
Dependence of the endemic equilibrium on evolutionary parameters. Coherently with Theorem 2 and Proposition 2, the numerical solutions presented in Figure 3 indicate that, if an endemic equilibrium is attained, the number of infected individuals I converges to the unique positive solution I E of the algebraic equation (32), while the size of the susceptible compartment S converges to the equilibrium value S E given by (31). Since I E is a monotonically decreasing function of ν and S E is a monotonically decreasing function of I E , we have that S E is a monotonically increasing function of ν. This is confirmed by the numerical solutions displayed in Figure 7. Moreover, the numerical results summarised by Figure 8 show that, for the choice of parameter values considered here, both I E and S E decrease with ζ.
The phenotypic distribution of susceptible individuals at endemic equilib-  (17). The values of the other model parameters are those reported in Table 1 with ν = 2.
rium is of the Gaussian type (29). The mean value of the equilibrium phenotypic distribution corresponds to the fittest phenotypic state ϕ(I E ) defined via (9), i.e.
The variance of the equilibrium phenotypic distribution is given by the ratio between β and the effective selection gradient h(I E ) defined via (9), i.e.
Hence, larger values of ζ and smaller values of ν -since I E is a monotonically decreasing function of ν, as mentioned earlier -will lead the peak of the equilibrium phenotype distribution of susceptible individuals to move from the phenotypic state characterised by the highest proliferative potential (i.e. x = 0) toward phenotypic states closer to the phenotypic state corresponding to the highest level of resistance to infection (i.e. x = 1). Moreover, higher ζ and smaller ν will correlate with a narrower equilibrium phenotype distribution of susceptible individuals (i.e. lower degrees of phenotypic heterogeneity in the susceptible compartment). This is confirmed by the numerical solutions presented in Figure 7 and Figure 8.  (29). The values of the other model parameters are those reported in Table 1.
These results support the idea that, if an endemic equilibrium is established, infections characterised by lower rates of death and exerting stronger selective pressures on susceptible individuals will lead to a smaller and less phenotypically diverse susceptible compartment, which will be mainly composed of individuals with lower proliferative potential and higher level of resistance to infection.

Conclusions and research perspectives
The results of our mathematical study disentangle the impact of different evolutionary parameters on the spread of infectious diseases and the consequent phenotypic adaption of susceptible individuals. In particular, the results obtained provide a theoretical basis for the observation that infectious diseases exerting stronger selective pressures on susceptible individuals and being characterised by higher infection rates are more likely to spread. Furthermore, these results indicate that heritable, spontaneous phenotypic changes in proliferative , where the dashed lines highlight the equilibrium population density function s E (x) given by (29). The values of the other model parameters are those reported in Table 1 with ν = 0.2.
potential and resistance to infection can either promote or prevent the spread of infectious diseases depending on the strength of selection acting on susceptible individuals prior to infection. Finally, we have demonstrated that, when an endemic equilibrium is established, higher levels of resistance to infection and lower degrees of phenotypic heterogeneity among susceptible individuals are to be expected in the presence of infections which are characterised by lower rates of death and exert stronger selective pressures on susceptible individuals. Note that such results have been obtained without making any smallness assumption on the rate of phenotypic changes. This strengthens the robustness of the biological inferences drawn from the mathematical results obtained. As a complement to the analytical part of our work, we chose some more specific biologically relevant definitions of the model functions and, by using a suitable ansatz and carrying out direct computations, we provided a fullyexplicit characterisation of the disease-free equilibrium and endemic equilibrium of the model. We believe that the results obtained showcase the application potential of the analyses undertaken in the present study.
We conclude with an outlook on possible research perspectives. From a mathematical point of view, it would be interesting to complement the asymptotic results established by Theorem 3 by proving convergence to the endemic equilibrium when condition R 0 > 1 is satisfied. This will require to find a suitable Lyapunov function for the PIDE-ODE system (3). Moreover, it would be interesting to explore cases where infection transmission is frequency-dependent, instead of being density-dependent, and cases where a phenotypic structure is introduced in the infective compartment. It would also be interesting to consider the more general case of SIR models, which will make it necessary to introduce a phenotypic structure both in the infected and recovered compartments, in order to take into account the effect of the reproduction of recovered individuals on the evolutionary dynamics of the system. All of these extensions would require further development of the methods of proof presented here in order to carry out similar analyses of evolutionary dynamics.
Since the analytical results we have obtained can accommodate parameter values for a wide range of SI systems and a variety of infectious diseases, from an application point of view it would be interesting to apply these results to specific datasets in order to test the validity of the explicit condition for the spread of infection that we have derived. Moreover, building upon previous work on the derivation of deterministic continuum models for the evolution of populations structured by phenotypic traits from stochastic individual-based models [3,15,16,17,52], another track to follow to further enrich the present study would be to develop a stochastic individual-based model corresponding to the deterministic continuum model presented here. This would make it possible to explore the impact of stochastic fluctuations in phenotypic properties of single individuals on the spread of infections.
Finally, given the fact that spatial interactions are central to the spread of several infectious diseases, an additional natural way of extending our study would be to introduce spatial structure in order to consider scenarios in which individuals move across space and/or migrate between different regions that occupy the nodes of a network.