# Modelling Vector Transmission and Epidemiology of Co-Infecting Plant Viruses

Linda J S Allen, Vrushali A Bokil, Nik J Cunniffe, Frédéric Marie Hamelin, Frank M Hilker, Michael Jeger

► To cite this version:

Linda J S Allen, Vrushali A Bokil, Nik J Cunniffe, Frédéric Marie Hamelin, Frank M Hilker, et al. Modelling Vector Transmission and Epidemiology of Co-Infecting Plant Viruses. Viruses, 2019, 11 (12), pp.1153. 10.3390/v11121153. hal-02735995

HAL Id: hal-02735995

https://institut-agro-rennes-angers.hal.science/hal-02735995

Submitted on 2 Jun 2020

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

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

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

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_2.png" /></figure>

Article

# Modelling Vector Transmission and Epidemiology of Co-Infecting Plant Viruses

Linda J. S. Allen <fnref n="1" />*, Vrushali A. Bokil <fnref n="2" />, Nik J. Cunniffe <fnref n="3" />, Frédéric M. Hamelin <fnref n="4" />✉, Frank M. Hilker <fnref n="5" /> and Michael J. Jeger <fnref n="6" />

Received: 5 November 2019; Accepted: 6 December 2019; Published: 13 December 2019

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

check for
updates

**Abstract:** Co-infection of plant hosts by two or more viruses is common in agricultural crops and natural plant communities. A variety of models have been used to investigate the dynamics of co-infection which track only the disease status of infected and co-infected plants, and which do not explicitly track the density of inoculative vectors. Much less attention has been paid to the role of vector transmission in co-infection, that is, acquisition and inoculation and their synergistic and antagonistic interactions. In this investigation, a general epidemiological model is formulated for one vector species and one plant species with potential co-infection in the host plant by two viruses. The basic reproduction number provides conditions for successful invasion of a single virus. We derive a new invasion threshold which provides conditions for successful invasion of a second virus. These two thresholds highlight some key epidemiological parameters important in vector transmission. To illustrate the flexibility of our model, we examine numerically two special cases of viral invasion. In the first case, one virus species depends on an autonomous virus for its successful transmission and in the second case, both viruses are unable to invade alone but can co-infect the host plant when prevalence is high.

**Keywords:** co-infection; invasion reproduction number; vector transmission

## 1. Introduction

Transmission is a key element in understanding the epidemiology of plant virus diseases, particularly those transmitted by arthropod vectors [1–4]. In general, four modes of transmission, non-persistent, semi-persistent, persistent-circulative and persistent-propagative, can be distinguished. Each of these modes has a characteristic time period for acquisition from infected plants, retention in the vector, and inoculation to healthy plants [5], although some virus groups such as the torradoviruses do not fit neatly into these categories [6]. Other aspects important for arthropod transmission include transovarial and transtadial transmission, and the “helper strategy” [3] in which a helper virus can be transmitted by the vector but the dependent virus can only be transmitted in the presence of the helper, a strategy modelled by Zhang et al. [7].

Co-infection of hosts by two or more plant viruses is common in both agricultural crops [8,9] and natural plant communities [10,11]. Because of this, the literature on co-infection by plant viruses

<footnote n="1">¹ Department of Mathematics and Statistics, Texas Tech University, Lubbock, TX 79409, USA</footnote>

<footnote n="2">² Department of Mathematics, Oregon State University, Corvallis, OR 97331, USA;
bokilv@math.oregonstate.edu</footnote>

<footnote n="3">³ Department of Plant Sciences, University of Cambridge, Cambridge CB2 3EA, UK; njc1001@cam.ac.uk</footnote>

<footnote n="4">⁴ IGEPP, Agrocampus Ouest, INRA, Université de Rennes 1, Université Bretagne-Loire, 35000 Rennes, France;
frédéric.hamelin@agrocampus-ouest.fr</footnote>

<footnote n="5">5 Institute of Environmental Systems Research, School of Mathematics/Computer Science, Osnabrück
University, 49069 Osnabrück, Germany; frank.hilker@uni-osnabrueck.de</footnote>

<footnote n="6">⁶ Centre for Environmental Policy, Imperial College London, Ascot SL5 7PY, UK; m.jeger@imperial.ac.uk</footnote>

<footnote>* Correspondence: linda.j.allen@ttu.edu</footnote>
is extensive, although often not related to transmission. Indeed, extending epidemiological models to go beyond a single pathogen species was relatively-recently highlighted as a key challenge in modeling plant diseases in [12] (challenge 4). Co-infection almost always leads to interactions between viruses during transmission and within-plant processes that can strongly influence disease development in individual plants and ultimately spread in a plant population. The strength and direction of interactions can vary with both negative and facilitating effects, involving within-cell processes, cell-to-cell movement, vector acquisition and inoculation, symptom development and virulence, and yield loss. Reports on cellular interactions have been the most prevalent, mostly for replication rates and virus titre, but some studies have shown clear interactions with vectors over short (epidemiological) and long-term (evolutionary) time scales [13].

Most experimental studies on the relationship between co-infection and transmission have been done for viruses with non-persistent transmission by aphids. Syller [14] reviewed the literature on "simultaneous" transmission of plant viruses by vectors, emphasizing the acquisition component of transmission with little consideration of inoculation. The use of the term "simultaneous" is ambiguous—what seems to be suggested is that two different virus particles can be acquired instantaneously by a vector during a single probe where there is spatial separation between two viruses. However, for non-persistent transmission, with a following probe [15] on the same or different plant, one of the viruses can become detached and no longer be available for inoculation. Transmission can present a real bottleneck in the virus life cycle [16]. In a subsequent review, Syller and Grupa [17] differentiate between simultaneous inoculation (which they call co-infection) and sequential inoculation (which they call super-infection). They claim that synergistic interactions within-plants most often arise between unrelated viruses. Synergism is defined as a facilitative effect in which accumulation of one or both viruses in the host plant increases; in the case of the effect on just one virus, it has been called asymmetric synergism [18]. Synergism has also been used to describe more severe disease symptoms than induced by either virus alone. Syller and Grupa [17] concentrate more on antagonistic effects, such as cross protection [19] or, as has been termed "super-infection exclusion" in which related viruses or virus strains are used preventively to exclude more virulent strains. Mascia and Gallitelli [20] note the contributions that mathematical modeling could make "in forecasting challenges deriving from the great variety of pathways of synergistic and antagonistic interactions" (p. 176).

Co-infection can cover scenarios ranging from two viruses (or virus strains)/one vector through to many viruses/many vectors, but with some nuances. There is an extensive literature on co-infection across this range. Some representative but not exhaustive publications are noted in Table 1, together with some key messages. Many publications acknowledge that there are several or many (in the case of aphids) vector species for a given plant virus, but the experiments reported only involve one vector. Similarly, the same virus and vector can infect more than one host (cucumber mosaic virus is an extreme example) and hence cause more than one disease. Co-infection with virus strains differing in virulence (or other characteristic) can lead to the same set of interactions and consequences as found with virus species. A good example of two strains of the same virus species with shared vector species is potato virus Y (PVY) [21–23] on potato and other hosts [15]. There are many examples of two co-infecting virus species with a shared vector species [24–26]. Similarly, there are many cases where two co-infecting viruses have quite different vectors taxonomically [27]. Co-infection is manifested in more complex situations with multiple viruses and vectors such as with grapevine leafroll disease [28,29] and sweet potato virus disease (SPVD) [30]. At an even higher level of complexity, the ecological networks formed by multiple co-infecting viruses and multiple hosts were analyzed by McLeish et al. [31].
<table><caption>**Table 1.** Summary of some representative experimental work on co-infection at different levels of complexity.</caption><thead><tr><th></th><th>Pathosystem</th><th>Authors</th><th>Key message</th></tr></thead><tbody><tr><td>Two virus strains</td><td>Potato virus Y (PYV) —<i>Myzus persicae</i> (host <i>Capsicum annuum</i>)</td><td>Moury et al. [15]</td><td>Two strains were equally transmissible and competition was studied to estimate the size of bottlenecks imposed by vector transmission. If there was a cost of virulence, modelling showed that virulent strains would go extinct.</td></tr><tr><td></td><td>PVY<sup>NTN</sup> compared with other strains—<i>Myzus persicae</i></td><td>Srinivasan et al. [21]</td><td>Previous work had suggested some specificity in transmission of strains. The rate of infection for PVY<sup>NTN</sup> was higher than for other strains, a vector-related outcome as this was not observed with mechanical transmission.</td></tr><tr><td></td><td>PVY<sup>NTN</sup> compared with PVY<sup>O</sup>—<i>Myzus persicae</i></td><td>Carroll et al. [22]</td><td>The necrotic strain was transmitted more efficiently than the wild-type. Co-infection would more likely result from inoculation by multiple aphids feeding on plants infected with the different strains rather than by single aphids feeding on multiple plants infected with the different strains.</td></tr><tr><td>Two virus species—common vector</td><td>Barley yellow dwarf virus/cereal yellow dwarf virus—<i>Rhopalosiphum padi</i></td><td>Lacroix et al. [24]</td><td>The co-inoculation of BYDV-PAV lowered the CYDV-RPV infection rate but only at low nutrient supply rates. Broader environmental and nutritional factors can alter co-infection interactions and outcomes.</td></tr><tr><td></td><td>Watermelon mosaic virus/zucchini yellow mosaic virus—<i>Aphis gossypii</i></td><td>Salvaudon et al. [25]</td><td>ZYMV accumulated at similar rates in single and mixed infections, whereas WMV was much reduced in the presence of ZYMV. ZYMV also induced host changes that gave strong vector preference for infected plants; whereas WMV did not, although it was still readily acquired from mixed infections.</td></tr><tr><td></td><td>Rice tungro spherical virus/rice tungro bacilliform virus—<i>Nephottetix virescens</i></td><td>Holt and Chancellor [26]</td><td>Infection by each virus alone results in less pronounced symptoms. RTBV is retained in the vector for a longer period. When a vector carries both viruses, co-inoculation is common. When inoculative with RSTV alone the infection probability is higher.</td></tr><tr><td>Two virus species—multiple vector species</td><td>Bean pod mottle virus—<i>Epilachna varivestis</i>/soybean mosaic virus—<i>Aphis glycines</i></td><td>Peñaflor et al. [27]</td><td>Singly-infected plants with either BPMV or SMV increased soybean palatability, potentially enhancing acquisition of BPMV from BPMV plants and secondary infection of BPMV from SMV plants. BPMV infection had little effect on <i>A. glycines</i>, whereas SMV infection reduced aphid population growth but increased the preference for infected plants. With co-infection, effects on population growth were reversed and aphids showed a preference for co-infected plants.</td></tr><tr><td>Multiple virus species—multiple vector species</td><td>Grapevine leafroll-associated viruses (GLRaVs)—mealybugs/scale insects</td><td>Naidu et al. [28]</td><td>The exact role of GLRaVs in disease etiology remains unclear. With mealybugs, transmission is of a semi-persistent manner with a lack of vector-virus specificity.</td></tr><tr><td></td><td>Sweet potato chlorotic stunt virus/sweet potato feathery mottle virus/multiple viruses—multiple vector species</td><td>Untiveros et al. [30]</td><td>Six viruses from the same or different virus families interacted synergistically with sweet potato virus disease, with increased disease symptoms, virus accumulation and movement in plants, and reduced yield of storage roots. All inoculations were made by grafting; no conclusions can be drawn on vector transmission effects.</td></tr><tr><td>Multiple virus species—multiple vector species, multiple hosts</td><td>Ecological networks formed by multiple co-infecting viruses in multiple hosts</td><td>McLeish et al. [31]</td><td>Co-infection networks were found to lead to strong non-random associations compared with single infections. Single infections were mostly related to habitat parameters, whereas co-infections were more related to ecological heterogeneity and ecosystem-level processes.</td></tr></tbody></table>
In this paper, we formulate a general epidemiological model for one vector species and one plant species that allows for co-infection of the host plant by two virus species or strains. The model is used to investigate the role of vector transmission on co-infection, specifically acquisition and inoculation, as well as antagonistic and synergistic interactions. The basic reproduction number provides a condition for invasion of a single virus infection. For co-infection, we derive a new invasion threshold. Given that a single virus can persist in the host plant, the invasion threshold highlights some key epidemiological parameters for successful co-infection. In addition, we investigate the roles of the vector acquisition and inoculation parameters when one virus depends on an autonomous virus for its successful transmission, or when both viruses are unable to invade alone but can facilitate co-infection if they occur in high enough prevalences. We also explicitly test when the simplifications—which are almost always left implicit—in models which do not explicitly include the infection status of vectors lead to potentially misleading results.

## 2. Materials and Methods

### 2.1. Modelling

The general epidemiological model for vectors and plants consists of a system of differential equations with either infection by a single virus, or co-infection by two viruses (species or strains). For simplicity, we refer to the two viruses as virus A and virus B. We make several simplifying assumptions. The two viruses are not transmitted vertically in the plant population (no transmission by seeds or other propagating material) nor in the vector (no transovarial transmission). The plants can be infected by a single virus or co-infected by both viruses. The vectors are only inoculative with a single virus; namely, we assume that acquisition of the first virus precludes inoculation by a second virus by the same vector for as long as the vector retains the first virus. This assumption applies to both modes of transmission. For non-persistent transmission, it may arise because stylet receptor sites are saturated by the first virus. For persistent-circulative transmission, it may be due to a phenomenon similar to "super-infection exclusion" taking place in the vector [16], or simply that the first virus moves back to the salivary glands faster than the second. In addition, the latent stages in the vector and plant are ignored. The following compartmental diagrams in Figure 1 illustrate the rates of change between the vector and plant stages. The vector and plant models are described in more detail below.

<figure><img src="image_4.png" /><figcaption>**Figure 1.** Compartmental diagrams (left) for the vector model and (right) the host plant model, described by the differential Equations (1) and (2), respectively. Meanings of the symbols for model parameters are summarised in Table 2.</figcaption></figure>

The vector model has three stages, X = the density of non-infective vectors, $Z_A$ = the density of infective vectors carrying virus A and $Z_B$ = the density of infective vectors carrying virus B (for vectors, "infective" more accurately means "inoculative"). The total vector density is $V = X + Z_A + Z_B$. The parameters $\Lambda$ = net vector birth rate, $c$ = per capita vector death rate, $\Phi$ = number of plants
visited per unit time by a single vector, $\alpha_i$ = probability a non-infective vector acquires a single virus from an infected plant $I_i$, $i = A, B, AB$, per plant visit and $\delta_i$ = per capita vector recovery rate from virus $i = A, B$. The parameters $\epsilon_A$ and $\epsilon_B$ multiplying the acquisition probability $\alpha_{AB}$ are the conditional probabilities that either virus $A$ or $B$ are acquired, given that acquisition is from a co-infected plant, $\epsilon_A + \epsilon_B = 1$. All of the parameters are non-negative. They are summarised in Table 2. With these assumptions, as well as the assumption of frequency-dependent transmission with $P$ being the total plant population density, the vector model takes the following form:

$$ \text{Vector} \quad \left\{ \begin{aligned} \frac{dX}{dt} &= \Lambda - \Phi X [\alpha_A I_A + \alpha_B I_B + \alpha_{AB} I_{AB}] / P - cX + \delta_A Z_A + \delta_B Z_B \\ \frac{dZ_A}{dt} &= \Phi X [\alpha_A I_A + \epsilon_A \alpha_{AB} I_{AB}] / P - (c + \delta_A) Z_A \\ \frac{dZ_B}{dt} &= \Phi X [\alpha_B I_B + \epsilon_B \alpha_{AB} I_{AB}] / P - (c + \delta_B) Z_B. \end{aligned} \right. \quad (1) $$

The plant model consists of four stages, $S$ = density of healthy plants, and three different classes for infection, $I_A$, $I_B$ and $I_{AB}$, equal to the density of infected plants carrying only virus $A$, virus $B$ or both viruses, respectively. The total plant population density is $P = S + I_A + I_B + I_{AB}$. The parameters for the plant model include $\sigma$ = net planting rate, $\mu$ = per capita rate of harvesting/mortality, $\theta_i$ = probability an infective vector $Z_i$ inoculates a healthy plant with virus $i = A, B$ per plant visit, $\omega$ = per capita recovery rate in a plant infected with a single virus, and $\omega_i$ = per capita viral recovery rate of virus $i$ in a co-infected plant $I_{AB}$, $i = A, B$. The parameters $\gamma_A$ and $\gamma_B$ denote the synergistic or antagonistic interactions between the viruses within the plant when infective vectors inoculate a plant carrying a different virus, i.e., account for the probability of a second infection being either increased or reduced relative to a healthy plant. The plant model is

$$ \text{Plant} \quad \left\{ \begin{aligned} \frac{dS}{dt} &= \sigma - \Phi S [\theta_A Z_A + \theta_B Z_B] / P - \mu S + \omega (I_A + I_B) \\ \frac{dI_A}{dt} &= \Phi \theta_A Z_A S / P - \Phi \gamma_B \theta_B Z_B I_A / P - (\mu + \omega) I_A + \omega_A I_{AB} \\ \frac{dI_B}{dt} &= \Phi \theta_B Z_B S / P - \Phi \gamma_A \theta_A Z_A I_B / P - (\mu + \omega) I_B + \omega_B I_{AB} \\ \frac{dI_{AB}}{dt} &= \Phi [\gamma_B \theta_B Z_B I_A + \gamma_A \theta_A Z_A I_B] / P - (\mu + \omega_A + \omega_B) I_{AB}. \end{aligned} \right. \quad (2) $$

The parameters $\epsilon_i$ and $\gamma_i$, $i = A, B$, reflect the fact that vector acquisition or inoculation with co-infection may differ from a single virus infection [24,25,30]. The acquisition of a single virus from a co-infected plant $I_{AB}$ may be greater or less than acquisition from a plant infected with a single virus (e.g., $\epsilon_A \alpha_{AB} > \alpha_A$ or $\epsilon_A \alpha_{AB} < \alpha_A$). In addition, the inoculation of a second virus into an infected plant may be greater or less than inoculation of a healthy plant (e.g., $\gamma_B \theta_B > \theta_B$ or $\gamma_B \theta_B < \theta_B$).
<table><caption>**Table 2.** Summary of model parameters for plants and vectors. All parameters are non-negative.</caption><thead><tr><th>Vector Parameters</th><th>Default Values Sections 3.1 and 3.2, Appendices B and C</th><th>Default Values Section 3.3</th></tr></thead><tbody><tr><td>Λ vector birth rate</td><td>10/month/area</td><td>1/month/area</td></tr><tr><td>c per capita vector natural death rate</td><td>1/month</td><td>1/month</td></tr><tr><td>Φ number of plants visited/time by a vector</td><td>1/day</td><td>8.33/day</td></tr><tr><td>δ<sub>A</sub> per capita infective vector recovery rate from virus A</td><td>3/month</td><td>0.66/day</td></tr><tr><td>δ<sub>B</sub> per capita infective vector recovery rate from virus B</td><td>3/month</td><td>1.33/day</td></tr><tr><td>α<sub>A</sub> probability non-infective vector acquires virus A from I<sub>A</sub> per plant visit</td><td>0.2</td><td>0.005</td></tr><tr><td>α<sub>B</sub> probability non-infective vector acquires virus B from I<sub>B</sub> per plant visit</td><td>0</td><td>0.005</td></tr><tr><td>α<sub>AB</sub> probability non-infective vector acquires a single virus, A or B, from I<sub>AB</sub> per plant visit</td><td>0.5</td><td>0.15</td></tr><tr><td>ε<sub>A</sub> from a co-infected plant I<sub>AB</sub>, given a successful acquisition</td><td>0.5</td><td>0.5</td></tr><tr><td>ε<sub>B</sub> conditional probability of acquiring virus B from co-infected plant I<sub>AB</sub>, given a successful acquisition (ε<sub>B</sub> = 1 - ε<sub>A</sub>)</td><td>0.5</td><td>0.5</td></tr></tbody></table><table><thead><tr><th>Plant Parameters</th><th>Default Values Sections 3.1 and 3.2, Appendices B and C</th><th>Default Values Section 3.3</th></tr></thead><tbody><tr><td>μ per capita mortality and or harvest of plants</td><td>1/year</td><td>1/year</td></tr><tr><td>σ seeding or planting rate</td><td>1000/year/area</td><td>100/year/area</td></tr><tr><td>θ<sub>A</sub> probability an infective vector with virus A inoculates a healthy plant per visit</td><td>0.8</td><td>0.8</td></tr><tr><td>θ<sub>B</sub> probability an infective vector with virus B inoculates a healthy plant per visit</td><td>0.5</td><td>0.5</td></tr><tr><td>γ<sub>A</sub> relative inoculation success of virus A (as compared to a heathy plant) in a plant I<sub>B</sub>, infected with a single virus B</td><td>0.9</td><td>0.5</td></tr><tr><td>γ<sub>B</sub> relative inoculation success of virus B (as compared to a healthy plant) in a plant I<sub>A</sub>, infected with a single virus A</td><td>0.25, 0.9</td><td>0.5</td></tr><tr><td>ω per-capita viral A or B loss rate in a plant infected with single virus</td><td>0</td><td>0.001/day</td></tr><tr><td>ω<sub>A</sub> per-capita viral B loss rate (A is retained) from a co-infected plant I<sub>AB</sub></td><td>0</td><td>0.001/day</td></tr><tr><td>ω<sub>B</sub> per-capita viral A loss rate (B is retained) from a co-infected plant I<sub>AB</sub></td><td>0</td><td>0.001/day</td></tr></tbody></table>

## 2.2. Invasion Thresholds

Two important disease threshold parameters are derived from models (1) and (2), the basic reproduction number and the invasion reproduction number. The density of non-infective vectors at the disease-free equilibrium (DFE) in the vector model (1) is $X = Λ/c = \bar{V}$ and the density of healthy plants at the DFE plant model (2) is $S = \sigma/\mu = \bar{P}$. The basic reproduction number can be computed from the next generation matrix approach [32–36] (Appendix A). Here it is defined as the maximum of two reproduction numbers,

$$ \mathcal{R}_0 = \max \left\{ \frac{\Phi \bar{V} \alpha_A \Phi \theta_A}{\bar{P}(c + \delta_A)(\mu + \omega)}, \frac{\Phi \bar{V} \alpha_B \Phi \theta_B}{\bar{P}(c + \delta_B)(\mu + \omega)} \right\} = \max\{\mathcal{R}_{0A}, \mathcal{R}_{0B}\}. \quad (3) $$

The two terms in the preceding definition are basic reproduction numbers corresponding to infection with either virus A or virus B, $\mathcal{R}_{0A}$ and $\mathcal{R}_{0B}$, respectively. If the basic reproduction number for virus A exceeds the value of one, $\mathcal{R}_{0A} > 1$, then virus A can invade the disease-free vector-plant system and if $\mathcal{R}_{0B} > 1$, then virus B can invade. An epidemiological interpretation of $\mathcal{R}_{0A}$ is that if one vector inoculative with virus A is introduced into a healthy vector-plant system, it will inoculate and infect plants at a rate of $\Phi\theta_A$ during the period of time $1/(c + \delta_A)$ the vector is infective. From an
infected plant, a non-infective vector will acquire the virus at a rate $\Phi\alpha_A\bar{V}/\bar{P}$ during the period of time $1/(\mu + \omega)$ the plant is infected. If the product of these two expressions exceeds the value of one, then one infective vector (or infected plant) will generate more than one infective vector (or infected plant), resulting in an epidemic. In general, if $R_0 > 1$, then either virus A or B can invade the vector-plant system.

An invasion reproduction number can be derived if the system is already infected with a single virus. We consider whether virus B can invade when the system is at the virus A equilibrium, i.e., $R_{0A} > 1$. The endemic equilibrium values for $(X, Z_A, S, I_A)$ with virus A are $(\bar{V} - za_{eq}, za_{eq}, P - ia_{eq}, ia_{eq})$. These endemic values can be expressed in terms of the basic reproduction number $R_{0A}$:

$$za_{eq} = \frac{(R_{0A} - 1)\bar{V}}{R_{0A} + \frac{\Phi\theta_A\bar{V}}{\bar{P}(\mu + \omega)}} \quad \text{and} \quad ia_{eq} = \frac{(R_{0A} - 1)\bar{P}}{R_{0A} + \frac{\Phi\alpha_A}{(c + \delta_A)}} \quad (4)$$

To determine whether virus B can invade, we apply the next generation matrix approach to derive an invasion matrix $M_{invB}$ from which an invasion reproduction number $R_{invB}$ can be computed (Appendix B). We assume that the vector-plant system is at the virus A equilibrium, defined in (4), with the remaining states $Z_B, I_B$, and $I_{AB}$ set equal to zero. An invasion matrix for the system (1)–(2) is defined as follows:

$$M_{invB} = \begin{pmatrix} 0 & \frac{\Phi \bar{X} \alpha_B}{\kappa_1} & \frac{\Phi \bar{X} \epsilon_B \alpha_{AB}}{\bar{P} \kappa_2} \\ \frac{\Phi \theta_B \bar{S}}{\bar{P}(c + \delta_B)} & 0 & \frac{\omega_B}{\kappa_2} \\ \frac{\Phi \gamma_B \theta_B ia_{eq}}{\bar{P}(c + \delta_B)} & \frac{\Phi \gamma_A \theta_A za_{eq}}{\kappa_1} & 0 \end{pmatrix}, \quad (5)$$

where $\bar{S} = \bar{P} - ia_{eq}$, $\bar{X} = \bar{V} - za_{eq}$, $\kappa_1 = (\mu + \omega)\bar{P} + \Phi\gamma_A\theta_Aza_{eq}$ and $\kappa_2 = \mu + \omega_A + \omega_B$. A mathematical definition of the invasion reproduction number is the spectral radius of the invasion matrix, that is, $R_{invB} = \rho(M_{invB})$. An epidemiological interpretation of $R_{invB}$ is the average number of new states, $Z_B, I_B$ or $I_{AB}$, that are produced after introduction of an average of one infective vector or infected plant containing virus B, $Z_B, I_B$ or $I_{AB}$, into the system infected with virus A. If $R_{invB} > 1$, then virus B can invade and if $R_{invB} < 1$, then virus B cannot invade. Each element in matrix $M_{invB}$ can be interpreted in terms of producing new infective vectors or infected plants, $Z_B, I_B$ or $I_{AB}$. For example, the entry in the second row and first column,

$$\frac{\Phi \theta_B \bar{S}}{\bar{P}(c + \delta_B)}$$

can be interpreted as the average number of new infected plants $I_B$ that are produced when one infective vector $Z_B$ is introduced into the vector-plant system (where virus A has already invaded).

Since the invasion matrix (5) has non-negative entries, an increase (or a decrease) in any matrix entry also increases (or decreases) the invasion reproduction number [37]. In particular, if the inoculation probability $\theta_B$ for virus B increases so does the invasion reproduction number. The direction of change for the invasion reproduction number is not as straightforward if the vector visitation rate $\Phi$ or parameters related to virus A are changed, as the equilibrium values $ia_{eq}$ and $za_{eq}$ also change.

## 2.3. Formal Reduction to a Model That Does Not Track Vectors Explicitly

The vector-plant system can be reduced to a simplified plant model, where the vectors are not explicitly included (Appendix C). If the vector recovery rates, $\delta_A$ and $\delta_B$, are sufficiently large, then the vector population dynamics occur on a faster time scale than the plant population dynamics.
The differential equations for the vector model (1) and the plant host model (2) can be approximated by another plant model without the vector variables:

$$ \frac{dS}{dt} = \sigma - \frac{\Phi^2 \bar{V}}{p^2} S (\hat{\theta}_A [\alpha_A I_A + \epsilon_A \alpha_{AB} I_{AB}] + \hat{\theta}_B [\alpha_B I_B + \epsilon_B \alpha_{AB} I_{AB}]) - \mu S + \omega (I_A + I_B) $$

$$ \frac{dI_A}{dt} = \frac{\Phi^2 \bar{V}}{p^2} (\hat{\theta}_A S [\alpha_A I_A + \epsilon_A \alpha_{AB} I_{AB}] - \gamma_B \hat{\theta}_B I_A [\alpha_B I_B + \epsilon_B \alpha_{AB} I_{AB}]) - (\mu + \omega) I_A + \omega_A I_{AB} \quad (6) $$

$$ \frac{dI_B}{dt} = \frac{\Phi^2 \bar{V}}{p^2} (\hat{\theta}_B S [\alpha_B I_B + \epsilon_B \alpha_{AB} I_{AB}] - \gamma_A \hat{\theta}_A I_B [\alpha_A I_A + \epsilon_A \alpha_{AB} I_{AB}]) - (\mu + \omega) I_B + \omega_B I_{AB} $$

$$ \frac{dI_{AB}}{dt} = \frac{\Phi^2 \bar{V}}{p^2} (\gamma_B \hat{\theta}_B I_A [\alpha_B I_B + \epsilon_B \alpha_{AB} I_{AB}] + \gamma_A \hat{\theta}_A I_B [\alpha_A I_A + \epsilon_A \alpha_{AB} I_{AB}]) - (\mu + \omega_A + \omega_B) I_{AB}, $$

where parameters $\hat{\theta}_i = \theta_i / (c + \delta_i)$ for $i = A, B$. To distinguish the two sets of models, we will refer to the simplified plant model (6) as the vector-implicit model and to the vector-plant model (1) and (2) as the vector-explicit model.

The vector-implicit model retains the inoculation and acquisition terms. At the DFE, the basic reproduction number is the same as in (3). But the equilibrium for virus A ($\bar{S}, \bar{I}_A$) differs from the vector-explicit model,

$$ \bar{S} = \frac{\bar{P}}{\mathcal{R}_{0A}} \quad \text{and} \quad \bar{I}_A = \bar{P} \left(1 - \frac{1}{\mathcal{R}_{0A}}\right). \qquad (7) $$

An invasion matrix for the vector-implicit model when $\mathcal{R}_{0A} > 1$ can be computed by a method similar to the computation of the invasion matrix in (5) (Appendix C). Invasion of virus B (via plants infected with virus B, $I_B$ or $I_{AB}$) is successful if the spectral radius of the following invasion matrix exceeds the value of one:

$$ M_{invB} = \left( \frac{\frac{\Phi^2 \bar{V} \alpha_B \hat{\theta}_B \bar{S}}{\eta}}{\frac{\Phi^2 \bar{V} \bar{I}_A [\gamma_A \alpha_A \hat{\theta}_A + \gamma_B \alpha_B \hat{\theta}_B]}{\eta}} \frac{\Phi^2 \bar{V} \hat{\theta}_B \alpha_{AB} \epsilon_B \bar{S} + \bar{P}^2 \omega_B}{\frac{\bar{P}^2 (\mu + \omega_A + \omega_B)}{\Phi^2 \bar{V} \gamma_B \hat{\theta}_B \alpha_{AB} \epsilon_B \bar{I}_A} \frac{\bar{P}^2 (\mu + \omega_A + \omega_B)}{\bar{P}^2 (\mu + \omega_A + \omega_B)}} \right), \quad (8) $$

$$ \text{where } \eta = \bar{P}^2 (\mu + \omega) + \Phi^2 \bar{V} \gamma_A \alpha_A \hat{\theta}_A \bar{I}_A. $$

### 3. Results

The effect of vector transmission on virus establishment in the vector-explicit model is examined in several numerical examples. In addition, the conditions for invasion of a second virus in the vector-explicit model are compared to the vector-implicit model.

#### 3.1. Only One Virus Can Invade in Absence of the Other

We assume that virus A can invade and persist in the vector-explicit model but not virus B. In particular, the acquisition probability of virus B is set to zero, $\alpha_B = 0$, and all other acquisition and inoculation probabilities are positive. We investigate the invasion reproduction number as a function of the remaining acquisition and inoculation parameters, $\alpha_A, \theta_A$ and $\theta_B$, whose values lie in $[0, 1]$. Table 2 is a summary of the default parameter values. The parameter values are chosen consistent with those summarised by Jeger et al. [5], but based on a time unit of one month. For each of the acquisition or inoculation parameters, Figure 2 shows a graph of the equilibrium prevalences for infective vectors and infected plants with virus A (with $z_a = za_{eq}/\bar{V}$ and $i_a = ia_{eq}/\bar{P}$) and the invasion reproduction curves. Two different invasion reproduction curves are graphed, one with $\gamma_B = 0.9$ (solid invasion curve) and the other with $\gamma_B = 0.25$ (dashed invasion curve). Parameter $\gamma_B$ denotes the inoculation success of virus B on a plant infected with virus A, relative to a healthy plant. Small values of $\gamma_B$ decrease the likelihood of a successful invasion of virus B.
<figure><img src="image_5.png" /><figcaption>**Figure 2.** For the vector-explicit model, the virus A equilibrium prevalences for the vector $z_a = za_{eq}/\bar{V}$ and for the host $i_a = ia_{eq}/\bar{P}$ are graphed in panels (A,C,E) and the corresponding invasion reproduction numbers $R_{invB}$ are graphed in panels (B,D,F) as a function of acquisition ($\alpha_A$) or inoculation ($\theta_A$ and $\theta_B$) parameters. Virus B may invade the virus A equilibrium if the invasion reproduction number is greater than the threshold value one. The parameters values that are not varied are fixed at their default values (see Table 2): $\Lambda = 10$ per month per unit area, $\Phi = 30$ per month, $c = 1$ per month, $\alpha_A = 0.2$, $\alpha_B = 0$, $\alpha_{AB} = 0.5$, $\epsilon_A = \epsilon_B = 0.5$, $\delta_A = \delta_B = 3$ per month, $\mu = 1/12$ per month, $\sigma = 1000\mu$ (1000 per year per unit area), $\theta_A = 0.8$, $\theta_B = 0.5$ and $\omega = \omega_A = \omega_B = 0$. Parameter values for the two invasion reproduction curves in (B,E,F) are $\gamma_B = 0.9$ for the solid curve and $\gamma_B = 0.25$ for the dashed curve. (B) The two black circles on the solid curve are at values of $\alpha_A = 0.1, 0.8$; (D) $\theta_A = 0.3, 0.8$; (F) $\theta_B = 0.2, 0.8$. The invasion reproduction number is $R_{invB} = \rho(M_{invB})$ as defined in (5).</figcaption></figure>

As the acquisition probability $\alpha_A$ increases (Figure 2A,B), there is initially an increase in the invasion reproduction number but then it decreases as $\alpha_A$ approaches one. From the invasion matrix $M_{invB}$ in (5), it can be seen that some of the matrix entries increase and some decrease with an increase in $\alpha_A$, resulting in potential increases or decreases in the invasion reproduction number. The equilibrium prevalences $ia_{eq}$ and $za_{eq}$ increase with $\alpha_A$ which provide more opportunities for an infective vector $Z_B$ to inoculate infected plants $ia_{eq}$ and more opportunities for an infected plant $I_B$ to become inoculated from infective vectors $za_{eq}$. But the density of healthy plants, $\bar{S}$ in matrix (5), and the density of non-infective vectors, $\bar{X}$ in matrix (5), decrease with increases in $\alpha_A$ which provide fewer opportunities for $Z_B$ to inoculate healthy plants or for non-infective vectors to acquire virus B from a co-infected plant $I_{AB}$. For a relative inoculation success ratio of $\gamma_B = 0.9$ (solid invasion curve), the invasion reproduction number exceeds the threshold value of one for a range of $\alpha_A$ values but for $\gamma_B = 0.25$ (dashed invasion curve), the invasion reproduction number never exceeds the value of one.

As the inoculation probability $\theta_A$ increases in Figure 2C,D, a similar effect on the invasion reproduction number could occur as for $\alpha_A$ but the decrease in the invasion reproduction number is not seen for these parameter values. The reduction in healthy plants and non-infective vectors is not as severe as for $\alpha_A$. Changes in the inoculation probability $\theta_B$ do not affect equilibrium prevalences $ia_{eq}$ and $za_{eq}$. Therefore, as $\theta_B$ increases, so does the invasion reproduction number (Figure 2E,F). Similarly, for large values of $\gamma_B$, the invasion reproduction number increases (Figure A1 in Appendix B).
Parameter values in Figure 2, indicated by the black circles on the invasion curves, represent either a successful ($R_{invB} > 1$) or an unsuccessful ($R_{invB} < 1$) invasion of virus B. For these parameter values, the infection prevalence of viruses A and B in the plant population are graphed as functions of time in Figure A2 (Appendix B).

## 3.2. Comparison with Results of Model That Does Not Track Vectors Explicitly

The dynamics of the vector-implicit model in Figure 3 are compared to the vector-explicit model when the same parameter values as in Figure 2 are applied. The differences in the invasion matrices and the virus A equilibrium values in the two models impact the invasion outcome. For the default parameter values in Table 2, the infected plant prevalences and the invasion reproduction values are larger in the vector-implicit model than in the vector-explicit model. The vector-implicit model predicts greater likelihood of co-infection than in the vector-explicit model. The differences may be attributed to the fact that the assumptions which led to the reduction to a plant model do not hold for all transmission classes. One of the assumptions in the vector-implicit model is that the infective vector life cycle is short, such as in non-persistent transmission. But in these examples, the parameter values for the vectors' infective period, $1/\delta_A$ and $1/\delta_B$, are equal to 10 days, characteristic of persistent-circulative transmission [5]. These differences between the vector-explicit model and the vector-implicit model illustrate the importance of considering the mode of vector transmission in models of vectored plant viruses.

In Appendix C, additional invasion reproduction curves are graphed for different values of the relative inoculation success of viruses A and B, $\gamma_A$ and $\gamma_B$ in Figure A3 and as a function of vector recovery $\delta_A$ in Figure A5. Also, the infection prevalence of viruses A and B in the plant population for parameter values indicated by the black circles in Figure 3 are graphed as functions of time in Figure A4.

<figure><img src="image_6.png" /><figcaption>**Figure 3.** For the vector-implicit model, the virus A host equilibrium prevalences $ia_{eq}/\bar{P}$ are graphed in panels (A,C,E) and the corresponding invasion reproduction numbers $R_{invB}$ are graphed in panels (B,D,F) as a function of the acquisition ($\theta_A$) and inoculation ($\theta_A$ and $\theta_B$) parameters. The same parameter values are applied as in the vector-explicit model in Figure 2. The three black circles on the solid curve in (B) are at values of $\alpha_A = 0.05, 0.1, 0.8$; in (D) $\theta_A = 0.2, 0.3, 0.8$; (F) $\theta_B = 0.05, 0.2, 0.8$. The invasion reproduction number equals $R_{invB} = \rho(M_{invB})$ as defined in (8).</figcaption></figure>
### 3.3. Neither Virus Can Invade in Absence of the Other

In this example, we assume neither virus can invade the vector-explicit model in the absence of the other virus, i.e., $R_{0A} < 1$ and $R_{0B} < 1$. There is only the disease-free equilibrium (DFE), which is always locally stable, and, for certain parameter values, we found numerically up to two co-infection equilibria. One of the latter is stable, so that the vector-explicit model is bistable. In the following analysis, we define $f_a$ and $f_b$ as the initial prevalences of virus A and virus B in the plant population, i.e., $f_a = (I_A + I_{AB})/P \in [0, 1]$ and $f_B = (I_B + I_{AB})/P \in [0, 1]$, respectively. For these prevalences, we will assume the following initial conditions of the plant population

$$S(0) = (1 - f_a)(1 - f_b)\bar{P}, \quad I_A(0) = f_a(1 - f_b)\bar{P}, \quad I_B(0) = f_b(1 - f_a)\bar{P}, \quad \text{and} \quad I_{AB}(0) = f_a f_b \bar{P}, \quad (9)$$

where $\bar{P} = \sigma/\mu$. For the initial conditions of the vector population, we assume

$$X(0) = \left( 1 - \frac{f_a^2 + f_b^2}{f_a + f_b} \right) \bar{V}, \quad Z_A(0) = \frac{f_a^2}{f_a + f_b} \bar{V}, \quad \text{and} \quad Z_B(0) = \frac{f_b^2}{f_a + f_b} \bar{V}, \quad (10)$$

where $\bar{V} = \Lambda/c$, and which accounts for the fact that the initial prevalences of the viruses in the vectors cannot add beyond one because there is no co-infection.

Figure 4A presents the bistable scenario and graphs the basins of attraction for the two stable equilibria as a function of the initial prevalences of virus A and virus B. In this graph, the origin represents the DFE, and the point representing the stable co-infection equilibrium indicates the virus prevalences at this equilibrium. The basin of attraction of an equilibrium is the set of initial conditions that will lead to that equilibrium in the long-run. As the system is bistable, there are two basins of attraction: one for the DFE which comprises smaller initial virus prevalences, and one for the co-infection state which comprises larger initial virus prevalences. The two basins of attraction are separated by a curve, the so-called separatrix. The time plots in Figure 4B illustrate that initial virus prevalences on one side of the separatrix approach the DFE in the long-run, while initial virus prevalences from the other side approach the co-infection state in the long-run.

Figure 5A shows the separatrix for varying values of the vector mortality, $c$. The basin of attraction of the co-infection equilibrium expands in size when vector mortality is decreased. The size of the basin of attraction is also a measure of the resilience of the corresponding equilibrium, because the larger the basin of attraction, the more resistant the equilibrium against perturbations. That is, decreasing vector mortality enhances the resilience of the co-infection equilibrium, whereas increasing vector mortality enhances the resilience of the DFE. Now consider a given initial condition marked by the star in Figure 5A. Depending on the value of vector mortality, the same initial condition leads either to the disappearance of both viruses from the system (for larger vector mortalities) or to co-infection (for smaller vector mortalities). This is shown in Figure 5B. Notably, there is no gradual transition between these different outcomes. Instead, there is a rather drastic establishment of co-infection at high prevalence levels once both viruses can persist.

Note that the separatrix between the two basins of attraction is asymmetric for the set of parameter values considered thus far. This is due to the different transmission potentials of the two viruses, as expressed by the different values of $R_{0A} = 0.95$ and $R_{0A} = 0.02$ (when $c = 1$). Figure 6 demonstrates that the 'asymmetry' in the basins of attraction can be 'reversed' when, for instance, the value of $\epsilon_A = 1 - \epsilon_B$ is decreased. This increases the conditional probability of vectors to acquire virus B from co-infected plants. That is, an advantage for virus B in the competition for vectors can compensate for the disadvantages of virus B assumed in the other parameter values.
<figure><img src="image_7.png" /><figcaption>**Figure 4.** (A) The basins of attraction of the disease-free equilibrium (DFE) (white region) and the co-infected equilibrium (shaded region) are graphed as a function of the initial frequency *f*<sub>*a*</sub> of virus *A* and frequency *f*<sub>*b*</sub> of virus *B* in the plant population. See Equations (9) and (10) for the initial conditions. The points mark the prevalences of virus *A* and *B* in the plant population at the DFE (at the origin) and the stable co-infection equilibrium (in the interior). The gray square marks the relative frequency of virus *A* and *B* in the vector population at the stable co-infection equilibrium. The blue cross and red asterisk indicate initial conditions in different basins of attraction, for which time plots show convergence either to the DFE ($I_{AB} = 0$) or to the co-infection equilibrium ($I_{AB}/P > 0$) in (B). Parameter values are $\Lambda = 1, \Phi = 250, c = 1, \alpha_A = 0.005, \alpha_B = 0.0005, \alpha_{AB} = 0.15, \epsilon_A = \epsilon_B = 0.5, \delta_A = 20, \delta_B = 40, \mu = 1/12, \sigma = 100\mu, \theta_A = 0.8, \theta_B = 0.4, \omega = \mu/2 = \omega_A = \omega_B, \gamma_A = \gamma_B = 0.5$ such that $R_{0A} = 0.95$ and $R_{0B} = 0.02$.</figcaption></figure>
<figure><img src="image_8.png" /><figcaption>**Figure 5.** (A) The changes in the size of the basins of attraction for the DFE and the co-infected equilibrium are graphed as a function of the initial frequency *f*<sub>*a*</sub> of virus *A* and *f*<sub>*b*</sub> of virus *B* when the death rate of vectors decreases from *c* = 1.0 to *c* = 0.7, with values of *c* = 0.7, 0.8, 0.9, 1.0. The boundary separating the two basins of attraction at *c* = 0.7 is closest to the origin=DFE and at *c* = 1.0 is furthest from the origin. Other parameter values are as in Figure 4. The star symbol marks the initial condition for the disease progress curves in (B), which shows the time series of the prevalence of co-infected plants for varying levels of per-capita plant mortality, *c*. Initial conditions are fixed at (9) and (10) with initial virus prevalences *f*<sub>*a*</sub> = *f*<sub>*b*</sub> = 0.1 in the plant population.</figcaption></figure>
<figure><img src="image_9.png" /><figcaption>**Figure 6.** Basins of attraction for the DFE and the co-infection equilibrium for different values of $\epsilon_A$ and $\epsilon_B = 1 - \epsilon_A$. The blue curve with $\epsilon = 0.5$ corresponds to the separatrix shown in Figure 4 with $R_{0A} = 0.95$, $R_{0B} = 0.02$. Decreasing $\epsilon_A$, i.e. increasing $\epsilon_B$, reverses the asymmetry in the basins of attraction.</figcaption></figure>

**4. Discussion**

Co-infections are pervasive in plant virus epidemiology; yet, mathematical models keeping track of co-infections often leave vector dynamics implicit. This is understandable since (i) keeping track of co-infections makes the models less tractable mathematically, (ii) modeling vector dynamics explicitly may not be necessary relative to the research question addressed, and (iii) the biological knowledge or the data may be too scarce to reasonably increase the model's complexity [38]. Moreover, it is well known [10,39,40] that vector-implicit models may reasonably approximate vector-explicit models when vector dynamics can be considered to be fast with respect to epidemiological dynamics in the plant host population (such as in non-persistent transmission). Although see [41] for a different model of non-persistent transmission of viruses that goes beyond differential equations to track the epidemiological effects of vector dynamics.

However, we have showed that vector-implicit and vector-explicit models may yield qualitatively different results (Figures 2 and 3). The default parameter set we considered were appropriate for persistent-circulative transmission with vector infective periods of about 10 days (Table 2), and the qualitative difference in the models' outcomes highlights a key biological feature of the vector-explicit model that is not accounted for in its vector-implicit analogue. This key feature is competition for vectors [16]. Competition occurs because we assumed one vector cannot be inoculative with two viruses at the same time at the point of inoculation. This assumption may reflect the fact that the first acquired virus saturates stylet receptor sites in non-persistent transmission, or that it precludes second inoculation of other viruses as long as the vector retains the first virus. The latter phenomenon would resemble the so-called "super-infection exclusion" usually associated with the host plant rather than with the vector [16]. Nevertheless, competitive exclusion (or pre-emptive competition) of one virus by another can also occur within vectors [42–45]. This way, infective vectors are only able to inoculate one virus at the point of inoculation. Therefore, viruses indirectly compete for vectors. However, the vector-implicit model is valid only if the vector recovery rate is sufficiently large, meaning that the vector loses the virus very quickly after acquisition (as is the case in non-persistent viruses).
Therefore, it is very unlikely that a virus cannot be transmitted because the vector had already acquired another virus. These results highlight the importance of the mode of vector transmission in models of vectored plant viruses.

To illustrate the implications of our model, we considered two distinct special cases with two viruses (or virus strains) denoted *A* and *B*:

i Virus *A* is able to complete a full infection cycle, invade a disease-free population and settle at an endemic equilibrium in the absence of virus *B*. By contrast, virus *B* can systematically infect a plant but cannot be acquired by a vector in the absence of co-infection in the host plant. Successful invasion of virus *B* depends on a complex relation between vector acquisition and inoculation rates, the relative inoculation success of viruses *A* and *B* and the prevalence of infected hosts and infective vectors carrying virus *A*.

ii Virus *A* and virus *B* can both complete their infection cycles but are unable to invade a disease-free population in the absence of the other virus. The only biologically feasible endemic equilibrium is a co-infection equilibrium. The system will approach the co-infection equilibrium if the initial prevalences of virus *A* and *B* are sufficiently large, i.e., if the initial condition is in the basin of attraction of the co-infection equilibrium. Otherwise, if the initial virus prevalences are so low that the initial condition is in the basin of attraction of the disease-free equilibrium, both viruses will disappear from the system. The separatrix between these two different outcomes represents a curve of tipping points. On either side of these tipping points, we have contrasting dynamics, namely a disease-free versus a co-infected system. We have seen that increased vector mortality rates (e.g., due to vector control programs) moves the separatrix by making the co-infection equilibrium less resilient. This could lead to an abrupt (rather than gradual) extinction of co-infection.

In this paper, we have formulated a general epidemiological model for potential co-infection in a host plant by two virus species or strains. We used the model to investigate the effects of vector transmission on co-infection, specifically acquisition and inoculation, as well as antagonistic and synergistic interactions between viruses. We showed that reducing the model to a vector-implicit model can lose some of the key features gained when the vector is included explicitly (such as competition for vector). We derived a new invasion threshold that determines whether or not a second virus can invade a host population in which a first virus had successfully established. The invasion threshold highlights the key epidemiological parameters that are important for successful co-infection. However, the invasion threshold only applies near the virus *A* equilibrium when the invading virus *B* species/strain is at low prevalence levels. The dynamics of the models away from the virus *A* equilibrium can be quite complex and exhibit bistability, where both the virus *A* equilibrium and the co-infection equilibrium are stable. In this case, the initial prevalences of viruses *A* and *B* determine whether virus *B* can successfully invade. We also investigated the potential for co-infection and the conditions that need to be satisfied when one virus depends on an autonomous virus for its successful transmission or when both viruses are unable to invade alone. For persistent-circulative transmission, competition between viruses/virus strains, either direct or indirect, as they move through the vector is identified as a key and challenging area for further research that would improve modeling attempts to predict the epidemiological consequences of co-infection.

**Author Contributions:** M.J.J. and F.M.H. (Frédéric M. Hamelin) designed the study; L.J.S.A., V.A.B., N.J.C., F.M.H. (Frédéric M. Hamelin), F.M.H. (Frank M. Hilker), and M.J.J. developed the models; L.J.S.A. formally analyzed the models; L.J.S.A., V.A.B., F.M.H. (Frédéric M. Hamelin), and F.M.H. (Frank M. Hilker) performed the visualization and model validation; L.J.S.A., V.A.B., N.J.C., F.M.H. (Frédéric M. Hamelin), F.M.H. (Frank M. Hilker), and M.J.J. drafted the manuscript; all authors approve the manuscript.

**Funding:** This work was initiated during the Multiscale Vectored Plant Viruses Working Group at the National Institute for Mathematical and Biological Synthesis, supported by the National Science Foundation through NSF Award #DBI-1300426, with additional support from The University of Tennessee, Knoxville. This material is based upon research supported by the Thomas Jefferson Fund of the Embassy of France in the United States and the FACE Foundation.
**Acknowledgments:** L.J.S.A. thanks Texas Tech University for support. N.J.C. thanks Girton College, University of Cambridge for support. The authors thank the reviewers for their helpful suggestions.

**Conflicts of Interest:** The authors declare no conflict of interest.

## Appendix A. Basic Reproduction Number

Linearization of the system of equations for the infected states ($Z_A, I_A, Z_B, I_B, I_{AB}$) in models (1) and (2) and evaluated at the DFE, $X = \Lambda/c = \bar{V}$, $S = \sigma/\mu = \bar{P}$, with the remaining states equal to zero, yields the following Jacobian matrix:

$$J = \begin{pmatrix} -c - \delta_A & \frac{\Phi \bar{V} \alpha_A}{\bar{P}} & 0 & 0 & \frac{\Phi \bar{V} \epsilon_A \alpha_{AB}}{\bar{P}} \\ \Phi \theta_A & -\mu - \omega & 0 & 0 & \omega_A \\ 0 & 0 & -c - \delta_B & \frac{\Phi \bar{V} \alpha_B}{\bar{P}} & \frac{\Phi \bar{V} \epsilon_B \alpha_{AB}}{\bar{P}} \\ 0 & 0 & \Phi \theta_B & -\mu - \omega & \omega_B \\ 0 & 0 & 0 & 0 & -\mu - \omega_A - \omega_B \end{pmatrix}.$$

The two 2 × 2 matrices on the diagonal are used to calculate $\mathcal{R}_{0A}$ and $\mathcal{R}_{0B}$. For example, the next generation matrix for virus A is

$$F_A V_A^{-1} = \begin{pmatrix} 0 & \frac{\Phi \bar{V} \alpha_A}{\bar{P}(\mu + \omega)} \\ \frac{\Phi \theta_A}{c + \delta_A} & 0 \end{pmatrix},$$

where

$$\mathcal{F}_A = \begin{pmatrix} 0 & \frac{\Phi \bar{V} \alpha_A}{\bar{P}} \\ \Phi \theta_A & 0 \end{pmatrix} \quad \text{and} \quad \mathcal{V}_A = \begin{pmatrix} c + \delta_A & 0 \\ 0 & \mu + \omega \end{pmatrix}.$$

Matrix $\mathcal{F}_A$ includes the rate of new infections and $\mathcal{V}_A$ includes all of the other transfer rates. We define the basic reproduction number as

$$\mathcal{R}_{0A} = \frac{\Phi \bar{V} \alpha_A \Phi \theta_A}{\bar{P}(c + \delta_A)(\mu + \omega)}, \quad (A1)$$

which is the product of the off-diagonal elements in the next generation matrix. The basic reproduction number is often defined as the spectral radius of the next generation matrix, the square root of this value, $\sqrt{\mathcal{R}_{0A}}$ [36]. Both definitions have been applied in vector-host epidemic models [32,34–36]. Roberts and Heesterbeek [34] call the expression in (A1) the type reproduction number, whereas Heffernan et al. [32] call it the basic reproduction number and van den Driessche [35] calls it the target reproduction number. The terms “type” or “target” reproduction numbers are used in relation to a particular stage or type, where control measures may be targeted. The two expressions, with or without the square root, are equal at the threshold value of one. We choose the definition in (A1) for its simplicity and its straightforward epidemiological interpretation. The two expressions $1/(\mu + \omega)$ and $1/(c + \delta_A)$ represent the average infective period of a plant or a vector, respectively. The entry in the first row and second column is the average number of infective vectors $Z_A$ produced by one infected plant $I_A$ during the plants’ infective period and the entry in the second row first column is the average number of infected plants $I_A$ produced by one infective vector $Z_A$ during the vectors’ infective period. Therefore,

$$\mathcal{R}_{0A} = \frac{\Phi \bar{V} \alpha_A}{\bar{P}(\mu + \omega)} \times \frac{\Phi \theta_A}{c + \delta_A}.$$
## Appendix B. Invasion Reproduction Number

The system of equations for the infected states ($Z_A, I_A, Z_B, I_B, I_{AB}$) are linearized about the virus A endemic equilibrium. This leads to the Jacobian matrix,

$$J = \begin{pmatrix} J_{\text{stab}} & C \\ O & J_{\text{invB}} \end{pmatrix}. \qquad (\text{A2})$$

Matrix $J_{\text{stab}}$ is a 2 × 2 stable matrix for the states ($Z_A, I_A$) (provided $\mathcal{R}_{0A} > 1$), C is a 2 × 3 matrix, O is a 3 × 2 zero matrix, and $J_{\text{invB}}$ is the invasion matrix for virus B. For simplicity, let $\bar{X} = \bar{V} - za_{eq}$ and $\bar{S} = \bar{P} - ia_{eq}$. The two diagonal block matrices in J are defined as

$$J_{\text{stab}} = \begin{pmatrix} -(c + \delta_A) - \frac{\Phi \alpha_A i a_{eq}}{\bar{P}} & \frac{\Phi \bar{X} \alpha_A}{\bar{P}} \\ \frac{\Phi \theta_A \bar{X}}{\bar{P}} & -(\mu + \omega) - \frac{\Phi \theta_A z a_{eq}}{\bar{P}} \end{pmatrix}$$

and

$$J_{\text{invB}} = \begin{pmatrix} -c - \delta_B & \frac{\Phi \bar{X} \alpha_B}{\bar{P}} & \frac{\Phi \bar{X} \epsilon_B \alpha_{AB}}{\bar{P}} \\ \frac{\Phi \theta_B \bar{S}}{\bar{P}} & -(\mu + \omega) - \frac{\Phi \gamma_A \theta_A z a_{eq}}{\bar{P}} & \omega_B \\ \frac{\Phi \gamma_B \theta_B i a_{eq}}{\bar{P}} & \frac{\Phi \gamma_A \theta_A z a_{eq}}{\bar{P}} & -(\mu + \omega_A + \omega_B) \end{pmatrix}.$$

Invasion of virus B can be determined by computing the eigenvalues of matrix $J_{\text{invB}}$. Invasion is successful if matrix $J_{\text{invB}}$ is unstable, i.e., if matrix $J_{\text{invB}}$ has an eigenvalue with positive real part. However, computation of the three eigenvalues yield complicated expressions in terms of all of the parameters and do not give a simple criterion for invasion. An alternative is to compute an invasion reproduction matrix.

To compute the invasion reproduction number, we define $J_{\text{invB}} = \mathcal{F}_{\text{invB}} - \mathcal{V}_{\text{invB}}$. An invasion reproduction number is defined as the spectral radius of the matrix $\mathcal{F}_{\text{invB}}\mathcal{V}_{\text{invB}}^{-1}$,

$$\mathcal{R}_{\text{invB}} = \rho(\mathcal{F}_{\text{invB}}\mathcal{V}_{\text{invB}}^{-1}). \qquad (\text{A3})$$

If $\mathcal{R}_{\text{invB}} > 1$, then virus B can invade the system. Equivalent thresholds can also be defined using the next generation matrix approach, such as type or target reproduction numbers [34–36].

We define matrix $\mathcal{F}_{\text{invB}}$ as the rate of new infections,

$$\mathcal{F}_{\text{invB}} = \begin{pmatrix} 0 & \frac{\Phi \bar{X} \alpha_B}{\bar{P}} & \frac{\Phi \bar{X} \epsilon_B \alpha_{AB}}{\bar{P}} \\ \frac{\Phi \theta_B \bar{S}}{\bar{P}} & 0 & \omega_B \\ \frac{\Phi \gamma_B \theta_B i a_{eq}}{\bar{P}} & \frac{\Phi \gamma_A \theta_A z a_{eq}}{\bar{P}} & 0 \end{pmatrix},$$

and matrix $\mathcal{V}_{\text{invB}}$ represents other transitions,

$$\mathcal{V}_{\text{invB}} = \begin{pmatrix} c + \delta_B & 0 & 0 \\ 0 & \mu + \omega + \frac{\Phi \gamma_A \theta_A z a_{eq}}{\bar{P}} & 0 \\ 0 & 0 & \mu + \omega_A + \omega_B \end{pmatrix}.$$
The recovery of virus B from a co-infected plant, $\omega_B$, is counted as a new infection in matrix $\mathcal{F}_{invB}$. Therefore, the invasion reproduction number is the spectral radius of the following matrix,

$$ M_{invB} = \mathcal{F}_{invB} V_{invB}^{-1} = \begin{pmatrix} 0 & \frac{\Phi \bar{X} \alpha_B}{\kappa_1} & \frac{\Phi \bar{X} \epsilon_B \alpha_{AB}}{\bar{P} \kappa_2} \\ \frac{\Phi \theta_B \bar{S}}{\bar{P}(c + \delta_B)} & 0 & \frac{\omega_B}{\kappa_2} \\ \frac{\Phi \gamma_B \theta_B i a_{eq}}{\bar{P}(c + \delta_B)} & \frac{\Phi \gamma_A \theta_A z a_{eq}}{\kappa_1} & 0 \end{pmatrix}, \quad (A4) $$

where $\kappa_1 = (\mu + \omega)\bar{P} + \Phi\gamma_A\theta_A z a_{eq}$ and $\kappa_2 = \mu + \omega_A + \omega_B$. The invasion reproduction number is the largest positive root of the characteristic equation of matrix $M_{invB}$, a cubic equation. Instead of writing the complex expression for the invasion reproduction number, we interpret how the nonzero entries in the invasion matrix $M_{invB}$ affect the invasion reproduction number.

Equilibrium prevalences of virus A in the plant population do not depend on the relative inoculation success of viruses A or B, parameters $\gamma_A$ or $\gamma_B$, but the invasion reproduction curves do depend on these parameters through the two matrix entries in (A4) in the last row. In Figure A1, three invasion reproduction curves are graphed for the relative inoculation success ratios $\gamma_A = 2$ and for $\gamma_B = 0.1, 0.5, 1$. As $\gamma_B$ increases, so do the invasion curves.

Figure A2 illustrates the time course of a successful invasion ($\mathcal{R}_{invB} > 1$) and an unsuccessful invasion ($\mathcal{R}_{invB} < 1$) of virus B, corresponding to the two sets of parameter values indicated by the black circles on the solid invasion curve in Figure 2.

<figure><img src="image_10.png" /><figcaption>**Figure A1.** For the vector-explicit model, virus A equilibrium prevalences for the vector $z_a = za_{eq}/\bar{V}$ and for the host $i_a = ia_{eq}/\bar{P}$ are graphed in (A,C,E) and the corresponding invasion reproduction curves are graphed in (B,D,F) as a function of acquisition and inoculation parameters. The parameter values are the same as in Figure 2 with the exception that in panels B, D and F parameters $\gamma_A = 2$ and $\gamma_B = 0.1$ (solid curve), 0.5 (dashed curve), 1 (dashed-dot curve).</figcaption></figure>
<figure><img src="image_11.png" /><figcaption>**Figure A2.** The time courses of one successful and one unsuccessful invasion of virus B are graphed for the vector-explicit model. The virus A host prevalences $I_A/P$ are graphed in panels (A,C,E) and the virus B host prevalences $I_B/P$ are graphed in panels (B,D,F) for the two sets of parameter values identified by black circles in Figure 2. The initial condition is perturbed slightly from the virus A equilibrium, ($\bar{V} - za_{eq} - 0.1, va_{eq}, 0.1, \bar{P} - ia_{eq} - 30, ia_{eq} - 10, 10, 20$). Parameter values for the graphs in panels (A–F) are the same as in Figure 2A–F with $\gamma_B = 0.9$.</figcaption></figure>

## Appendix C. Reduction to a Vector-Implicit Model

The vector dynamics generally occur on a faster time scale than the plant dynamics. We use this assumption to make a quasi-steady-state approximation for vector abundances [39]. Consider the vector model (1), where the total vector population is $V = X + Z_A + Z_B$. Assume that the vector recovery rates, $\delta_A$ and $\delta_B$, are large (short infective vector life cycle). Summing the differential equations for $X$, $Z_A$ and $Z_B$ leads to the following differential equation for $V$:

$$\frac{dV}{dt} = \Lambda - cV.$$

Thus, $V(t) = \Lambda/c + (V(0) - \Lambda/c)e^{-ct}$. For $V(0) \approx \Lambda/c$, it follows that $V(t) \approx \Lambda/c = \bar{V}$ is approximately constant, an assumption applied in our model. Replace $X$ by $\bar{V} - Z_A - Z_B$ in the differential equations for $Z_A$ and $Z_B$ and divide by $c + \delta_A$ or $c + \delta_B$:

$$\epsilon_1 \frac{dZ_A}{dt} = \frac{\Phi}{(c + \delta_A)P} (\bar{V} - Z_A - Z_B) [\alpha_A I_A + \epsilon_A \alpha_{AB} I_{AB}] - Z_A \\ \approx \frac{\Phi}{(c + \delta_A)P} \bar{V} [\alpha_A I_A + \epsilon_A \alpha_{AB} I_{AB}] - Z_A \approx 0$$

$$\epsilon_2 \frac{dZ_B}{dt} = \frac{\Phi}{(c + \delta_B)P} (\bar{V} - Z_A - Z_B) [\alpha_B I_B + \epsilon_B \alpha_{AB} I_{AB}] - Z_B \\ \approx \frac{\Phi}{(c + \delta_B)P} \bar{V} [\alpha_B I_B + \epsilon_B \alpha_{AB} I_{AB}] - Z_B \approx 0.$$
In particular, $\epsilon_1 = 1/(c + \delta_A) \ll 1$ and $\epsilon_2 = 1/(c + \delta_B) \ll 1$ due to the short life cycle of infective vectors. In addition, we assumed $\bar{V} \ge Z_A + Z_B$. This latter assumption holds if the vector visitation rate is much smaller than the vector death rate plus recovery rate $\Phi \ll (c + \delta_i)$, $i = A, B$. Solving the right-hand sides of the differential equations for $Z_A$ and $Z_B$ leads to

$$ 
\begin{aligned}
Z_A &\approx \frac{\Phi \bar{V}}{(c + \delta_A)P} [\alpha_A I_A + \epsilon_A \alpha_{AB} I_{AB}] \\
Z_B &\approx \frac{\Phi \bar{V}}{(c + \delta_B)P} [\alpha_B I_B + \epsilon_B \alpha_{AB} I_{AB}].
\end{aligned}
\quad (A5) $$

Let $\hat{\theta}_i = \theta_i / (c + \delta_i)$ for $i = A, B$. The vector-implicit model (Equation (6) in the main text) is

$$ \frac{dS}{dt} = \sigma - \frac{\Phi^2 \bar{V}}{P^2} S (\hat{\theta}_A [\alpha_A I_A + \epsilon_A \alpha_{AB} I_{AB}] + \hat{\theta}_B [\alpha_B I_B + \epsilon_B \alpha_{AB} I_{AB}]) - \mu S + \omega (I_A + I_B) $$

$$ \frac{dI_A}{dt} = \frac{\Phi^2 \bar{V}}{P^2} (\hat{\theta}_A S [\alpha_A I_A + \epsilon_A \alpha_{AB} I_{AB}] - \gamma_B \hat{\theta}_B I_A [\alpha_B I_B + \epsilon_B \alpha_{AB} I_{AB}]) - (\mu + \omega) I_A + \omega_A I_{AB} $$

$$ \frac{dI_B}{dt} = \frac{\Phi^2 \bar{V}}{P^2} (\hat{\theta}_B S [\alpha_B I_B + \epsilon_B \alpha_{AB} I_{AB}] - \gamma_A \hat{\theta}_A I_B [\alpha_A I_A + \epsilon_A \alpha_{AB} I_{AB}]) - (\mu + \omega) I_B + \omega_B I_{AB} $$

$$ \frac{dI_{AB}}{dt} = \frac{\Phi^2 \bar{V}}{P^2} (\gamma_B \hat{\theta}_B I_A [\alpha_B I_B + \epsilon_B \alpha_{AB} I_{AB}] + \gamma_A \hat{\theta}_A I_B [\alpha_A I_A + \epsilon_A \alpha_{AB} I_{AB}]) - (\mu + \omega_A + \omega_B) I_{AB} $$

The vector parameters appear in the transmission coefficients. For example, the transmission coefficients from a co-infected plant $I_{AB}$ when $P = \bar{P}$ are

$$ \beta_i = \frac{\Phi^2 \bar{V}}{\bar{P}^2} \hat{\theta}_i \alpha_{AB}, \quad i = A, B. \quad (A6) $$

The basic reproduction number and the invasion matrix for the plant model are calculated by a method similar to that used for the vector-explicit model. The basic reproduction number for the vector-implicit model is the same as in (3). The virus A equilibrium ($\bar{S}, \bar{I}_A$) is distinct from the equilibrium for the vector-explicit model in (4),

$$ \bar{S} = \frac{\bar{P}}{\mathcal{R}_{0A}} \quad \text{and} \quad \bar{I}_A = \bar{P} \left(1 - \frac{1}{\mathcal{R}_{0A}}\right). \quad (A7) $$

The invasion matrix for virus B, given $\mathcal{R}_{0A} > 1$, has the following form (matrix (8) in the main text):

$$ \mathcal{M}_{invB} = \begin{pmatrix} \frac{\Phi^2 \bar{V} \alpha_B \hat{\theta}_B \bar{S}}{\eta} & \frac{\Phi^2 \bar{V} \hat{\theta}_B \alpha_{AB} \epsilon_B \bar{S} + \bar{P}^2 \omega_B}{\bar{P}^2 (\mu + \omega_A + \omega_B)} \\ \frac{\Phi^2 \bar{V} \bar{I}_A [\gamma_A \alpha_A \hat{\theta}_A + \gamma_B \alpha_B \hat{\theta}_B]}{\eta} & \frac{\Phi^2 \bar{V} \gamma_B \hat{\theta}_B \alpha_{AB} \epsilon_B \bar{I}_A}{\bar{P}^2 (\mu + \omega_A + \omega_B)} \end{pmatrix}, $$

where $\eta = \bar{P}^2(\mu + \omega) + \Phi^2 \bar{V} \gamma_A \alpha_A \hat{\theta}_A \bar{I}_A$. The invasion reproduction number $\mathcal{R}_{invB} = \rho(\mathcal{M}_{invB})$ is the largest positive root of the quadratic characteristic equation of matrix $\mathcal{M}_{invB}$. Alternately, a successful invasion can be assessed directly from the trace and determinant of the non-negative matrix $\mathcal{M}_{invB}$ [46]. That is,

$$ \operatorname{trace}(\mathcal{M}_{invB}) - \operatorname{det}(\mathcal{M}_{invB}) > 1 \quad (A8) $$

if and only if $\mathcal{R}_{invB} > 1$.

Condition (A8) for invasion of virus B simplifies in the case that virus B cannot invade alone and there is no plant recovery, i.e., if $\alpha_B = 0$ ($\mathcal{R}_{0B} = 0$) and in the plant population, $\omega = 0 = \omega_A = \omega_B$. Under these assumptions, condition (A8) can be expressed as
$$ \frac{\hat{\theta}_B \epsilon_B \alpha_{AB}}{\alpha_A \hat{\theta}_A} \left( \gamma_B + \frac{\gamma_A}{1 + (\mathcal{R}_{0A} - 1)\gamma_A} \right) (\mathcal{R}_{0A} - 1) > 1. $$

The preceding condition shows that the invasion success of virus B increases if the parameters related to virus B are larger than those related to virus A, $\hat{\theta}_B \epsilon_B \alpha_{AB} > \alpha_A \hat{\theta}_A$, or if $\gamma_B$ is large.

For the parameter values $\alpha_B = 0$ and $\omega = 0 = \omega_A = \omega_B$ and the same parameter values as in Figure A1, the equilibrium prevalences for plants infected with virus A and three invasion reproduction curves are graphed for the vector-implicit model in Figure A3. In Figure A4, the time courses of two successful invasions ($\mathcal{R}_{invB} > 1$) and one unsuccessful invasion ($\mathcal{R}_{invB} < 1$) are graphed for the vector-implicit model and for parameter values as in Figure 3.

If the vector recovery rate $\delta_A$ is increased (infective period of vector is shortened) there is better agreement between the vector-implicit and vector-explicit models. That is, if the vector-explicit model predicts no invasion of virus B so does the vector-implicit model and vice versa. In addition, the values of the host equilibrium prevalences for the two models are in closer agreement for large values of $\delta_A$. We can see this if we fix the composite parameter $\hat{\theta}_A = \theta_A/(c + \delta_A) = K$ =constant in the vector-implicit model but vary $\theta_A$ and $\delta_A$. The virus A prevalence and the invasion reproduction number in the vector-implicit model does not change as $\delta_A$ and $\theta_A$ change but this is not the case for the virus A prevalence and the invasion reproduction number in the vector-explicit model (see Figure A5). The invasion reproduction number in the vector-implicit model predicts species B invasion success $\mathcal{R}_{invB} > 1$ but in the vector-explicit model, a successful invasion occurs only if $\delta_A$ is sufficiently large. In the example in Figure A5, $K = 0.05$, $\gamma_B = 0.9$, and $\alpha_A = 0.3$ with other parameter values as in Table 2. Panels A and B in Figure A5 are the results for the vector-explicit model and panels C and D are for the vector-implicit model. The invasion curve where $\mathcal{R}_{invB} > 1$ in the vector-explicit model predicts a successful invasion only if $\delta_A > 7$.

<figure><img src="image_12.png" /><figcaption>**Figure A3.** For the vector-implicit model, the virus A host equilibrium prevalences $i_a = ia_{eq}/\bar{P}$ are graphed in panels (A,C,E) and the corresponding invasion reproduction curves are graphed in panels (B,D,F) as a function of acquisition and inoculation parameters. Parameter values are as in Figure 3A–F with the exception that in panels B, D and F the parameter values for $\gamma_A = 2$ and for $\gamma_B = 0.1$ (solid curve), 0.5 (dashed curve), 1 (dashed-dot curve).</figcaption></figure>
<figure><img src="image_13.png" /><figcaption>**Figure A4.** The time courses of two successful invasions and one unsuccessful invasion of virus *B* are graphed for the vector-implicit model. The virus *A* host prevalences *I*<sub>A</sub>/P are graphed in panels (**A**,**C**,**E**) and the virus *B* host prevalences *I*<sub>B</sub>/P are graphed in panels (**B**,**D**,**F**) for the three sets of parameter values identified by the black circles in Figure 3. The initial condition is perturbed slightly from the virus *A* equilibrium, (P̄ − 30, Ī<sub>A</sub> − 10, 10, 20). Parameter values are as in Figure 3A–F with γ<sub>B</sub> = 0.9.</figcaption></figure>

<figure><img src="image_14.png" /><figcaption>**Figure A5.** Comparison of the virus *A* equilibrium prevalences and invasion reproduction numbers in the vector-explicit model in (**A**,**B**) to the vector-implicit model in (**C**,**D**) as a function of vector recovery rate *δ*<sub>A</sub>. Parameter values are *θ*̂<sub>A</sub> = *θ*<sub>A</sub>/(*c* + *δ*<sub>A</sub>) = 0.05, *α*<sub>A</sub> = 0.3, *γ*<sub>B</sub> = 0.9, *δ*<sub>A</sub> ∈ [0, 19], *θ*<sub>A</sub> = 0.05(*c* + *δ*<sub>A</sub>) ∈ [0.05, 1] and other parameter values as in Table 2.</figcaption></figure>
References

1. Gray, S.M.; Banerjee, N. Mechanisms of arthropod transmission of plant and animal viruses. *Microbiol. Mol. Biol. Rev.* **1999**, *63*, 128–148. [PubMed]

2. Blanc, S. Virus transmission–getting out and in. In *Viral Transport in Plants*; Springer: Berlin/Heidelberg, Germany, 2007; pp. 1–28.

3. Brault, V.; Uzest, M.; Monsion, B.; Jacquot, E.; Blanc, S. Aphids as transport devices for plant viruses. *C. R. Biol.* **2010**, *333*, 524–538. [CrossRef] [PubMed]

4. Bragard, C.; Caciagli, P.; Lemaire, O.; Lopez-Moya, J.; MacFarlane, S.; Peters, D.; Susi, P.; Torrance, L. Status and prospects of plant virus control through interference with vector transmission. *Annu. Rev. Phytopathol.* **2013**, *51*, 177–201. [CrossRef] [PubMed]

5. Jeger, M.; Van Den Bosch, F.; Madden, L.; Holt, J. A model for analysing plant-virus transmission characteristics and epidemic development. *Math. Med. Biol. A J. IMA* **1998**, *15*, 1–18. [CrossRef]

6. Verbeek, M.; van Bekkum, P.J.; Dullemans, A.M.; van der Vlugt, R.A. Torradoviruses are transmitted in a semi-persistent and stylet-borne manner by three whitefly vectors. *Virus Res.* **2014**, *186*, 55–60. [CrossRef] [PubMed]

7. Zhang, X.S.; Holt, J.; Colvin, J. Mathematical models of host plant infection by helper-dependent virus complexes: Why are helper viruses always avirulent? *Phytopathology* **2000**, *90*, 85–93. [CrossRef]

8. Mahuku, G.; Lockhart, B.E.; Wanjala, B.; Jones, M.W.; Kimunye, J.N.; Stewart, L.R.; Cassone, B.J.; Sevgan, S.; Nyasani, J.O.; Kusia, E.; et al. Maize lethal necrosis (MLN), an emerging threat to maize-based food security in sub-Saharan Africa. *Phytopathology* **2015**, *105*, 956–965. [CrossRef]

9. Legg, J.; Jeremiah, S.; Obiero, H.; Maruthi, M.; Ndyetabula, I.; Okao-Okuja, G.; Bouwmeester, H.; Bigirimana, S.; Tata-Hangy, W.; Gashaka, G.; et al. Comparing the regional epidemiology of the cassava mosaic and cassava brown streak virus pandemics in Africa. *Virus Res.* **2011**, *159*, 161–170. [CrossRef]

10. Seabloom, E.W.; Hosseini, P.R.; Power, A.G.; Borer, E.T. Diversity and composition of viral communities: Coinfection of barley and cereal yellow dwarf viruses in California grasslands. *Am. Nat.* **2009**, *173*, E79–E98. [CrossRef]

11. Susi, H.; Barrès, B.; Vale, P.F.; Laine, A.L. Co-infection alters population dynamics of infectious disease. *Nat. Commun.* **2015**, *6*, 5975. [CrossRef]

12. Cunniffe, N.J.; Koskella, B.; Metcalf, C.J.E.; Parnell, S.; Gottwald, T.R.; Gilligan, C.A. Thirteen challenges in modelling plant diseases. *Epidemics* **2015**, *10*, 6–10. [CrossRef] [PubMed]

13. Escriu, F.; Fraile, A.; García-Arenal, F. The evolution of virulence in a plant virus. *Evolution* **2003**, *57*, 755–765. [CrossRef]

14. Syller, J. Biological and molecular events associated with simultaneous transmission of plant viruses by invertebrate and fungal vectors. *Mol. Plant Pathol.* **2014**, *15*, 417–426. [CrossRef] [PubMed]

15. Moury, B.; Fabre, F.; Senoussi, R. Estimation of the number of virus particles transmitted by an insect vector. *Proc. Natl. Acad. Sci. USA* **2007**, *104*, 17891–17896. [CrossRef] [PubMed]

16. Gallet, R.; Michalakis, Y.; Blanc, S. Vector-transmission of plant viruses and constraints imposed by virus-vector interactions. *Curr. Opin. Virol.* **2018**, *33*, 144–150. [CrossRef] [PubMed]

17. Syller, J.; Grupa, A. Antagonistic within-host interactions between plant viruses: Molecular basis and impact on viral and host fitness. *Mol. Plant Pathol.* **2016**, *17*, 769–782. [CrossRef] [PubMed]

18. Zhang, X.S.; Holt, J.; Colvin, J. Synergism between plant viruses: A mathematical analysis of the epidemiological implications. *Plant Pathol.* **2001**, *50*, 732–746. [CrossRef]

19. Zhang, X.S.; Holt, J. Mathematical models of cross protection in the epidemiology of plant-virus diseases. *Phytopathology* **2001**, *91*, 924–934. [CrossRef]

20. Mascia, T.; Gallitelli, D. Synergies and antagonisms in virus interactions. *Plant Sci.* **2016**, *252*, 176–192. [CrossRef]

21. Srinivasan, R.; Hall, D.G.; Cervantes, F.A.; Alvarez, J.M.; Whitworth, J.L. Strain specificity and simultaneous transmission of closely related strains of a Potyvirus by Myzus persicae. *J. Econ. Entomol.* **2012**, *105*, 783–791. [CrossRef]
22. Carroll, J.; Smith, D.; Gray, S. Preferential acquisition and inoculation of PVYNTN over PVYO in potato by the green peach aphid Myzus persicae (Sulzer). *J. Gen. Virol.* **2016**, *97*, 797. [CrossRef] [PubMed]

23. Mondal, S.; Wenninger, E.J.; Hutchinson, P.J.; Whitworth, J.L.; Shrestha, D.; Eigenbrode, S.D.; Bosque-Pérez, N.A. Comparison of transmission efficiency of various isolates of Potato virus Y among three aphid vectors. *Entomol. Exp. Et Appl.* **2016**, *158*, 258–268. [CrossRef]

24. Lacroix, C.; Seabloom, E.W.; Borer, E.T. Environmental nutrient supply alters prevalence and weakens competitive interactions among coinfecting viruses. *New Phytol.* **2014**, *204*, 424–433. [CrossRef] [PubMed]

25. Salvaudon, L.; De Moraes, C.M.; Mescher, M.C. Outcomes of co-infection by two potyviruses: Implications for the evolution of manipulative strategies. *Proc. R. Soc. B Biol. Sci.* **2013**, *280*, 20122959. [CrossRef]

26. Holt, J.; Chancellor, T. Simulation modelling of the spread of rice tungro virus disease: The potential for management by roguing. *J. Appl. Ecol.* **1996**, *33*, 927–936. [CrossRef]

27. Peñaflor, M.F.G.; Mauck, K.E.; Alves, K.J.; De Moraes, C.M.; Mescher, M.C. Effects of single and mixed infections of Bean pod mottle virus and Soybean mosaic virus on host-plant chemistry and host-vector interactions. *Funct. Ecol.* **2016**, *30*, 1648–1659. [CrossRef]

28. Naidu, R.A.; Maree, H.J.; Burger, J.T. Grapevine leafroll disease and associated viruses: A unique pathosystem. *Annu. Rev. Phytopathol.* **2015**, *53*, 613–634. [CrossRef]

29. Blaisdell, G.K.; Cooper, M.L.; Kuhn, E.J.; Taylor, K.A.; Daane, K.M.; Almeida, R.P. Disease progression of vector-mediated Grapevine leafroll-associated virus 3 infection of mature plants under commercial vineyard conditions. *Eur. J. Plant Pathol.* **2016**, *146*, 105–116. [CrossRef]

30. Untiveros, M.; Fuentes, S.; Salazar, L.F. Synergistic interaction of Sweet potato chlorotic stunt virus (Crinivirus) with carla-, cucumo-, ipomo-, and potyviruses infecting sweet potato. *Plant Dis.* **2007**, *91*, 669–676. [CrossRef]

31. McLeish, M.; Sacristán, S.; Fraile, A.; García-Arenal, F. Coinfection Organizes Epidemiological Networks of Viruses and Hosts and Reveals Hubs of Transmission. *Phytopathology* **2019**, *109*, 1003–1010. [CrossRef]

32. Heffernan, J.M.; Smith, R.J.; Wahl, L.M. Perspectives on the basic reproductive ratio. *J. R. Soc. Interface* **2005**, *2*, 281–293. [CrossRef] [PubMed]

33. Diekmann, O.; Heesterbeek, J.A.P.; Metz, J.A. On the definition and the computation of the basic reproduction ratio R 0 in models for infectious diseases in heterogeneous populations. *J. Math. Biol.* **1990**, *28*, 365–382. [CrossRef] [PubMed]

34. Roberts, M.; Heesterbeek, J. A new method for estimating the effort required to control an infectious disease. *Proc. R. Soc. Lond. Ser. B Biol. Sci.* **2003**, *270*, 1359–1364. [CrossRef] [PubMed]

35. van den Driessche, P. Reproduction numbers of infectious disease models. *Infect. Dis. Model.* **2017**, *2*, 288–303. [CrossRef] [PubMed]

36. van den Driessche, P.; Watmough, J. Reproduction numbers and sub-threshold endemic equilibria for compartmental models of disease transmission. *Math. Biosci.* **2002**, *180*, 29–48. [CrossRef]

37. Ortega, J.M. *Matrix Theory: A Second Course*; Plenum Press: New York, NY, USA; London, UK, 1987.

38. Hilker, F.M.; Allen, L.J.; Bokil, V.A.; Briggs, C.J.; Feng, Z.; Garrett, K.A.; Gross, L.J.; Hamelin, F.M.; Jeger, M.J.; Manore, C.A.; et al. Modeling virus coinfection to inform management of maize lethal necrosis in Kenya. *Phytopathology* **2017**, *107*, 1095–1108. [CrossRef]

39. Keeling, M.J.; Rohani, P. *Modeling Infectious Diseases in Humans and Animals*; Princeton University Press: Princeton, NJ, USA, 2008.

40. Hamelin, F.M.; Allen, L.J.; Prendeville, H.R.; Hajimorad, M.R.; Jeger, M.J. The evolution of plant virus transmission pathways. *J. Theor. Biol.* **2016**, *396*, 75–89. [CrossRef]

41. Donnelly, R.; Cunniffe, N.J.; Carr, J.P.; Gilligan, C.A. Pathogenic modification of plants enhances long-distance dispersal of nonpersistently transmitted viruses to new hosts. *Ecology* **2019**, e02725. [CrossRef]

42. Gildow, F.E.; Rochow, W.F. Transmission interference between two isolates of barley yellow dwarf virus in Macrosiphum avenae. *Phytopathology* **1980**, *70*, 122–126. [CrossRef]

43. Cohen, S.; Duffus, J.E.; Liu, H.Y. Acquisition, interference, and retention of cucurbit leaf curl viruses in whiteflies. *Phytopathology* **1989**, *79*, 109–113. [CrossRef]

44. Ng, J.C.; Chen, A.Y. Acquisition of Lettuce infectious yellows virus by Bemisia tabaci perturbs the transmission of Lettuce chlorosis virus. *Virus Res.* **2011**, *156*, 64–71. [CrossRef] [PubMed]
45. Seabloom, E.W.; Borer, E.T.; Gross, K.; Kendig, A.E.; Lacroix, C.; Mitchell, C.E.; Mordecai, E.A.; Power, A.G. The community ecology of pathogens: Coinfection, coexistence and community composition. *Ecol. Lett.* **2015**, *18*, 401–415. [CrossRef] [PubMed]

46. Allen, L.J. *Introduction to Mathematical Biology*; Pearson/Prentice Hall: Upper Saddle River, NJ, USA, 2007.

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

© 2019 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).