How heterogeneity in density dependence affects disease spread: when lifestyle matters

Abstract

People's lifestyles play a major role in disease risk. Some employment sectors and transport modes involve fixed exposures regardless of community size, while in other settings exposure tracks with population density. MERS-CoV, a coronavirus discovered in Saudi Arabia in 2012 closely related to those causing SARS and COVID-19, appears to need extended contact time for transmission, making some segments of a community at greater risk than others. We model mathematically how heterogeneity in contact rate structure impacts disease spread, using as a case study a MERS outbreak in two Saudi Arabian communities. We divide the at-risk population into segments with exposure rates either independent of population density or density-dependent. Analysis shows disease spread is minimized for intermediate size populations with a limited proportion of individuals in the density-independent group. In the case study, the high proportion of density-independent exposure may explain the historical outbreak's extinction in the larger city.

Keywords:

• Middle East respiratory syndrome
• contact rate
• population density
• epidemic outbreak

1. Introduction

In the twenty-first century, new outbreaks of coronaviruses have been a major source of concern, beginning with Severe Acute Respiratory Syndrome (SARS) in 2003, followed by Middle East Respiratory Syndrome (MERS) in 2012, and most recently Coronavirus disease (COVID-19) in 2019. MERS is a viral respiratory illness caused by the MERS-CoV virus. First reported in Saudi Arabia in 2012, it has since spread to several other countries in the Middle East and Asia. MERS can cause severe respiratory illness, including fever, cough, and shortness of breath. In some cases, it can also lead to pneumonia and kidney failure. MERS is primarily spread through extended close contact with infected individuals, such as through respiratory droplets or close contact with surfaces contaminated with the virus. Notably, such extended contacts have occurred in spillover outbreak events from dromedary camels, and in nosocomial settings. It is considered a severe public health threat, and efforts are ongoing to better understand MERS's transmission dynamics and develop effective treatments and preventative measures [5, 15].

Growth in global infrastructure and international travel has amplified the transmission risk of viruses such as SARS-CoV-1, SARS-CoV-2, MERS-CoV, and H1N1 influenza (swine flu) [1, 8]. Preventative measures, such as restrictions on international travel, can help contain the spread of such diseases and prevent a global epidemic. However, even an isolated outbreak in one country can have widespread negative impacts on the rest of the world [13], making local contact structures equally important. MERS shares many characteristics with SARS and COVID-19, including significant mortality rates. However, its community transmission has been limited due to the requirement for extended periods of close contact for the virus to be transmitted [3]. Such extended contacts may occur in greatest number in work or public settings (household exposure is important but limited to household size), but even within a single community there may be significant heterogeneity in the amounts of time that individuals spend in such settings. The COVID-19 pandemic caused major changes in people's everyday contact patterns, some of which have become permanent as certain employment sectors allow working remotely. This observation inspires the study of a population in which some individuals (such as those who work in large groups or take long trips via public transportation) have density-dependent contacts, while others (such as those who work from home and travel alone) have density-independent contacts, and yet both groups still have contact with each other. This study examines how population density impacts disease transmission in such heterogeneous populations to gain insight into ways to minimize disease spread.

Two classical models have been used to describe disease transmission dynamics: mass action incidence and standard incidence. The assumption underlying mass action (density-dependent) incidence is that the [per capita] number of daily encounters adequate for transmission (here, involving extended exposure time) is proportional to population size. This model is typically suitable for small populations below a threshold value at which infectious contacts saturate. Standard (density-independent) incidence, on the other hand, depends instead on the proportion (rather than absolute number) of infected hosts in the population, under an assumption that the overall contact rate is already saturated. The infectivity threshold depends on the average number of contacts per infected individual, as measured by the basic reproduction number. The disease persists if infectives spend enough time in contact with susceptible individuals to cause more than one new infection per infective [7]. The difference between mass action and standard incidence is observed primarily under significant changes or differences in population density.

Although MERS outbreaks have been amplified most significantly in healthcare settings, where extended contacts are facilitated, in this study we consider MERS as an example of an infectious disease where populations are highly heterogeneous as to the rates at which individuals make contacts adequate for transmission. Several modelling studies have considered aspects of MERS transmission, but we found none which considered heterogeneity in density dependence of contacts. However, Sardar et al. compared three different incidence functions for fit to weekly incidence data of MERS-CoV from outbreaks in three major provinces of Saudi Arabia: Riyadh, Makkah, and Madinah [11]. They considered infection rates described by mass action (density dependence), Holling type II saturation in the number of infectives (per capita contacts $\text{β}\left(I\right)={\text{β}}_0/(1+\alpha I)$), and a kind of hypersaturation ($\text{β}\left(I\right)={\text{β}}_0/(1+\alpha I^2)$) which leads to a total infection rate that peaks at a limited outbreak size. The study used three measures: the timing and amplitude of the outbreak peak, and the total number of cases in a year. Of the three forms, the study found that simple density dependence best explained the observed data except in Riyadh, where two distinct strains of MERS were cocirculating and (perhaps due to interactions between the two strains) the hypersaturation incidence function best captured the multiple peaks observed.

Building upon a preliminary investigation [12], this study aims to identify how heterogeneity in the infectious contact structure affects overall community transmission in the population. In the context of MERS, how does heterogeneity in extended exposure through the work and transport sectors affect the spread of MERS in a community? As a first approximation, we consider two groups within the population: those who work on-site in large groups and/or use public transportation to and from work, thus spending significant periods of time in density-dependent settings; and those who work from home (or with limited exposure to others) and use their own vehicles to travel, thus primarily in density-independent settings. Mass action and standard incidence, respectively, capture transmission dynamics within each type of setting; which type is more conducive to the spread of infection depends on the size of each group, which may be affected by economic or political factors. The prolonged contacts required to transmit MERS may increase when more of the population experiences a primarily density-dependent exposure, unlike in settings with density-independent exposure.

This study uses a compartmental dynamical systems model, with qualitative analysis to identify overall trends and numerical analysis to highlight specific cases and features. To illustrate the impact of population density and contact structure heterogeneity on disease transmission in a hypothetical MERS outbreak, we consider two cities with different population densities: Buraydah and Ash Shimasiyah, in Saudi Arabia, where the disease was first diagnosed. Buraydah was selected due to its high population density and a 2016 MERS outbreak reported there. Ash Shimasiyah is a smaller community approximately 37 km from Buraydah.

2. Model development

In the classic SIRS epidemiological model, the total population is classified by infection status as susceptible (S), infected (I), or recovered (R) individuals. The model is cyclical because susceptible individuals who become infected can then recover from the disease, and eventually return to susceptible status again after a period of temporary immunity. This paper modifies the simple SIRS model to gain deeper insight into how heterogeneity in infectious contact structure impacts disease dynamics while retaining the fundamental principles of the model. The main modification is to split the population into two groups (see Figure 1): one whose contact rate is not affected by population density (Group 1, $N_1\left(t\right)=S_1\left(t\right)+I_1\left(t\right)+R_1\left(t\right)$) and one whose contacts with others are influenced by population density (Group 2, $N_2\left(t\right)=S_2\left(t\right)+I_2\left(t\right)+R_2\left(t\right)$). Because the dynamics of the two groups are indistinguishable except for this dependence in their infectious contact rates, we can attribute any variations between those dynamics more precisely to this population density factor.

Figure 1. Flow diagram reflecting infection dynamics for groups 1 and 2, interconnected by the infection rates via the weighted average proportion of contacts made with infectives, $q=\frac{{\text{β}}_1{I}_1+{\text{β}}_2{I}_{2}N}{{\text{β}}_1{N}_1+{\text{β}}_2{N}_{2}N}$.

This diagram shows new arrivals entering either the class S1 or the class S2. There are flows from S1 to I1, from I1 to R1, from R1 back to S1, and from all three of these classes leaving the study population altogether. Similar flows exist among S2, I2, and R2.

We therefore begin by describing the contact structure which drives the transmission dynamics for both groups. Individuals in group 1 are assumed to have the kinds of extended contacts with others necessary to transmit a virus like MERS-CoV, at a per capita rate independent of population size. These individuals spend time in primarily density-independent settings such as small office facilities and private vehicles (whether their own or rideshare) which offer a fixed number of exposures per day. This leads to a constant per capita contact rate, say ${\text{β}}_1$, as in standard incidence. On the other hand, individuals in group 2 spend time primarily in density-dependent settings like large work sites (such as factories) and public transportation, where the number of exposures per day is affected by population size. This leads to a contact rate of the form ${\text{β}}_{2}N$ consistent with mass-action incidence. We assume that individuals' group membership is determined by their profession and lifestyle, which are not likely to change over the course of an outbreak, and therefore remains fixed throughout the epidemic regardless of infection status. Note, however, that cross-group contacts occur regularly: in public locations such as markets, stores, and places of worship, and in private locations such as work sites (whether corporate or domestic). Our simple model explores only the effects of heterogeneity in contact rate, and does not attempt to capture the contact network in detail. Note also that in order to manage model complexity and focus on the role lifestyle plays in community transmission, we do not also consider healthcare settings.

The infection rate for each group is given by the product of the per capita contact rate, the size of the corresponding susceptible class, and the proportion of potentially infectious contacts which are made with infectious individuals. In a simple model with a homogeneous contact structure, this product has the form $\text{β}\left(N\right)\cdot S\cdot \frac{I}{N}$. In our model $\text{β}\left(N\right)$ varies by group, and the proportion is a weighted average of the contacts made by infectives in both groups: (1) $\frac{{\text{β}}_1{N}_1}{{\text{β}}_1{N}_1+{\text{β}}_2{N}_{2}N}\left(\frac{I_1}{N_1}\right)+\frac{{\text{β}}_2{N}_{2}N}{{\text{β}}_1{N}_1+{\text{β}}_2{N}_{2}N}\left(\frac{I_2}{N_2}\right)=\frac{{\text{β}}_1{I}_1+{\text{β}}_2{I}_{2}N}{{\text{β}}_1{N}_1+{\text{β}}_2{N}_{2}N},$(1) The proportions contributed by each group come from the total number of contacts made by group 1 members in unit time, ${\text{β}}_1{N}_1$, and the total for group 2, $\left({\text{β}}_{2}N\right)N_2$. The total infection rate in group 1 is then given by (2) ${\text{β}}_1\left(\frac{{\text{β}}_1{I}_1+{\text{β}}_2{I}_{2}N}{{\text{β}}_1{N}_1+{\text{β}}_2{N}_{2}N}\right)S_1,$(2) while the total group 2 infection rate is given by (3) ${\text{β}}_{2}N\left(\frac{{\text{β}}_1{I}_1+{\text{β}}_2{I}_{2}N}{{\text{β}}_1{N}_1+{\text{β}}_2{N}_{2}N}\right)S_2.$(3) Thus individuals in each group are exposed to infectious individuals in both groups, but in different proportions, where the proportion of exposure to group 2 infectives increases with population size $N\left(t\right)$ (as does group 2's total exposure to each group). In this way the model depicts a contact rate which is a hybrid of standard and mass-action incidence.

Some epidemiological models do not account for births and deaths in a population undergoing an outbreak. However, demographic renewal rates can be significant on the timescale of an epidemic, especially disease-related deaths, as seen in the recent COVID-19 pandemic, and are thus appropriate to include. We use a modified logistic growth term to account for population growth due to births, based on the total population size $N=N_1+N_2$, the region's carrying capacity K, the maximum per capita reproductive rate r, and the proportion p (assumed fixed) of births into group 1. Correspondingly, a proportion $(1-p)$ of the total births and new members is assigned to group 2. We assume that infection and immunity cannot be passed from mother to offspring, so all new arrivals enter the susceptible class of the given group.

The remaining transitions between classes are assumed to be linear (i.e. driven by constant per capita rates) and independent of group membership (1 or 2). The natural death rate d is also independent of infection status. Infected individuals recover (or enter treatment and isolation) at per capita rate γ and experience additional mortality due to the disease at rate δ, making their total per capita death rate $d+\delta$. Recovered individuals lose their immunity at rate θ, returning to the susceptible classes.

The model is given by the following dynamical system and by Figure 1; state variables and parameters are summarized in Table 1. Disease-induced death requires us to consider N, $N_1$, and $N_2$ as functions of time. (4) $\begin{array}{rl}\frac{d{S}_1}{dt}& =prN(1-\frac{N}{K})+\theta R_1-{\text{β}}_1\left(\frac{{\text{β}}_1{I}_1+{\text{β}}_2{I}_{2}N}{{\text{β}}_1{N}_1+{\text{β}}_2{N}_{2}N}\right)S_1-d{S}_1,\\ \frac{d{S}_2}{dt}& =(1-p)rN(1-\frac{N}{K})+\theta R_2-{\text{β}}_{2}N\left(\frac{{\text{β}}_1{I}_1+{\text{β}}_2{I}_{2}N}{{\text{β}}_1{N}_1+{\text{β}}_2{N}_{2}N}\right)S_2-d{S}_2,\\ \frac{d{I}_1}{dt}& ={\text{β}}_1\left(\frac{{\text{β}}_1{I}_1+{\text{β}}_2{I}_{2}N}{{\text{β}}_1{N}_1+{\text{β}}_2{N}_{2}N}\right)S_1-(\gamma +d+\delta )I_1,\\ \frac{d{I}_2}{dt}& ={\text{β}}_{2}N\left(\frac{{\text{β}}_1{I}_1+{\text{β}}_2{I}_{2}N}{{\text{β}}_1{N}_1+{\text{β}}_2{N}_{2}N}\right)S_2-(\gamma +d+\delta )I_2,\\ \frac{d{R}_1}{dt}& =\gamma I_1-(\theta +d)R_1,\\ \frac{d{R}_2}{dt}& =\gamma I_2-(\theta +d)R_2.\end{array}$(4)  

Table 1. State variable and parameter definitions and their units.

	Notation	Definition
State variables	$S_i\left(t\right)$	Number of susceptible individuals in group i at time t
	$I_i\left(t\right)$	Number of infected individuals in group i at time t
	$R_i\left(t\right)$	Number of recovered individuals in group i at time t
Demographic parameters	p	Fraction of new arrivals assigned to group 1 (dimensionless)
	$(1-p)$	Fraction of new arrivals assigned to group 2 (dimensionless)
	r	Maximum growth rate ($1/\mathrm{t}\mathrm{i}\mathrm{m}\mathrm{e}$)
	K	Carrying capacity ($\mathrm{\#}\mathrm{p}\mathrm{e}\mathrm{o}\mathrm{p}\mathrm{l}\mathrm{e}/{\mathrm{k}\mathrm{m}}^2$)
	d	Natural death rate ($1/\mathrm{t}\mathrm{i}\mathrm{m}\mathrm{e}$)
Epidemiological parameters	${\text{β}}_1$	per capita infectious contact rate for group 1 ($1/\mathrm{t}\mathrm{i}\mathrm{m}\mathrm{e}$)
	${\text{β}}_{2}N\left(t\right)$	per capita infectious contact rate for group 2 ($1/\mathrm{t}\mathrm{i}\mathrm{m}\mathrm{e}$)
	δ	Infection induced death rate ($1/\mathrm{t}\mathrm{i}\mathrm{m}\mathrm{e}$)
	γ	Recovery rate ($1/\mathrm{t}\mathrm{i}\mathrm{m}\mathrm{e}$)
	θ	Immunity waning rate for recovered individuals ($1/\mathrm{t}\mathrm{i}\mathrm{m}\mathrm{e}$)

3. Qualitative analysis

The state of the population in the absence of disease is given by the disease-free equilibrium (DFE). Here the infected and recovered classes are empty, leaving groups 1 and 2 entirely in the susceptible classes, making the total population $N=S_1+S_2$ follow a logistic law, $\frac{dN}{dt}=rN(1-\frac{N}{K})-dN.$Thus as long as $N\left(0\right)>0$, the total population approaches $N^{\ast }=K(1-\frac{d}{r})$, with $S_1^{\ast }=p{N}^{\ast }$, $S_2^{\ast }=(1-p)N^{\ast }$.

The canonical measure of an infection's ability to invade a population is given by its basic reproduction number ${\mathcal{R}}_0$, defined in terms of the DFE: it gives the average number of new cases produced by one infected individual in a fully susceptible population. The next generation operator approach developed by Diekmann et al. [7] calculates ${\mathcal{R}}_0$ as the dominant eigenvalue of the next generation matrix, which is drawn from the rates of change of the infected classes. For the model of this study, the method yields ${\mathcal{R}}_0=\frac{{\text{β}}_1\left({\text{β}}_{1}p\right)+{\text{β}}_2{N}^{\ast }\left({\text{β}}_2{N}^{\ast }\right(1-p\left)\right)}{{\text{β}}_{1}p+{\text{β}}_2{N}^{\ast }(1-p)}\cdot \left(\frac{1}{\delta +\gamma +\theta }\right).$This expression gives the average infection rate weighted by each group's number of contacts, multiplied by the average duration of infectivity. Averages are taken across both groups. The expression can also be rewritten to show each group's contribution to disease spread, of which ${\mathcal{R}}_0$ again becomes a weighted average: ${\mathcal{R}}_0=\left(\frac{{\text{β}}_1}{\delta +\gamma +\theta }\right)\left[\frac{{\text{β}}_{1}p}{{\text{β}}_{1}p+{\text{β}}_2{N}^{\ast }(1-p)}\right]+\left(\frac{{\text{β}}_2{N}^{\ast }}{\delta +\gamma +\theta }\right)\left[\frac{{\text{β}}_2{N}^{\ast }(1-p)}{{\text{β}}_{1}p+{\text{β}}_2{N}^{\ast }(1-p)}\right].$Since the expressions in square brackets are just contact-based weighting terms for each group, the expressions $\left(\frac{{\text{β}}_1}{\delta +\gamma +\theta }\right)$ and $\left(\frac{{\text{β}}_2{N}^{\ast }}{\delta +\gamma +\theta }\right)$ can be seen to symbolize the contributions ${\mathcal{R}}_{01}$ and ${\mathcal{R}}_{02}$, respectively, of each group.

The next generation approach ties the local asymptotic stability of the DFE to the condition ${\mathcal{R}}_0<1$. Analytical identification of endemic equilibria is more complex; however, graphical and numerical evaluations provide evidence for the existence of a unique endemic equilibrium when ${\mathcal{R}}_0>1$.

4. Numerical analysis

Numerical analysis of a hypothetical case study illustrates the impact of heterogeneity in contact structure and of population size on MERS transmission. Estimates for model parameter values were derived from existing demographic and epidemiological data and are summarized in Table 2.

Table 2. Summary of estimated model parameters.

Parameter (unit)	Value
d (${\mathrm{d}\mathrm{a}\mathrm{y}\mathrm{s}}^{-1}$)	7.283$\times {10}^{-6}$
r (${\mathrm{d}\mathrm{a}\mathrm{y}\mathrm{s}}^{-1}$)	3.128$\times {10}^{-4}$
K (people)	12,600; 794,600
δ (${\mathrm{d}\mathrm{a}\mathrm{y}\mathrm{s}}^{-1}$)	0.0833
γ (${\mathrm{d}\mathrm{a}\mathrm{y}\mathrm{s}}^{-1}$)	0.1548
θ (${\mathrm{d}\mathrm{a}\mathrm{y}\mathrm{s}}^{-1}$)	0.000913
${\text{β}}_1$ (${\mathrm{d}\mathrm{a}\mathrm{y}\mathrm{s}}^{-1}$)	0.05
${\text{β}}_2$ (${\mathrm{p}\mathrm{e}\mathrm{o}\mathrm{p}\mathrm{l}\mathrm{e}}^{-1}{\mathrm{d}\mathrm{a}\mathrm{y}\mathrm{s}}^{-1}$)	$5\times {10}^{-7}$
p (dimensionless)	0–1

Demographic data for the country of Saudi Arabia were taken from the World Bank [14]. The per capita death rate d was approximately 0.00266/yr immediately before the COVID-19 pandemic, at which point the per capita birth rate (given in our model as $r(1-N/K)$) was 0.01845/yr. We fit a logistic model to the birth rates and population sizes for 2016–2020 to calculate the best-fit value for the carrying capacity K (about 42.7M) and used the proportion of the most recent population estimate for Saudi Arabia (35.95M, 2021) to K, to obtain a scale factor (0.8414) which we then applied to the populations of Buraydah and Ash Shimasiyah [10], to estimate the values of K for those communities. We also used the national value of K to back-calculate r from the value and expression above, as 0.1142/yr.

The MERS case fatality ratio (CFR) is given by the CDC [5] and other organizations as about 35% and occurring approximately 12 days from onset of illness. Infection may last only a few days for mild cases and a month or more for severe hospitalized cases. A recovery rate $\gamma =0.1548/$day representing a mean duration of about 6.5 days together with a disease-induced death rate $\delta =1/$(12 days) replicates the 35% CFR. The rate of loss of immunity (θ) was set to the inverse of 3 years, reflecting the persistence of MERS-specific antibodies and immune response factors [9].

The infection rate for MERS has not been measured in literature; however, given the extended nature of the contact assumed necessary for transmission between two humans, it should be much lower than that for other coronaviruses. Chowell et al. estimated an infection rate of 0.25/day for SARS in Hong Kong based on data fitting [6]. In considering the effect of population density through the examples of Buraydah and Ash Shimasiyah, we consider a value for group 1's ${\text{β}}_1$ of one-fifth that figure, or 0.05/day, and set the total per capita rate for group 2, ${\text{β}}_{2}N$, equal to that when $N={10}^5$, making ${\text{β}}_2=5\times {10}^{-7}$/day. (We will consider other values for specific scenarios below.)

The value of p, the proportion of the incoming population assigned to Group 1, whose daily contact patterns lead to a density-independent extended contact rate, depends on the socioeconomic structure of the community, and the proportion of individuals in different sectors of employment. In order to examine the effect of heterogeneity in contact structure, in our analysis we vary p between 0 and 1, while noting for reference that in the Al-Qassim region where the cities are located, the proportion of the workforce using private cars to travel to work and the proportion working in government jobs are both around 80% [10].

Figure 2 shows ${\mathcal{R}}_0$ as a function of the proportion of individuals in group 1 (vs. 2) for Buraydah and Ash Shimasiyah. It is evident that population size plays the biggest role when most of the population lives in a density-dependent environment (group 2). (Note that for a purely density-independent population, p = 1, ${\mathcal{R}}_0$ is independent of the population size.) For small populations like Ash Shimasiyah, the basic reproduction number increases with p because the disease spreads better in group 1 when population density is low. For larger populations like Buraydah, however, the high population density makes the disease spread better in group 2. In Buraydah MERS can spread except when most people live a density-independent lifestyle, $p\ge 80\mathrm{\%}$, which is approximately the proportion suggested by the employment and transport data cited above. This is consistent with the limited outbreak seen there in 2016.

Figure 2. The basic reproduction number ${\mathcal{R}}_0$ as a function of p for high (Buraydah, solid curve) and low (Ash Shimasiyah, dashed curve) values of $N^{\ast }$.

Motivated by the distinct trends observed in these two case studies, we consider two additional scenarios. In one, we adjust our infection rate estimates slightly (${\text{β}}_1=0.15$/day, ${\text{β}}_2=3\times {10}^{-6}$/day, making ${\text{β}}_1={\text{β}}_{2}N$ when $N=\mathrm{50,000}$) so that ${\mathcal{R}}_{01}$ is slightly less than 1. In this scenario, ${\mathcal{R}}_0>1$ only for sufficiently large population sizes; otherwise an outbreak dies out regardless of p (Figure 3). As seen in Figure 4, however, p still plays a role in ${\mathcal{R}}_0$ – in fact, for each value of p there is a different population density $N^{\ast }$ which minimizes disease spread, graphed in Figure 5. These critical values are small when nearly all of the population lives a density-dependent lifestyle, but when most of the population lives in density-independent settings, ${\mathcal{R}}_0$ is minimized at population sizes larger than some of the municipalities in the Al-Qassim region, including Ash Shimasiyah. This optimal relationship between population density and the level of heterogeneity in the contact structure suggests that disease spread is minimized when the proportion of the population living a density-independent lifestyle grows only once the population exceeds a given size.

Figure 3. Variation in ${\mathcal{R}}_0$ as a function of p and $N^{\ast }$, using ${\text{β}}_1=0.15$/day, ${\text{β}}_2=3\times {10}^{-6}$/day.

Figure 4. Effects of population density $N^{\ast }$ on ${\mathcal{R}}_0$ for different values of p. ${\text{β}}_1,{\text{β}}_2$ as in Figure 3.

Figure 5. The population size which minimizes ${\mathcal{R}}_0$, as a function of p. ${\text{β}}_1,{\text{β}}_2$ as in Figure 3.

In the other scenario that we consider here, ${\mathcal{R}}_{01}$ is slightly above 1, achieved by making ${\text{β}}_1=0.25$/day (the same as for SARS) and ${\text{β}}_2=5\times {10}^{-6}$/day (again ${\text{β}}_1={\text{β}}_{2}N$ at $N=\mathrm{50,000}$). In this case, the only way to prevent the spread of an outbreak (${\mathcal{R}}_0<1$) is to have a relatively small population with a predominantly density-dependent lifestyle, as seen in Figures 6 and 7. However, these figures also show the same phenomenon for very small populations observed in the first scenario: as $N\to 0$, ${\mathcal{R}}_0\to {\mathcal{R}}_{01}$ since nearly all contacts involve group 1. Thus for very small populations as well as large populations, ${\mathcal{R}}_0>1$ regardless of contact structure heterogeneity (p). Here again an intermediate population size minimizes disease spread and, in this case, limits the size of an outbreak.

Figure 6. Variation in ${\mathcal{R}}_0$ as a function of p and $N^{\ast }$, using ${\text{β}}_1=0.25$/day, ${\text{β}}_2=5\times {10}^{-6}$/day.

Figure 7. The dark shaded region shows combinations of p and $N^{\ast }$ for which ${\mathcal{R}}_0<1$, using ${\text{β}}_1=0.25$/day, ${\text{β}}_2=5\times {10}^{-6}$/day.

5. Conclusion

This study uses the basic reproduction number as a lens through which to examine the interplay between population density and contact structure heterogeneity on the one hand, and transmission of viruses like MERS-CoV which require extended contact on the other hand. In a simple model where each individual follows either a density-independent (fixed exposure through limited contacts at work and in transportation) or density-dependent lifestyle, disease spread is minimized for small to intermediate size populations with only a limited proportion of individuals in the density-independent group (p). The optimal population size for a given heterogeneity level increases (from 0) with p, rising quickly at first and then slowing down past a certain level. The non-monotone dependence of disease spread on population size is caused by the hybrid nature of the infection rate, where the portion represented by mass-action incidence increases with population size while the portion represented by standard incidence causes an asymptotic rise for very small populations.

A rough estimate of MERS community transmission in two neighbouring municipalities in Saudi Arabia which have seen a MERS outbreak indicates that any outbreaks in the smaller community are self-limiting regardless of contact structure, whereas in the larger city outbreaks may be stymied only because of the high proportion of individuals with a density-independent lifestyle. This case study highlights the impact of socioeconomic factors in shaping disease transmission in ways that can amplify or limit an outbreak.

The heterogeneity in contact structure depicted here (using two extremes) is just a sketch and could be fleshed out using data for specific cities on employment levels in different sectors, proportions of the day spent in extended contact with members of each group, and even a full contact network if desired. For MERS in particular, outbreaks are believed to have originated through human-to-camel contact, so data on where in the contact structure those kinds of contact occur would clarify how the outbreak spreads. Likewise the addition of a healthcare setting could distinguish the role of nosocomial transmission. A key quantity for comparing disease spread in different sectors, which would also require detailed data to determine, is the population size at which the two groups have equal infectious contact rates. In addition, this study has focused on the infection's baseline ability to spread as measured by ${\mathcal{R}}_0$; time-series data tracking the disease's spread in each segment of the population separately offers an additional measure of the effectiveness of preventive measures within each group. To the extent that lifestyle affects birth and death rates, those could also be differentiated.

This structure may also be applicable to the study of other infectious diseases known to require extended contacts for transmission such as tuberculosis, for which it has already been shown that public transportation networks can be an important source of infection for workers who travel long distances every day [4]. Alternatively, such a model could be used to study the effects on disease transmission of temporary but significant events like mass gatherings, such as the Hajj pilgrimage, where extremely high population densities can drive outbreaks in one group which then spread to another part of the population [2].

Acknowledgement

The authors acknowledge the preliminary work done by the authors of [12], whose research was supported by an NSF UBM-Institutional grant, DUE-0827136, as part of the UTTER program at UT Arlington.

Disclosure statement

No potential conflict of interest was reported by the author(s).