Towards a thermodynamic mechanistic model for the effect of temperature on arthropod vector competence for transmission of arboviruses.

Arboviruses such as West Nile virus (WNV), bluetongue virus (BTV), dengue virus (DENV) and chikungunya virus (CHIKV) infect their arthropod vectors over a range of average temperatures depending on the ambient temperature. How the transmission efficiency of an arbovirus (i.e. vector competence) varies with temperature influences not only the short term risk of arbovirus outbreaks in humans and livestock but also the long term impact of climate change on the geographical range of the virus. The strength of the interaction between viral surface (glyco)protein (GP) and the host cell receptor (Cr) on binding of virus to host cell is defined by the thermodynamic dissociation constant Kd_receptor which is assumed to equal 10-3 M (at 37 °C) for binding of a sialic acid (SA) on the arthropod midgut epithelial cell surface to a SA-binding site on the surface of BTV, for example. Here virus binding affinity is modelled with increasing number of GP/Cr contacts at temperatures from 10 °C to 35 °C taking into account the change in entropy on immobilization of the whole virus on binding (ΔSa_immob). Based on published data, three thermodynamic GP/Cr binding scenarios, namely enthalpy-driven, entropy-assisted and entropy-driven, are shown to affect the temperature sensitivity of virus binding in different ways. Thus for enthalpy-driven GP/Cr binding, viruses bind host cells much more strongly at 10 °C than 35 °C. A mechanistic model is developed for the number of arthropod midgut cells with bound virus and by building in a kinetic component for the rate of arbovirus replication and subsequent spread to the arthropod salivary glands, a model for the effect of temperature on vector competence is developed. The model separates the opposing effects of temperature on midgut cell binding affinity from the kinetic component of virogenesis. It successfully accommodates both increases in vector competence with temperature as for DENV and WNV in mosquitoes and decreases as for the CHIKV 2010-1909 strain in various populations of Aedes albopictus mosquitoes. Enhanced cell binding at lower temperatures through enthalpy-driven GP/Cr binding compensates for the lower replication rate to some degree such that some transmission can still occur at lower temperatures. In contrast, the strength of entropy-driven GP/Cr binding diminishes at low temperatures although there is no minimum temperature threshold for transmission efficiency. The magnitude of ΔSa_immob is an important data gap. It is concluded that thermodynamic and kinetic data obtained at the molecular level will prove important in modelling vector competence with temperature.

Introduction

Arboviruses are so named because they are arthropod-borne and cycle between an arthropod vector (e.g. mosquito, biting midge, tick) and a vertebrate host. Arboviruses may cause disease in humans or livestock or both and examples include West Nile virus (WNV), bluetongue virus (BTV), dengue virus (DENV) and chikungunya virus (CHIKV). Temperature affects arbovirus transmission through the geographic range and global distribution of the arthropod vector species and through the local density of arthropods and the rates at which the arthropods feed on and bite vertebrate hosts (Guis et al., 2012, Wittmann and Baylis, 2000, Carpenter et al., 2011, Turner et al., 2013). Temperature is also a major factor in arthropod vector competence which is the ability of an arthropod species or sub-population of species to transmit an arbovirus to a subsequent host individual through feeding, after itself being infected through a blood-meal from a previous viraemic vertebrate host (Mullens et al., 2004, Wittmann et al., 2002). Arthropods are poikilothermic and arboviruses thus have to infect the arthropod vector at the ambient temperature, the average of which may vary considerably depending on climate, season and local conditions (micro-climate). In this respect, understanding how vector competence varies with temperature is important for assessing not only the short term risk of a BTV epizootic in livestock for example but also the long term impact of climate change on the geographical range of an arbovirus (Guis et al., 2012, Wittmann and Baylis, 2000). The vector competence considered here is restricted to a given life cycle stage (namely the adult midge or mosquito) of the arthropod and does not include transstadial or transovarial transmission of the virus within the arthropod.

Vector competence is quantified as arbovirus transmission efficiency. In general arthropod vector competence increases with temperature as in the case of BTV in Culicoides biting midges (Wittmann and Baylis, 2000, Carpenter et al., 2011) and flaviviruses such as DENV and WNV in mosquitoes (Vogels et al., 2016, Liu et al., 2017). However, there are also a few specific arbovirus/vector combination examples for which the opposite effect has been observed (Samuel et al., 2016) with enhanced vector competence for arbovirus being associated with lower temperatures. For example, in a study with CHIKV in Ae. albopictus, Zouache et al. (2014) reared all the mosquito larval stages at 26 °C before testing the competence of the adult females to virus at 20 °C and 28 °C; 20 °C representing the low-temperature threshold recorded when local CHIKV transmission occurred during the Italian and French epidemics and 28 °C the typical average temperature in tropical regions. The 2010–1909 CHIKV strain in general showed higher transmission efficiencies in Brazilian Ae. albopictus mosquito populations at 20 °C compared to 28 °C while the 06–021 CHIKV strain in general showed increased competence at the higher temperature (Zouache et al., 2014).

The effect of temperature on vector competence is fundamentally determined by molecular events involving proteins and other macromolecules of the arbovirus and the arthropod host during host cell infection. In physical biochemistry, temperature has two independent effects. First, it affects the strength of the binding affinity between a virus and its host cell through the thermodynamic equilibrium constant Ka (Gale, 2017) according to Equation 1 (see Table 1 ) and second it affects the rate of biochemical reaction according to the Arrhenius equation (Equation 2). The overall objective of the work here was to develop a mechanistic model for the effect of temperature on arthropod vector competence by linking Equation 1 and Equation 2. The main parameters of the model, in addition to temperature, would therefore be the thermodynamic terms, namely the changes in the enthalpy (ΔHa_virus) and entropy (ΔSa_virus) on virus binding to the arthropod midgut epithelial cells in Equation 1 and the kinetic terms, namely the activation energy (EA) and a “rate” constant (pcomplete283) in Equation 2.

Table 1

Equations used. See methods for derivation.

Equation 1	K_a=e^−ΔH_a_virusRT+ΔS_a_virusR
Equation 2	p_completeT=p_complete283×e^(E_AR×(1283−1T))
Equation 3	ΔH_a_virus=n×ΔH_a_receptor
Equation 4	ΔS_a_virus=ΔS_a_immob+n×ΔS_a_receptor
Equation 5	F_c=[V_free][V_free]+1/K_a
Equation 6	C.V_T=F_c×C_total_midgut
Equation 7	p_transmissionT=1−(1−p_completeT)^C.V_T
Equation 8	ΔG_a_receptor=RTlnK_d_receptor
Equation 9	ΔG_a_receptor=ΔH_a_receptor−TΔS_a_receptor
Equation 10	ΔS_a_receptor=ΔS_solvent+ΔSconfGP+ΔSconfCr
Equation 11	ΔH_a_receptor=TΔS_a_receptor+RTlnK_d_receptor
Equation 12	C_free+V_free↔C.V
Equation 13	K_a=[C.V][V_free][C_free]
Equation 14	ΔS_a_immob=ΔS_trans+ΔS_rot+ΔS_mem
Equation 15	K_a_T2K_a_T1=e^−nΔH_a_receptorR(1T_2−1T_1)
Equation 16	K_a=e^−nΔH_a_receptorRT×e^(ΔS_a_immob+nΔS_a_receptor)R
Equation 17	K_a=1(K_d_receptor)^n×e^ΔS_a_immobR
Equation 18	lnK_a=−nΔG_a_receptorRT+ΔS_a_immobR

The advantage of this physical biochemistry-based approach is that three of these parameters, namely ΔHa_virus and ΔSa_virus (or some of their component terms according to Equation 3 and Equation 4 respectively) and EA can be measured in biochemical “test tube” experiments (Carneiro et al., 2002, Fei et al., 2015, Scholtessik and Rott, 1969). In theory therefore it may be possible to partly parameterize a model for the effect of temperature on vector competence based on biochemical data obtained from molecular components isolated from the arbovirus and arthropod host cells. Furthermore, the model developed here identifies the data gaps, namely the pcomplete283 “rate” term in Equation 2 and of great interest, the change in entropy on immobilization of the whole virus on its binding to the host cell, ΔSa_immob, in Equation 4. To the author's knowledge the ΔSa_immob term identified here through application of a thermodynamic analysis to virus binding has never been defined before let alone quantified. Its importance is therefore unknown but is potentially huge as discussed here. It is anticipated that through application of the model developed here to experimental virus binding data and to empirical vector competence data it may be possible to back calculate ΔSa_immob and pcomplete283. However, even if the mechanistic model developed here is never sufficient to replace insectary-based vector competence experiments measuring infection rates in live arthropods, it offers a novel approach that complements current statistical approaches (Carpenter et al., 2011) to relate experimental vector competence data to temperature by providing a thermodynamic/kinetic framework that naturally accommodates temperature through Equation 1 and Equation 2, the parameters of which are directly related to the molecular mechanisms of the various infection stages. Furthermore changing these parameters in the model e.g. on the basis of molecular structure data (Gale, 2018), could aid the interpretation and understanding of experimental vector competence data (Wittmann et al., 2002) in terms of molecular mechanisms, highlighting the direction of further experiments at the molecular level. For example, it is shown here that the thermodynamic nature of virus binding could determine whether significant transmission occurs at low temperatures or whether instead there is greatly diminished transmission. This is important for assessing the seasonal risks of BTV transmission in northern Europe, for example. The model here could also be used to develop a dose-response for infection of an arthropod at constant temperature although the probability of getting through the pathogen-recognition receptors of the innate immune system in the arthropod midgut would need to be modelled, perhaps using the thermodynamic approach developed for virus binding to mucins in the human intestinal tract model (Gale, 2018). Although not discussed here, this could provide useful information for arbovirus transmission in an arthropod at low challenge doses and hence separate “dead end” vertebrate host species from potential reservoir vertebrate host species on the basis of measured viral titres in their blood and hence in the blood-meal ingested by an individual arthropod. Insectary experiments to measure low transmission efficiencies at low virus challenge doses would require individual analysis of huge numbers of arthropods. The model could be used to model transmission efficiency in terms of two independent variables, namely arbovirus dose and vector temperature. A unique advantage of the model developed here is that the predicted differences in virus infectivity and vector competence are due to differences in the input parameters ΔHa_virus, ΔSa_virus, EA and pcomplete283 which in the arbovirus/arthropod infection process will be related to the gene sequences of specific viral and host proteins (see Gale, 2018). Thus this thermodynamic treatment may provide the link to relate genetic sequence data to virus infectivity and vector competence which is a major goal of current research.

The steps in the mechanistic model developed here for infection of an arthropod vector are shown in Fig. 1 . The thermodynamic component of the model here focuses on arbovirus binding to the microvillae surface (the so-called brush borders) of the arthropod midgut epithelial cell which is the first step in infection of the arthropod (Forrester et al., 2012, Franz et al., 2015) and involves specific interactions between a receptor molecule (Cr) on the brush border surface and a (glyco)protein (GP) on the arbovirus coat. In the case of BTV for example, the Cr may be a sialic acid (SA) on the arthropod midgut epithelial cell surface (Zhang et al., 2010). Here three thermodynamic GP/Cr binding scenarios, namely enthalpy-driven (Fei et al., 2015), entropy-driven (Casasnovas and Springer, 1995) and entropy-assisted (Caro et al., 2017) are considered in terms of their impact in modulating the effect of temperature on the strength of virus binding to a midgut cell as represented by Ka. In enthalpy-driven binding, it is the large negative ΔHa_virus term in Equation 1 that dominates in maximizing Ka and overcoming an unfavorable negative ΔSa_virus term (see Table 2 ), while in entropy-driven binding it is the large positive ΔSa_virus term that overcomes an unfavorable positive ΔHa_virus term. In entropy-assisted binding a small positive ΔSa_virus term adds to a small negative ΔHa_virus term in maximizing Ka (Table 2). The magnitude of Ka together with the virus dose (∼Vfree) in the midgut dictate the number (C.VT) of arthropod midgut cells with bound virus at a given temperature (Gale, 2018) according to Equation 5 and Equation 6. In addition to the strength of the binding of GP to Cr as represented by the thermodynamic dissociation constant, Kd_receptor (Xiong et al., 2013), two other factors affect Ka, namely the number, n, of GP/Cr contacts made on binding of the virus to a cell and the unknown entropy term ΔSa_immob (see above). The kinetic component as represented by pcompleteT in Equation 2 is the probability of spread of infectivity through the arthropod to its salivary glands at a given temperature once virus in the initial challenge dose in the ingested blood-meal has bound to a midgut cell brush border. This represents the rate of virus replication and has to be completed within the lifetime of the arthropod so that it can infect another vertebrate host through feeding before it dies. In this respect, the probability pcompleteT has an implicit time dimension. The thermodynamic output, C.VT, and kinetic output, pcompleteT, are linked through the standard “dose-response” relationship in Equation 7 to give the transmission efficiency (vector competence) at temperature T as represented by ptransmissionT.

Fig. 1

Mechanistic model for infection of an arthropod vector by an arbovirus at temperature T. The output probability, ptransmissionT, of successful infection of the arthropod salivary glands given oral exposure is calculated by combining the number of midgut cells with bound virus (C.VT) and pcompleteT through Equation 7. Stage 1 is the thermodynamic binding component and Stage 2 is kinetic.

Table 2

Published thermodynamic parameters for GP/Cr interactions.

GP/Cr system studied	ΔHa_receptor (kJ/mol)	ΔSa_receptor (J/mol/K)	Kd_receptor (M)	Thermodynamic scenario
*AIV (HA) binding to soluble α2–3 sialyllactose (Fei et al., 2015)	−173.2	−550.4	2.6 × 10−2 (25 °C)	Enthalpy-driven
*AIV (HA) binding to soluble α2–6 sialyllactose (Fei et al., 2015)	−134.7	−424.0	3.4 × 10−2 (25 °C)	Enthalpy-driven
*Low affinity HRV3 binding to soluble ICAM-1 (Casasnovas and Sringer, 1995)	+26.4	+184.1	1.2 × 10−5 (20 °C)	Entropy-driven
*High affinity HRV3 binding to soluble ICAM-1 (Casasnovas and Sringer, 1995)	+14.6	+167.4	0.73 × 10−6 (20 °C)	Entropy-driven
Protein/protein interaction between calmodulin and neuronal nitric oxide synthase (Caro et al., 2017)	−39.5	+15.6	2.8 × 10−8 (35 °C)	Entropy-assisted

Although these studies used immobilized whole virus which may present multiple GPs, the Cr used in each experiment was soluble and thus the binding measured can be interpreted as for a single GP/Cr interaction and results therefore represent Kd_receptor. Thus although AIV HA is trimeric, the Kd_receptor represents the HA monomer. Values presented for ΔSa_receptor may be too negative for use in the model here because Cr is soluble (see text).

The specific objectives are:

To derive realistic thermodynamic parameters to represent enthalpy-driven, entropy-assisted and entropy-driven GP/Cr binding scenarios each with a Kd_receptor of 10−3 M (37 °C) for use in the vector competence model. This is important to understand how different GP/Cr binding scenarios may control or modulate the effect of temperature on Ka and hence vector competence;

To calculate the biologically significant range of Ka for arbovirus binding in the arthropod midgut model. This is important to understand when C.VT will not be affected by changes in temperature despite changing Ka;

To understand how the magnitude of Ka depends on the temperature, the number of GP/Cr contacts for the three thermodynamic GP/Cr binding scenarios and the magnitude of ΔSa_immob. This provides a basis for understanding the outputs from the vector competence model;

To develop, using the best available data, a two-stage mechanistic model based on virus binding affinity and virus replication rate which can accommodate both the increases and decreases in vector competence with temperature observed experimentally; and

To identify mechanisms by which the effect of temperature on virus binding may be modulated and hence explain the temperature paradox. The temperature paradox is that according to thermodynamics, strong enthalpy-driven virus binding should be much weaker at higher temperatures.

The Methods, Results and Discussion sections of this paper are structured and numbered according to these five specific objectives.

Methods

The outputs produced here using the various equations in Table 1 are referred to as models because they allow the effect on Ka or vector competence of changing the variables n and temperature together with the thermodynamic/kinetic parameters to be assessed. The model outputs were calculated with Microsoft Excel. The basic thermodynamic and kinetic equations can be obtained from most advanced Biochemistry or Physical Chemistry text books.

Nomenclature: The terminology for the “receptor” and “virus” subscripts used in the thermodynamic parameters (e.g. ΔHa_receptor and ΔHa_virus) is adopted from Xiong et al. (2013) where “receptor” represents the interaction between GP and Cr and “virus” represents the interaction between whole virus and host cell. Thus for example, although Equation 3 does not specifically refer to GP and host cell, ΔHa_virus would formally be defined as the difference in enthalpy between a host cell with a bound whole virus (product) and the sum of the enthalpies of a free whole virus and a free host cell (reactants) while ΔHa_receptor is defined as the difference in enthalpy between the GP.Cr complex (product) and the sum of the enthalpies of an unbound GP and an unbound Cr (reactants). In thermodynamics it is the differences in enthalpies as represented by ΔHa_receptor and ΔHa_virus that are quantified and used rather than the absolute enthalpies of the “reactants” and “products” which are generally unknown. The same applies to the entropy and the Gibbs free energy.

To derive realistic thermodynamic parameters to represent enthalpy-driven, entropy-assisted and entropy-driven GP/Cr binding scenarios each with a Kd_receptor of 10−3 M (37 °C) for use in the vector competence model

These parameters allow the model to test how different GP/Cr binding scenarios may control or modulate the effect of temperature on Ka and hence vector competence.

Strength of binding of GP to Cr

The strength of the binding affinity between the GP and Cr molecules is quantified by the dissociation constant, Kd_receptor (Xiong et al., 2013) which is related to the change in the Gibbs free energy (ΔGa_receptor) for association of GP and Cr according to Equation 8. The magnitude of ΔGa_receptor is related to the changes in enthalpy (ΔHa_receptor) and entropy (ΔSa_receptor) in the GP and Cr proteins on molecular interaction during binding according to Equation 9 (Gale, 2017, Gale, 2018). The ΔSa_receptor term may be broken down into a number of components which act additively (Bostrom et al., 2011, Mammen et al., 1998, Caro et al., 2017, Gale, 2018) according to Equation 10. These components are the change in entropy of the water solvent molecules on the GP and Cr surfaces (ΔSsolvent) and the changes in the internal conformational freedoms within the GP and Cr proteins as represented by and respectively. Any rotational and translational entropy components are ignored here because GP and Cr are already immobilized in the virus and host cell membranes, respectively.

Available thermodynamic data for binding of GP to Cr

While there are many studies reporting Kd_receptor values, there are relatively few papers detailing the thermodynamic parameters (i.e. ΔHa_receptor and ΔSa_receptor) for GP/Cr interactions. The ΔHa_receptor and ΔSa_receptor data for avian influenza virus (AIV) haemagglutinin (HA) monomers binding to soluble SA glycans (sialyllactose) (Fei et al., 2015) and for human rhinovirus serotype 3 (HRV3) binding to soluble intercellular adhesion molecule-1 (ICAM-1) (Casasnovas and Springer, 1995) are presented in Table 2 and are examples of enthalpy-driven and entropy-driven GP/Cr binding respectively. However, the Kd_receptor values for HA monomer/sialyllactose binding at 25 °C (Table 2) are much higher than the 10−3 M generally associated with HA monomer/SA binding (Taylor and Drickamer, 2006, Xiong et al., 2013). This may be because the sialyllactose (Cr) used in Fei et al. (2015) being soluble has a high degree of motional freedom such that the ΔSa_receptor values observed by Fei et al. (2015) (Table 2) are too negative to represent SA glycans in Cr immobilized on the host cell surface and thus the Kd_receptor values in Table 2 may be unrealistically high for the virus/host cell interaction. Indeed Caro et al. (2017) reported values of ΔSa_receptor for protein/ligand binding ranging from −243 J/mol/K to +263 J/mol/K reflecting negative and positive values of and (Equation 10). For the model here it is desired to have realistic values of ΔSa_receptor and values of −125 J/mol/K, +40 J/mol/K and +125 J/mol/K are therefore chosen for the enthalpy-driven, entropy-assisted and entropy-driven GP/Cr binding scenarios respectively, (Table 3 ) as these are well within the range reported by Caro et al. (2017). Rearranging Equation 8 and Equation 9, the magnitude of ΔHa_receptor required to achieve the specified Kd_receptor value of 10−3 M at 37 °C given an assumed or known value of ΔSa_receptor is calculated using Equation 11.

Table 3

Thermodynamic parameters to achieve Kd_receptor of 10−3 M (37 °C) for enthalpy-driven, entropy-assisted and entropy-driven GP/Cr binding scenarios as used in models here for arbovirus attaching to arthropod midgut brush borders.

Scenario	ΔHa_receptor (kJ/mol)a	ΔSa_receptor (J/mol/K)b	Kd_receptor (M) (37 °C)	Molecular description
Enthalpy-driven	−56.545	−125	10−3	Good fit with immobilization
Entropy-assisted	−5.395	+40	10−3	Weaker fit with some increase in mobility
Entropy-driven	+20.954	+125	10−3	Unfavourable fit overcome with large increase in mobility

The ΔHa_receptor value required to give a Kd_receptor of 10−3 M in a single GP/Cr interaction at 37 °C (Xiong et al., 2013, Taylor and Drickamer, 2006) is calculated from ΔSa_receptor according to Equation 11.

See text for origin.

To calculate the biologically significant range of Ka for arbovirus binding in the arthropod midgut model

This is important to explain why the number of midgut cells with bound virus (C.VT) may not be affected by changes in temperature despite changing Ka. In a prototype mechanistic dose-response model for infection of a human intestine (Gale 2018), it was shown that increasing Ka above 1015 M−1 had no effect on Fc (Equation 5) while for Ka < 105 M−1, unrealistically high virus challenge doses are required to give a reasonable probability of infection. The main conclusion was that the value of Ka over the range 105 to 1015 M−1 is critical in affecting Fc and hence the risk of infection. Here, Equation 5 is used to calculate the Ka values for which Fc equals 0.001, 0.5, and 0.999 representing 1 midgut cell with bound virus, half the midgut cells with bound virus and nearly all the midgut cells with bound virus, respectively (see below for description of midgut model).

To understand how the magnitude of Ka depends on the temperature, the number of GP/Cr contacts for the three thermodynamic GP/Cr binding scenarios and the magnitude of ΔSa_immob

This provides a basis for understanding the outputs from the vector competence model. The central role of Ka is apparent from Fig. 1.

Strength of binding of virus to cell as represented by Ka

Binding of the virus to its host cell may involve a multiple number (n) of GP/Cr contacts depending on the geometry of the virus (Xiong et al., 2013) and the mobility of the Cr molecules within the host cell membrane (Kuhmann et al., 2000, Harada et al., 2004). The binding of a virus to a host cell is a dynamic equilibrium (Handel et al., 2014, Gale, 2017, Gale, 2018) such that within the given volume of the human intestine, lung airway fluid, body fluid or arthropod midgut where infection takes place, free virus (Vfree) binds reversibly to a host cell with no bound virus (Cfree) to give a cell with bound virus (C.V) as represented by Equation 12. The strength of interaction is quantified by the association constant, Ka, which can be expressed both in terms of the concentrations in moles dm−3 (M) (represented by square brackets, []) of C.V, Vfree and Cfree in Equation 13 (Gale, 2017, Gale, 2018) and also in terms of the changes in enthalpy and entropy, ΔHa_virus and ΔSa_virus respectively, of the virus/host cell system on association and the temperature T in Equation 1 (Gale, 2018). The larger Ka is in magnitude, the stronger the binding.

Since enthalpy is a state function, the total enthalpy change on binding of a virus to a cell, ΔHa_virus, is the sum of the enthalpy changes for the n GP/Cr contacts made during binding (Equation 3). Entropy is also a state function, and the same applies for the GP/Cr entropy term (ΔSa_receptor) in Equation 4. In addition, the change in entropy, ΔSa_immob, of the whole virus as it is immobilized on the host cell surface has to be added in separately. Thus, overall the change in entropy on virus binding to a host cell, ΔSa_virus, can be separated into the whole virus term and the GP/Cr receptor term according to Equation 4. The ΔSa_immob term is itself comprised of several entropy terms (Equation 14). These include the changes in translational and rotational entropy of the whole virus (ΔStrans and ΔSrot respectively) and an entropic pressure (ΔSmem) associated with bringing two membranes close together, as for example, when the virus envelope approaches the host cell membrane (Sharma, 2013).

Data for ΔSa_immob: Changes in rotational and translational entropy on binding of whole virus to host cells

To the author's knowledge Equation 4 represents the first time ΔSa_immob has been formally considered for virus binding and not surprisingly there are no data. Therefore, some assumptions are made here to enable the significance of the magnitude of ΔSa_immob to be assessed. Both ΔStrans and ΔSrot are negative for soluble protein-protein associations and the ΔStrans + ΔSrot term for binding of the enzyme acetylcholinesterase to its snake toxin inhibitor, fasciculin-2, is ∼−125 J/mol/K (Minh et al., 2005). ΔStrans and ΔSrot and hence ΔSa_immob would be expected to be much more negative for binding of a whole virus to a host cell. Thus, according to Mammen et al. (1998) the translational entropy of a molecule is related to the natural logarithm of its mass. On the basis of the molecular weights of Bromegrass mosaic virus at 4.6 × 106 Da (Bockstahler and Kaesberg, 1962) and Semliki Forest virus at 6.5 × 107 Da (Houk et al., 1990) ΔStrans would be 6.5-fold and 9-fold higher, respectively, than that for fasciculin-2 at 6735 Da (Sigma-Aldrich, 2018). For the purpose of demonstration of the effect of ΔSa_immob here, ΔSa_immob is assumed to be 6-fold that for fasciculin-2 binding (Minh et al., 2005), i.e. 6 x −125 J/mol/K = −750 J/mol/K.

Data on number of GP/Cr contacts made during virus/host cell binding

There are multiple Cr molecules on the cell surface of each susceptible host cell (Houk et al., 1990, Nunes-Correia et al., 1999, Schlegel et al., 1982) allowing the possibility of a single virion making multiple GP/Cr contacts on binding depending on their spatial proximity. The number of GP/Cr contacts per virion also depends on the number of GP molecules per virion and their spatial geometry. Xiong et al. (2013) estimate that 5.5 AIV HA monomers can bind to sugars from the geometry of the virus particles. Attachment of bluetongue virus (BTV) is mediated by the VP2 coat protein (Zhang et al., 2010) which is a trimer and presents three SA-binding pockets for interaction with SA glycans on the host cell surface (Zhang et al., 2010). Cell binding of BTV could involve more than one VP2 protein as there are 60 VP2 proteins per BTV virion (Zhang et al., 2010) and thus more than 3 GP/Cr contacts could be possible. Kuhmann et al. (2000) demonstrated that approximately four to six Cr molecules assemble around the HIV virion to form a complex needed for infection. Carneiro et al. (2002) estimated 7 GP/Cr contacts for vesicular stomatitis virus (VSV) binding, albeit to phospholipid bilayers. Overall it is concluded here that the number of GP/Cr contacts could range from 3 to 7 for spherical viruses. However, some viruses such as AIV and Zaire ebolavirus (EBOV) form filaments (Seladi-Schulman et al., 2013, Dadonaite et al., 2016) which could make many more GP/Cr contacts during infection than spherical forms. Thus here up to 20 GP/Cr contacts are modelled for the purpose of demonstration.

Modelling the effect of temperature on Ka taking into account the number of molecular contacts

The binding affinity of viruses to host cells is predicted over a range of temperatures using Equation 1. The change in Ka with temperature is related to the magnitudes of ΔHa_receptor and n according to the Van't Hoff Isochore (Gale, 2017) which on integrating with respect to temperature gives the relative change in Ka over the temperature range T1 to T2 (Equation 15).

Modelling the effect of n and ΔSa_immob on Ka at constant temperature for each of the three thermodynamic GP/Cr binding scenarios

Substituting ΔHa_virus and ΔSa_virus in Equation 1 with Equation 3 and Equation 4 respectively allows Ka to be calculated in terms of the number, n, of GP/Cr contacts (Equation 16) together with the thermodynamic parameters for GP/Cr binding and ΔS. Equation 16 can be rearranged into Equation 17 and Equation 18 which express Ka in terms of Kd_receptor and ΔGa_receptor, respectively, in addition to n and ΔSa_immob, thus simplifying interpretation.

To develop, using the best available data, a two-stage mechanistic model based on virus binding affinity and virus replication rate which can accommodate both the increases and decreases in vector competence with temperature observed experimentally

Overview of the model

Forrester et al. (2012) define three main infectious transitions or stages during mosquito infection after taking a viraemic blood-meal from a vertebrate host: 1) midgut epithelial cells infection, 2) midgut escape, when the virus enters the haemocoel to infect secondary tissues, and 3) salivary gland infection resulting in release of virus with saliva for horizontal transmission to an uninfected vertebrate host (Smith et al., 2008, Franz et al., 2015). The spread of arbovirus through the midgut and salivary glands of a mosquito is shown schematically in Franz et al. (2015). For simplicity and to demonstrate proof of principle, the mechanistic model here (Fig. 1) is based on just two stages namely: 1) binding of virus in the blood-meal to the midgut epithelial cell brush borders through GP/Cr interactions, and 2) subsequent infection of the midgut cell followed by midgut escape and salivary gland infection. The first is the thermodynamic component, the second is the kinetic component. The output of the model is the probability of transmission, ptransmissionT, of an arbovirus by an arthropod at temperature T. This is defined as the probability of completion of infection of the arthropod salivary glands within the lifetime of the arthropod given the arthropod has ingested a blood-meal containing arbovirus. The time between initial oral challenge of the arthropod and the ability to infect a new vertebrate host is called the extrinsic incubation period (EIP) and decreases with increasing temperature for orbiviruses in biting midges (Carpenter et al., 2011, Mullens et al., 2004, Wittmann et al., 2002) and mosquitoes (Liu et al., 2017), as does the lifetime of the arthropod (Wittmann and Baylis, 2000, Turner et al., 2013). In effect the transmission probability is the probability that the EIP is less than the lifetime of the arthropod. In the case of Venezuelan equine encephalitis virus (VEEV) in Aedes taeniorhynchus mosquitoes, a single infected midgut epithelial cell may sometimes be sufficient to establish infection of the whole mosquito including its salivary glands (Smith et al., 2008). Indeed Smith et al. (2008) describe probabilistic events where ∼23 infected midgut cells are needed to obtain infection in ∼100% of mosquitoes, and ∼1 infected midgut cell has a 50% probability of leading to infection of the whole mosquito. The model here builds on this probabilistic concept for the infected midgut cell by going back to the stage of initial binding of virus to the midgut epithelial cell surface and defining pcompleteT as the probability that, given a virion has bound to the surface of a midgut cell, that midgut cell becomes infected and its progeny viruses go on to infect the salivary gland so completing the arthropod infection process within the life time of the arthropod at temperature T. The model output, ptransmissionT, is related to the number of midgut cells (C.VT) with bound virus at temperature T according to Equation 7 and is the probability of at least one of those cells successfully infecting the rest of the arthropod at temperature T within its lifetime. Thus the more midgut cells (C.VT) with bound virus then the greater the chance that infection will reach the salivary glands from at least one of those midgut cells within the lifetime of the arthropod.

Thermodynamic component of Equation 7: Modelling C.VT and the binding of arbovirus to epithelial cells in the arthropod midgut

This is Stage 1 in Fig. 1. The thermodynamic equilibrium component is based on the mechanistic dose-response model developed for faecal-oral viruses in a human intestinal tract (Gale, 2018) and is applied to the midgut of the arthropod. The approach requires an estimate of the arthropod midgut volume, the number of epithelial cells lining the midgut (to which virus can bind) and the challenge dose in virions. The volume of blood ingested by an arthropod vector is related to the midgut volume and varies depending on the life cycle stage and the species of arthropod. Thus an adult female Culicoides variipennis biting midge is able to ingest ∼10−4 ml (10−10 m3) (Fu, 1995) and larval tick stages ingest ∼1 µl (10−9 m3) while adult ticks may ingest ∼1 ml (10−6 m3) (Koch and Sauer, 1984). The luminal surface area of the mosquito midgut is 3 mm2 (3.0 × 10−6 m2) (Houk et al., 1990) and assuming a spherical shape, the volume of the midgut in an engorged female mosquito would be 0.49 mm3 (4.9 × 10−10 m3). Houk et al. (1990) estimated 10,000 epithelial cells in the mosquito midgut, although Smith et al. (2008) estimated that, of these, only 100 midgut cells in Ae. taeniorhyncus mosquitoes were actually susceptible to infection by VEEV even with high-titred blood-meals. For the purpose of the model here, the total number, Ctotal_midgut, of arthropod midgut epithelial cells able to bind virus is assumed to be 103 and the arthropod midgut volume is assumed to be 10−6 dm3 (10−9 m3). The arbovirus challenge dose in the midgut, Vtotal, used here is relatively high and fixed at 105 virions which is representative of vector competence experiments (Liu et al., 2017, Vogels et al., 2016, Zouache et al., 2014) which use a high dose in the blood-meal (typically 107.5 pfu cm−3) to ensure a high chance of infecting a high proportion of the arthropods at optimum temperature. Here it is assumed that all 105 virions have reached the midgut and are contained within the midgut for long enough to reach dynamic equilibrium with epithelial cells according to Equation 12. This Vtotal of 105 virions is realistic based on Kay and Jennings (2002) feeding individual mosquitoes doses of 20,000 TCID50 of Ross River virus and consistent with Smith et al. (2008) feeding ∼10−3 cm3 volume blood-meals of up to 108 VEEV pfu/cm3 when taking into account that the TCID50 and pfu typically comprise many virions in the case of RNA viruses (Rodriguez-Calvo et al., 2011). Lower doses are not considered here. The advantage of using a high Vtotal that greatly exceeds Ctotal_midgut is that as previously demonstrated by Gale (2018) the fraction (Fc) of midgut cells with bound virus at temperature T can be calculated from Ka and the concentration of free virus in the midgut ([Vfree]) using Equation 5. Equation 5 simplifies the calculation and avoids the need for tedious difference equation approaches at low virus challenge doses (Gale, 2018). Furthermore since Vtotal at 105 virions greatly exceeds Ctotal_midgut at 103 cells, Vtotal may be used as an approximation for Vfree in Equation 5. This is analogous to the Michaelis equation in enzyme kinetics where [substrate] >>[enzyme]. Thus even when all the midgut cells have a single bound virus, Vfree = 105 – 103 = 99,000 ∼ Vtotal. Using the Avogadro number, L = 6.02 × 1023 particles/mole, and a midgut volume of 10−6 dm3, [Vfree] is approximately constant at 1.7 × 10−13 M. For each temperature, T, Ka is calculated from Equation 16 and used to calculate Fc from Equation 5. Finally, C.VT is calculated from Equation 6 and used in Equation 7.

Kinetic component of Equation 7: Modelling pcompletet

This is Stage 2 in Fig. 1. The probability pcompleteT depends on rates of biochemical processes including rate of entry of bound virus to the midgut cell, rate of virogenesis and budding and rate of spread through the mosquito to the salivary glands and is assumed here to be a kinetic process related to temperature over the 10 °C to 35 °C range at least. In chemical kinetics, the rate of reaction, k, is modelled by the Arrhenius equation which is based on the activation energy EA. EA is often positive (representing an energy barrier to be overcome) and k increases exponentially with increasing temperature as more molecules acquire the kinetic energy to overcome EA. Rates of viral protein synthesis and nucleic acid replication increase with temperature (below protein denaturation temperatures). For example, Scholtessik and Rott (1969) demonstrated that the rates of synthesis of AIV HA and neuraminidase (NA) increased with temperature from 28 °C to 34 °C and the enzymatic activity of the AIV RNA dependent RNA polymerase (RdRp) increased from 24 °C to 36 °C, peaking at 38 °C and then falling by 40 °C as the enzyme was completely inactivated by heat at 41 °C. Here, the kinetic component of infection of an arthropod by an arbovirus is modelled using the Arrhenius equation where k is replaced by the probability, pcompleteT (Equation 2) and pcomplete283 is the probability that a single midgut cell with bound virus results in successful infection of the salivary glands within the lifetime of the arthropod at a temperature of 283 K i.e. 10 °C. The magnitude of pcomplete283 covers a number of factors in the host pathogen response including the ability of the virus to enter and then replicate in the host cell once attached and then bud (Gale, 2018) and also the immunity of the arthropod to the virus. Thus the more immune the arthropod is to an arbovirus then the lower the magnitude of pcomplete283 in Equation 2 and hence the lower the magnitude of pCompleteT as the probability of spread of the virus through the arthropod decreases. It should be noted that the units of pcompleteT are not those of rate (s−1) because pcomplete283 is a dimensionless probability in Equation 2. It was checked that pcompleteT did not exceed 1 at the temperature of 35 °C as it is a probability.

To identify mechanisms by which the effect of temperature on virus binding may be modulated and hence explain the temperature paradox

This draws heavily on the model results and outputs produced here and also focuses on considering the role of ΔSa_virus on virus/cell binding together with molecular mechanisms which either increase n or change ΔSa_receptor for example through intrinsically disordered regions (IDRs) of the GP/Cr proteins on interaction.

Results

The Kd_receptor for a single SA binding site such as that on an AIV HA monomer is ∼10−3 M (Taylor and Drickamer, 2006, Xiong et al., 2013) and it is assumed that each SA binding pocket on the BTV surface (Zhang et al., 2010), for example, has a Kd_receptor = 10−3 M (37 °C) on binding to a SA on the arthropod midgut epithelial cell surface. For the purposes of modelling, the ΔHa_receptor values required to achieve a Kd_receptor of 10−3 M at 37 °C are calculated in Table 3 using Equation 11 for three scenarios, namely enthalpy-driven, entropy-assisted and entropy-driven as based on their respective ΔSa_receptor assumptions (see Methods).

For the arthropod midgut model developed here the Ka at which half of the 1000 midgut cells have bound virus may be calculated from Equation 5 as 6.0 × 1012 M−1 based on [Vfree] = 1.7 × 10−13 M (i.e. 105 virions per 10−6 dm3 midgut volume). Similarly at the same virus concentration, 999 of the 1000 midgut cells have bound virus with Ka = 6.0 × 1015 M−1. This confirms that for this midgut model Ka ∼1015 M−1 is the upper limit of biological significance for virus binding. If an upper limit of a biologically realistic arbovirus challenge dose is assumed to be 109 virions per arthropod per blood-meal, then [Vfree] = 1.7 × 10−9 M in the model here. For Fc to equal 0.001 representing virus bound to just one of the 1000 midgut cells, Ka must be greater than 6.0 × 105 M−1 with [Vfree] = 1.7 × 10−9 M according to Equation 5. It is concluded that the biologically significant range of Ka is between 106 M−1 and 1015 M−1 for the individual arthropod. This range is marked as horizontal dotted lines on the Ka axis of Figs. 2 and 3 .

Fig. 2

Predicted variation of Ka according to Equation 1 over the biologically relevant temperature range of 10 °C to 40 °C, using ΔHa_virus = −4958 kJ/mol and ΔSa_virus = −16,062 J/mol/K as reported by Carneiro et al. (2002) for binding of vesicular stomatitis virus (VSV) to phospholipid bilayers (solid line). The range of the biologically significant Ka (horizontal dotted lines) relative to the temperature range is controlled by the entropy with ΔSa_virus = −16,562 J/mol/K (dashed line) and −15,562 J/mol/K (dotted line).

Fig. 3

Predicted effect of increasing number of GP/Cr contacts on the binding affinity (Ka) of the virus to its host cell at temperatures of 10 °C (o) and 35 °C (x) according to Equation 16 for enthalpy-driven (a, b, c), entropy-assisted (d) and entropy-driven (e, f) GP/Cr binding scenarios. For graphs a), d) and e) GP/Cr Kd_receptor = 10−3 M at 37 °C with parameters from Table 3. For graphs b) AIV HA/α2–3 sialyllactose contacts and c) AIV HA/α2–6 sialyllactose contacts with parameters from Table 2 (Fei et al., 2015). For graph f) HRV3/ICAM-1 contacts with parameters from Table 2 (Casasnovas and Springer, 1995). Dotted horizontal lines enclose Ka values of biological significance (see text). In each graph, ΔSa_immob = 0 J/mol/K for upper pair of lines (solid) and −750 J/mol/K for lower pairs of lines (dashed).

Predicted effect of temperature on virus binding affinity

The magnitude of Ka calculated with Equation 1 for VSV binding to phospholipid membranes decreases hugely with temperature (Fig. 2) as expected with the large negative ΔHa_virus of −4958 kJ/mol reported by Carneiro et al. (2002). In Table 5, ΔHa_virus values required to achieve a range of Kd_virus values at 37 °C for three different ΔSa_virus scenarios are calculated using Equation 11. For each ΔHa_virus the decrease in Ka on increasing the temperature from 10 °C to 35 °C is calculated using Equation 15. Kd_virus values reported for viruses are in the 10−9 to 10−15 M range (Xiong et al., 2013), for which Ka is decreased by six-fold to 20-fold between 10 °C and 35 °C if ΔSa_virus is negligible. To compensate the huge ΔSa_virus of −16,062 J/mol/K from Carneiro et al. (2002) requires ΔHa_virus in the range of −5000 kJ/mol according to Equation 11 such that the decrease in Ka over the 10 °C to 35 °C temperature range is ∼1076–fold even for the weakest binding i.e. Kd_virus of 10−3 M (at 37 °C). Even using a more realistic value of ΔSa_virus = −750 J/mol/K (derived in Methods) gives decreases in Ka of 19,000-fold to 66,000-fold over the 10 °C to 35 °C temperature range for viruses with Kd_virus of 10−9 to 10−15 M at 37 °C (Table 5).

Table 5

The importance of understanding ΔSa_virus for enthalpy-driven virus binding. Effect of ΔHa_virus on temperature sensitivity of Ka (Equation 15) for a range of Kd_virus values for three ΔSa_virus scenarios.

ΔSa_virus	0 J/mol/K		−750 J/mol/K (See Methods)		−16,062 J/mol/K (from Carneiro et al., 2002)
Kd_virus (M) at 37 °C (310 K)	ΔHa_virus (kJ/mol)*	Decrease in Ka on increasing temperature from 10 °C to 35 °C	ΔHa_virus (kJ/mol)*	Decrease in Ka on increasing temperature from 10 °C to 35 °C	ΔHa_virus (kJ/mol)*	Decrease in Ka on increasing temperature from 10 °C to 35 °C
10−3	−17.8	1.8	−250.3	5646	−4997	0.8 × 1075
10−6	−35.6	3.4	−268.1	10,435	−5015	1.5 × 1075
10−9	−53.4	6.3	−285.9	19,287	−5033	2.7 × 1075
10−12	−71.2	11.7	−303.7	35,645	−5050	5.0 × 1075
10−15	−89.0	21.6	−321.5	65,878	−5068	0.9 × 1076
10−18	−106.8	39.9	−339.3	121,753	−5086	1.7 × 1076

ΔHa_virus required to achieve each Kd_virus at 310 K calculated from Equation 11 given specified ΔSa_virus.

Predicted effect of increasing number of GP/Cr contacts on virus/host cell binding affinity at constant temperature

The logarithm of Ka calculated with Equation 16 increases linearly with the number, n, of GP/Cr contacts made during virus/host cell binding at a given temperature as shown for the enthalpy-driven (Fig. 3a, b, c), entropy-assisted (Fig. 3d) and entropy-driven (Fig. 3e, f) GP/Cr scenarios. It can be seen that the sensitivity of Ka to increasing n depends on the temperature and that this sensitivity is different for the enthalpy-driven (Fig. 3a), entropy-assisted (Fig. 3d) and entropy-driven (Fig. 3e) GP/Cr binding scenarios despite each having the same Kd_receptor of 10−3 M at 37 °C. Thus, the effect of increasing the number of GP/Cr contacts on enhancing virus binding is much more marked at the lower temperature of 10 °C compared to 35 °C for the enthalpy-driven GP/Cr scenarios in Fig. 3a, b, and c with the difference between Ka at 10 °C compared to 35 °C increasing with n according to Equation 15. This is due to the negative ΔHa_receptor in enthalpy-driven binding (Tables 2 and 3) such that binding affinity is increased at lower temperatures as for all exothermic processes. In contrast, the sensitivity of Ka to changing the number of the entropy-assisted GP/Cr contacts in Fig. 3d is less affected by temperature reflecting the smaller magnitude of ΔHa_receptor (Table 3). Indeed, with ΔSa_receptor = +57.4 J/mol/K, a Kd_receptor of 10−3 M (37 °C) is achieved with ΔHa_receptor = 0 kJ/mol (Equation 11) and there is then no temperature effect on Ka (Equation 15) with the lines for all temperatures being overlaid (not shown). While the entropy-assisted binding of GP to Cr minimizes temperature sensitivity of increasing n, the entropy-driven binding of GP to Cr actually reverses the temperature sensitivity such that the effect of increasing the number of GP/Cr contacts on enhancing virus binding is more marked at higher temperatures (Fig. 3e, f). This is due to the positive ΔHa_receptor (Tables 2 and 3) such that binding affinity is increased at higher temperatures as for all endothermic processes. It is assumed that ΔHa_receptor is constant over the 10 °C to 35 °C temperature range which is acceptable.

Non-specific binding of virus (n = 0)

For n = 0, there are no GP/Cr interactions between the virus and the host cell and ln Ka is related to the magnitude of ΔSa_immob according to Equation 18. This represents a virus interacting with an inert surface and with ΔSa_immob set to 0 J/mol/K, Ka = 1 M−1 (from Equation 18) such that there is neither attraction nor repulsion between virus and host cell. However the fraction, Fc, of arthropod midgut cells with bound virus according to Equation 5 is so small at 1.7 × 10−13 that none of the 1000 midgut cells would have bound virus. This simply reflects the very low concentration (1.7 × 10−13 M) of virus at 105 virions per midgut. Indeed mathematically with Ka = 1 M−1 in Equation 5, Fc → [Vfree] which approximates [Vtotal]. In contrast with ΔSa_immob set to −750 J/mol/K, Ka = 10−40 M−1 (Fig. 3a) and there is a very strong repulsion between virus and cell, driven by the unfavourable entropy of immobilization of the virus. This repulsion represents an entropic force as for ΔSmem (Sharma, 2013) and can only be overcome by a large number of GP/Cr contacts (Fig. 3).

The magnitude of ΔSa_immob affects the number of GP/Cr contacts required for Ka to be biologically significant

For ΔSa_immob = 0 J/mol/K, single GP/Cr contacts of Kd_receptor = 10−3 M (37 °C) give low affinity binding of the virus while five GP/Cr contacts give very high affinity binding with Ka in the region of 1015 M−1 depending on the temperature and GP/Cr binding scenario (Fig. 3a, d, e). For ΔSa_immob = −750 J/mol/K, the virus experiences repulsion below n = 10–14 with Ka < 1 M−1 (Fig. 3a, d, e). However, increasing the number of GP/Cr contacts by approximately five to n = 14–20 increases the Ka to 1015 M−1 above which there is no additional increase in infectivity.

Biologically significant binding of AIV HA to soluble SA glycans occurs at the lower temperature of 10 °C with fewer than five GP/Cr contacts if ΔSa_immob is negligible while at 35 °C, the biological significance of binding increases over a wide range of GP/Cr contacts from 10 to 20 for both the α2–3 and α2–6 SA receptors (Fig. 3b, c). However, with the more negative ΔSa_immob of −750 J/mol/K binding is only biologically significant at the lower temperature with more than 15 GP/Cr contacts while at the higher temperature of 35 °C biologically significant binding cannot be achieved with fewer than 70 contacts. This may reflect the unrealistically negative values for ΔSa_receptor (Table 2) addressed above.

For HRV3/ICAM-1 binding, a single GP/Cr contact is sufficient for biologically significant binding assuming the ΔSa_immob is negligible (Fig. 3f), and increasing the number of GP/Cr contacts above three has no significant biological benefit on infectivity. To overcome the larger ΔSa_immob of −750 J/mol/K requires 10 to 11 GP/Cr contacts.

Quantifying arbovirus transmission efficiency

In vector competence experiments, arbovirus transmission efficiency at a given temperature is calculated by dividing the number of female mosquitoes with positive salivary glands (or infected saliva) by the total number of females which were fed an infectious blood-meal 6 to 14 days previously (Liu et al., 2017, Vogels et al., 2016, Zouache et al., 2014). The experimental transmission efficiency may therefore be interpreted as the probability of transmission and modelled as ptransmission_T in Equation 7.

DENV in Aedes albopictus

Liu et al. (2017) fed DENV to Ae. albopictus mosquitoes in blood-meals and showed that the transmission efficiency increased over the temperatures 23 °C, 28 °C and 32 °C (Fig. 4 a). The GP/Cr parameters used in the model here for calculating C.VT are those for the enthalpy-driven scenario in Table 3 with n = 4 GP/Cr contacts of Kd_receptor = 10−3 M (37 °C) and the values of ΔHa_virus and ΔSa_virus are presented in Table 4 . For the kinetic component of Equation 7 (Table 4), a pcomplete283 of 8.0 × 10−7 reflecting the very low probability of virus replication in an arthropod at 10 °C (Wittmann et al., 2002, Carpenter et al., 2011, Mullens et al., 2004) together with an EA of +270 kJ/mol in Equation 2 gave a good approximation of ptransmission_T to the transmission efficiency data of Liu et al. (2017) over the three temperature points as shown in Fig. 4a (solid line).

Fig. 4

Output of thermodynamic/kinetic model in Equation 7 for the effect of temperature on the transmission efficiency of a) DENV by Aedes albopictus mosquitoes (experimental data (■) from Liu et al. (2017)) and c) WNV by Culex pipiens biotype pipiens mosquitoes (experimental data (□) from Vogels et al. (2016)) after oral infection. In a) parameters in Table 4 with number, n, of enthalpy-driven GP/Cr contacts of Kd_receptor 10−3 M (37 °C) equal to 3 (x), 4 (solid line) and 5 (dotted line). In c) the n = 4 enthalpy-driven GP/Cr contacts (solid line) are replaced with n = 4 entropy-assisted (dashed line) or n = 4 entropy-driven (dotted line) GP/Cr contacts of Kd_receptor 10−3 M (37 °C) with parameters in Table 3. Plots b) and d) represent a) and c) respectively with y-axis logarithm-transformed to visualize low transmission probabilities.

Table 4

Parameters to model temperature dependence of transmission efficiency of arboviruses in mosquitoes. Model output Ka values also presented.

Scenario	DENV in Aedes albopictus (Liu et al., 2017)	WNV in Culex pipiens biotype pipiens (Vogels et al., 2016)	CHIKV 2010–1909 strain in Aedes albopictus populations from Brazil (Zouache et al., 2014)
Replication kinetic parameters for Equation 2
pcomplete283 (10 °C)	8.0 × 10−7	6.0 × 10−7	8.2 × 10−4
EA (kJ/mol)	270.0	270.0	105.0
GP/Cr binding parameters	Enthalpy-driven (Table 3)	Enthalpy-driven (Table 3)	Enthalpy-driven (Table 3)
ΔHa_receptor (kJ/mol)	−56.545	−56.545	−56.545
ΔSa_receptor (J/mol/K)	−125	−125	−125
Kd_receptor at 37 °C	10−3 M	10−3 M	10−3 M
Whole virus binding parameters
ΔSa_immob (J/mol/K)	0	0	−32
Number of GP/Cr contacts (n)	4	4	4
ΔHa_virus (kJ/mol)	−226.2	−226.2	−226.2
ΔSa_virus (J/mol/K)	−500	−500	−532
Ka at 10 °C	4.3 × 1015 M−1	4.3 × 1015 M−1	9.2 × 1013 M−1
Ka at 20 °C	1.6 × 1014 M−1	1.6 × 1014 M−1	*3.5 × 1012 M−1
Ka at 37 °C	1.0 × 1012 M−1	1.0 × 1012 M−1	2.1 × 1010 M−1

Reasonable agreement with Ka of 2.2 × 1011 M−1 (20 °C) for WEEV binding to BBFs from susceptible Culex tarsalis mosquitoes (Houk et al., 1990).

The number of GP/Cr contacts up to five has a marked effect on vector competence

The solid line with n = 4 GP/Cr contacts in Fig. 4a gives the best fit to the vector competence data. With n = 3 GP/Cr contacts (crosses), the vector competence is greatly decreased with only 0.15% transmission at 32 °C. That low level transmission does indeed increase with temperature even with just three GP/Cr contacts is more apparent when the transmission probabilities are log-transformed (Fig. 4b). With n = 5 GP/Cr contacts transmission is improved further over n = 4 at higher temperatures with 95% transmission efficiency predicted at 32 °C (dotted line Fig. 4a). Below 20 °C the fifth GP/Cr contact offers no advantage over four (Fig. 4b). Increasing the number of GP/Cr contacts to n > 5 offers no further improvement according to the model (not shown).

WNV in Culex pipiens

The transmission efficiency for WNV by Cx. pipiens biotype pipiens mosquitoes after oral infection increases with temperature between 18 °C and 28 °C (Vogels et al., 2016) and a model differing from that for DENV in Ae. albopictus (Fig. 4a solid line) only in the value of pcomplete283 for the kinetic component (Table 4) fits the three temperature points well (Fig. 4c solid line).

The importance of enthalpy-driven binding for the GP/Cr contacts in promoting vector competence at lower temperatures

Replacing the n = 4 enthalpy-driven GP/Cr contacts in Fig. 4c with n = 4 entropy-assisted or n = 4 entropy-driven GP/Cr binding scenarios (each with Kd_receptor = 10−3 M at 37 °C and parameters in Table 3) not only shifts the predicted transmission efficiency curve to higher temperatures but also steepens the curves such that transmission at lower temperatures is greatly diminished. The low temperature differences are more apparent when the transmission probabilities are log-transformed (Fig. 4d). Thus the lower the temperature, the more the transmission efficiency with entropy-driven GP/Cr binding is diminished relative to that of enthalpy-driven binding although there is no minimum temperature threshold. Thus replacing n = 4 enthalpy-driven GP/Cr contacts with n = 4 entropy-driven contacts decreases transmission by 135-fold at 10 °C but only by 1.4-fold at 35 °C. All three curves in Fig. 4c,d intercept at 37 °C because they each have the same Kd_receptor of 10−3 M at 37 °C. The predicted changes in C.VT with temperature for the three thermodynamic scenarios for WNV binding to Cx. pipiens biotype pipiens midgut epithelial cells are compared in Fig. 5 . Thus on decreasing the temperature, C.VT increases for the enthalpy-driven scenario but markedly decreases for the entropy-driven scenario and is little affected for the entropy-assisted scenario.

Fig. 5

Effect of temperature on the number (C.VT) of midgut cells with bound virus for n = 4 enthalpy-driven (solid line), n = 4 entropy-assisted (dashed line) and n = 4 entropy-driven (dotted line) GP/Cr contacts of Kd_receptor 10−3 M (37 °C) with parameters in Table 3 and ΔSa_immob = 0 J/mol/K. Parameters representative of models for DENV binding in Aedes albopictus mosquitoes and WNV binding in Culex pipiens biotype pipiens mosquitoes (Table 4).

CHIKV in Aedes albopictus: Modelling decreasing vector competence with temperature

Transmission efficiencies for the CHIKV 2010–1909 strain in Ae. albopictus females from a Brazilian mosquito population were higher (statistically significant) at 20 °C than at 28 °C (Zouache et al., 2014) suggesting vector competence decreases with temperature as shown in Fig. 6 a. The kinetic component parameters used for CHIKV were very different from those for DENV and WNV in their respective mosquitoes (Table 4) with CHIKV pcomplete283 some three orders of magnitude higher and CHIKV EA less than half the value of those for DENV and WNV. The GP/Cr parameters used in the CHIKV model here for calculating C.VT are those for the enthalpy-driven scenario in Table 3 with n = 4 GP/Cr contacts of Kd_receptor = 10−3 M (37 °C). The best fit for ptransmission_T to the data for CHIKV peaks at 18 °C and then decreases with temperature up to 35 °C at least and is achieved with ΔSa_immob = −32 J/mol/K (Fig. 6a).

Fig. 6

Output of thermodynamic/kinetic model in Equation 7 for the effect of temperature on the transmission efficiency of CHIKV 2010–1909 strain by Brazilian Aedes albopictus mosquitoes. Experimental data (●) from Zouache et al. (2014). Model parameters in Table 4 with; a) n = 4 enthalpy-driven GP/Cr interactions (Kd_receptor = 10−3 M (37 °C), ΔHa_receptor = −56.5 kJ/mol; ΔSa_receptor = −125 J/mol/K) with ΔSa_immob = −10 J/mol/K (o); −20 J/mol/K (Δ); −32 J/mol/K (x); −40 J/mol/K (□); and −50 J/mol/K (+); b) with number, n, of enthalpy-driven GP/Cr interactions (Kd_receptor = 10−3 M (37 °C), ΔHa_receptor = −56.5 kJ/mol; ΔSa_receptor = −125 J/mol/K) and ΔSa_immob as follows: 4, −32 J/mol/K (x); 5, −97 J/mol/K (□); 6, −162 J/mol/K (Δ); and 7, −227 J/mol/K (o).

The importance of ΔSa_immob on temperature dependence of vector competence

Decreasing ΔSa_immob from −10 J/mol/K to −50 J/mol/K not only decreases the maximum transmission efficiency predicted but also shifts the maximum transmission peak to lower temperatures (Fig. 6a). Thus with ΔSa_immob = −10 J/mol/K, highly efficient transmission (>95%) occurs between 19 °C and 31 °C while with ΔSa_immob = −50 J/mol/K, maximum transmission efficiency is only 40% peaking at 11 °C.

Increasing the number of GP/Cr contacts compensates for more negative values of ΔSa_immob

The baseline model uses n = 4 enthalpy-driven GP/Cr contacts (Table 4) which is sufficient to compensate a relatively small ΔSa_immob of −32 J/mol/K. Each successive decrease in ΔSa_immob by 65 J/mol/K can be compensated by an additional GP/Cr contact up to n = 7 as shown in Fig. 6b, albeit with a steeper effect of temperature on transmission efficiency.

How the opposing effects of temperature on virus binding affinity and virus replication rate may affect vector competence

The opposing effects of temperature on the number (C.VT) of cells with bound virus and the kinetic component (pcompleteT) are compared for the DENV model and the CHIKV model in Fig. 7 . The value of pcompleteT increases more strongly with temperature for DENV reflecting the larger value of EA (Table 4) while pcompleteT for CHIKV is higher for all temperatures reflecting the larger value of pcomplete283. C.VT decreases by three orders of magnitude with increasing temperature for CHIKV but is less affected for DENV.

Fig. 7

Opposing effects of temperature on the kinetic component as represented by pcompleteT (open symbols) and the GP/Cr binding component as represented by the number (C.VT) of midgut cells with bound virus (solid symbols) used in Equation 7 for models for DENV in Aedes albopictus (triangles) and CHIKV 2010–1909 in Brazilian Aedes albopictus (circles). The differences in C.VT are entirely due to the difference of 32 J/mol/K in ΔSa_immob (Table 4).

The magnitude of ΔSa_virus aligns the biologically significant Ka range relative to the biologically relevant temperature range

Taking logarithms of Equation 1, it can be seen that while the magnitude of ΔSa_virus does not affect the sensitivity of Ka to temperature it does control the temperature range over which the change in magnitude of Ka is biologically significant, i.e. between 106 and 1015 M−1. Thus the ΔSa_virus /R term is the y-axis intercept for lnKa when plotted against 1/T. This is apparent in Fig. 2 where with the ΔSa_virus of −16,062 J/mol/K as reported by Carneiro et al. (2002), the change in Ka for binding of VSV to phospholipid bilayers is biologically significant (i.e. 106 M−1 Sa_virus by 500 J/mol/K to −15,562 J/mol/K shifts the range of biologically significant Ka to temperatures above the biologically relevant temperature range, i.e. >40 °C (Fig. 2), such that at all temperatures below 40 °C, the Ka is >1015 M−1. Conversely decreasing ΔSa_virus by 500 J/mol/K to −16,562 J/mol/K hugely decreases the binding affinity (Fig. 2), such that at temperatures above 26 °C the Ka is less than 1 M−1 and there is a “repulsive” force between virus and host cell.

Discussion

The terms enthalpy-driven, entropy-assisted and entropy-driven reflect the relative contributions of the ΔHa_receptor and TΔSa_receptor terms to ΔGa_receptor in Equation 9 as shown in Tables 2 and 3 and are well documented in protein interactions (Bostrom et al., 2011). The enthalpy-driven and entropy-driven scenarios are supported by some experimental data for virus GP/Cr interactions (Table 2), but are modified for modelling in Table 3, firstly so that each GP/Cr interaction scenario has the same Kd_receptor of 10−3 M (37 °C) to represent a single BTV VP2 SA binding pocket interaction with an SA and secondly so the ΔSa_receptor values are not too extreme. X-ray structure data for MERS-CoV and EBOV GPs’ binding to their respective Crs (Lu et al., 2013, Wang et al., 2016, Yuan et al., 2015) suggest many virus GP/Cr interactions are enthalpy-driven although HRV3/ICAM-1 is one example of entropy-driven virus binding as reported by Casasnovas and Springer (1995) (Table 2). Indeed, Carneiro et al. (2002) reported a huge ΔHa_virus of −4958 kJ/mol driving the binding of VSV, and concluded that virus binding was enthalpy-driven. To the author's knowledge entropy-assisted binding has not been documented for virus GP/Cr systems, although it has been demonstrated for protein/protein interactions (Table 2) and for about a quarter of the 28 protein/ligand complexes studied by Caro et al. (2017) the conformational entropy change assisted the driving of specific binding. There is much indirect evidence to support the important contribution of entropy in GP/Cr binding. Thus according to De Graaf and Fouchier (2014) among the major difficulties of X-ray crystallography structural studies of AIV HA interacting with SA receptors are the flexible nature of both the HA protein and the glycan structure, which is largely undetectable in crystallography. This motional freedom is consistent with the importance of and in the GP/Cr interaction, as is the presence of IDRs in many viral proteins (Dolan et al., 2015) and their receptors (Fong et al., 2009).

The Ka values vary by orders of magnitude depending on temperature and n (Figs. 2 and 3). However, the biologically significant range of Ka is from 106 to 1015 M−1 such that for Ka >1015 M−1, temperature induced changes in Ka have no biological effect on infectivity either because at high virus challenge doses all the 1000 arthropod midgut cells have bound virus, or because at low virus challenge doses all the virus is bound to midgut cells (Gale, 2018). Similarly, Kas below 106 M−1 require unachievably high virus challenge doses for a realistic chance of infection at the individual arthropod level. It should be noted that a Ka of <106 M−1 does not represent a threshold just because Fc is less than 0.001 (albeit with an unrealistically high challenge dose of 109 virions), and it should be stressed that Fc decreases linearly with decreasing Ka (Gale, 2018). Indeed at the arthropod population level where infection of 1 in 1000 arthropods may be epidemiologically significant in disease transmission, the low Ka values may be biologically important. Thus a Ka of 6.0 × 106 M−1 would give virus binding in 1 in 106 midgut cells, i.e. 1 in 1000 midgut cells in just one of a 1000 individual arthropods (each receiving a blood-meal of 105 virions which is realistic). It is suggested that Ka ∼106 M−1 may be epidemiologically significant in arbovirus transmission depending on the number of arthropods and their biting rates. The biologically significant range of Ka (106 M−1 to 1015 M−1) for the individual arthropod midgut is very similar to that for the human intestine model (Gale, 2018), suggesting the effective range of Ka is not affected by the volume of the midgut (10−6 dm3) or the intestine (0.314 dm3; Gale (2018)), and would also apply to other virus routes such as AIV in the fluid in a lung for example. This is not surprising given that Ka is dimensionless.

The results in Table 5 and Fig. 2 highlight the importance of understanding the magnitude of ΔSa_virus and hence ΔSa_immob in terms of effect of temperature on Ka. Thus the more negative ΔSa_virus is, then the more negative ΔHa_virus has to be to overcome it and the greater the sensitivity of Ka to temperature according to Equation 15. The results in Table 5 and Fig. 2 demonstrate that Ka may be decreased by orders of magnitude for enthalpy-driven GP/Cr binding over relatively small temperature increases depending on the magnitudes of n and ΔHa_receptor. Since the “n x ΔHa_receptor” term (and hence ΔHa_virus) is constant, the change in Ka with temperature according to Equation 15 does not depend on whether there is one single interaction of enthalpy ΔHa_virus kJ/mol or multiple (n) interactions each of ΔHa_receptor kJ/mol which when summed give ΔHa_virus according to Equation 3. This is important because it might be thought that having multiple interactions each of very small ΔHa_receptor could alleviate the effect of temperature on Ka compared to a single large ΔHa_virus for example from a very high affinity single GP/Cr interaction. Therefore, information on the number, n, of GP/Cr contacts per virus/host cell interaction is as important as the magnitude of ΔHa_receptor itself.

That Ka according to Equation 16 increases logarithmically with n at constant temperature in Fig. 3 is apparent from Equation 18 which separates the GP/Cr receptor binding component (ΔGa_receptor) from the whole virus particle entropy component (ΔSa_immob) and gives a straight line when ln Ka is plotted against n at constant temperature with the y-axis intercept being related to ΔSa_immob. For strong binding to occur between GP and Cr, Kd_receptor must be much less than 1 M, and therefore ΔGa_receptor must be negative (<0) according to Equation 8. For example, for Kd_receptor = 10−3 M, ΔGa_receptor = −17.8 kJ/mol at 310 K (37 °C). Thus despite the negative n term in Equation 18, the logarithm of Ka increases linearly with n at constant temperature for all three thermodynamic GP/Cr binding scenarios.

Using the three “realistic” thermodynamic GP/Cr binding scenarios in Table 3 the results in Fig. 3a,d,e show that it is the number of GP/Cr contacts and the magnitude of ΔSa_immob which are important in affecting Ka and in particular in achieving a biologically significant Ka as indicated by the horizontal dashed lines in Fig. 3. The relative lack of importance of temperature between 10 °C and 35 °C reflects the relatively small values of ΔHa_receptor (Table 3) according to Equation 15. However, for AIV HA/siallylactose binding (Fig. 3b,c) the much more negative values of ΔHa_receptor (Table 2) result in much lower Ka values at the higher temperature (35 °C) so that many more GP/Cr contacts are needed to achieve a biologically significant Ka. Thus to achieve biological significance at higher temperatures with only a few GP/Cr contacts, it is beneficial for ΔHa_receptor not to be too negative, hence the potential benefit of entropy-assisted or entropy-driven GP/Cr binding.

Here a two stage thermodynamic/kinetic approach is developed to model the effect of temperature on the competence of arthropod vectors to transmit arboviruses. Through Equation 7, the model separates the thermodynamic equilibrium component as represented by C.VT from the kinetic component as represented by pcompleteT. The advantage of this is that it allows the opposing effects of temperature on each component (Fig. 7) to be accommodated and treated separately in calculating the transmission efficiency, ptransmissionT. A notable success is its ability to accommodate both increasing (Fig. 4) and decreasing (Fig. 6) vector competence with increasing temperature as observed in laboratory mosquito infection studies. The model also allows changes in biophysical parameters that can be measured by experiments to be tested. For example, Fig. 4a shows there is no gain in overall vector competence on increasing n above 5 according to the model here. Interestingly the upper limit of biological significance for Ka of 1015 M−1 is equivalent to the maximum of 10−15 M reported for Kd_virus by Xiong et al. (2013) for AIV binding involving 5 GP/Cr contacts each of Kd_receptor = 10−3 M. In evolutionary terms there may be no point increasing n above 5, although this depends on the magnitude of ΔSa_immob which is currently unknown.

Justification of the approach used

Although the second (kinetic) stage also involves cell binding by the virus (e.g. to the salivary glands) this is ignored here in effect assuming that midgut cell brush border binding in the first stage is the most important binding process. This is justified because differences in mosquito vector competence are in part due to arbovirus-specific midgut epithelial Cr receptors on the brush borders (Franz et al., 2015). Indeed Western equine encephalitis virus (WEEV) bound with high affinity to brush border fragments (BBFs) from the midguts of susceptible Cx. tarsalis mosquitoes but not to those from refractory mosquitoes, confirming the importance of the appropriate Cr molecules (Houk et al., 1990). The second stage lumps together all the infection stages of the arthropod following midgut cell surface binding and is entirely kinetic being based on the Arrhenius equation (Equation 2). Justification of this is as follows. First, infection of just one midgut cell is sufficient for a 50% probability of subsequent infection of the salivary glands in the case of VEEV in Ae. taeniorhynchus (Smith et al., 2008). This suggests the midgut escape, dissemination and salivary gland stages are more efficient than midgut infection and are not major barriers and can therefore be lumped together into a single kinetic component reflecting the time required for viral replication and spread within the arthropod. This is supported by the results of Vogels et al. (2016) showing that while oral challenge of Cx. pipiens mosquitoes with WNV resulted in salivary gland infection in only 2% – 35% of mosquitoes depending on temperature (Fig. 4c), intrathoracic injection of WNV achieved salivary gland infection in 95% – 100% mosquitoes at all temperatures (not shown). Intrathoracic injection by-passes the midgut entry, infection and escape stages. Second, Scholtessik and Rott (1969) concluded that the activity of RdRp is rate-limiting in AIV infection of cell cultures. Thus by analogy all that is needed for the kinetic component of the vector competence model here are the EA and pcomplete283 values for the rate-limiting step of the midgut cell infection, dissemination and salivary glands stages. This rate-limiting step being the slowest of these will limit pcompleteT and hence overall vector competence in Equation 7.

Validation and appropriateness of input parameters

A formal validation cannot be performed here because of the lack of the required thermodynamic/kinetic parameters for arboviruses. For the next step in the development of this mechanistic model, data for the kinetic parameters, pcomplete283 and EA, and the thermodynamic parameters, n, ΔHa_receptor, ΔSa_receptor and ΔSa_immob for arbovirus molecular components set out in Table 4 should be obtained experimentally and the model outputs validated against competence data from insectary experiments as shown for Figs. 4 and 6. Parameterization of ΔHa_receptor and ΔSa_receptor for arboviruses as Fei et al. (2015) presented for AIV HA may be achieved in the near future and furthermore changes in ΔHa_receptor may be predictable from molecular data as described previously (Gale, 2018) particularly as more detailed information becomes available on arbovirus structure and Cr binding (Zhang et al., 2010). Due to the current paucity of Kd_receptor data for arboviruses, the models here are based on multiple integer GP/Cr contacts each of Kd_receptor = 10−3 M (37 °C) to represent a typical SA binding site on an arbovirus for the purpose of demonstration. That cell surface SA is important for BTV binding was demonstrated by the finding that blocking the SA on the cell surface reduced production of BTV in cell culture by ∼20-fold (Zhang et al., 2010). The thermodynamic parameters used here in Table 4 are believed to be realistic and reassuringly the model outputs give good approximations to the experimental vector competence data (Figs. 4 and 6). It is also reassuring to note that the Ka of 3.5 × 1012 M−1 predicted for CHIKV binding to Ae. albopictus at 20 °C in the model (Table 4) is within two orders of magnitude of that of 2.2 × 1011 M−1 measured experimentally for the arbovirus WEEV binding to Cx. tarsalis BBFs at 20 °C (Houk et al., 1990). This offers some validation of the binding component of the model at least. For the kinetic component, while no data are available for pcomplete283, Scholtessik and Rott (1969) presented EA values for AIV RdRp of +83.7 kJ/mol and +150.6 kJ/mol, which are remarkably similar to the values of +105 kJ/mol and +270 kJ/mol used in the models (Table 4). In the absence of data for pcomplete283 the values used here in the kinetic component were chosen to give good approximations to the data in Figs. 4 and 6 to demonstrate proof of principle. The parameter pcompleteT is related to ν in Turner et al. (2013) which is the rate at which latently-infected Culicoides midge vectors (i.e. those exposed to an effective BTV dose through feeding) become infectious and is the reciprocal of the EID in units of per day. Turner et al. (2013) calculate ν for a given temperature T (in °C) using ν = 0.0003T(T – 10.4) from Mullens et al. (2004). In effect, ν is the probability that the EIP is 1 day and at 11 °C for example has a value of 0.002 according to the equation of Mullens et al. (2004). The value of pcompleteT at 11 °C for CHIKV in Fig. 7 is very similar at 0.001, although pcompleteT is not in units of per day but is by definition the probability over the lifetime of the vector in the model here. Differences in pcomplete283 and EA at the biological level could represent the ability of virus to replicate in the midgut or salivary gland cells, for example, with the magnitude of EA representing how the efficiencies in these processes increase with temperature.

The importance of getting data on ΔSa_immob

The importance of ΔSa_immob is apparent from the results that the differences in the temperature sensitivities of C.VT in DENV and CHIKV in Fig. 7 are solely due to the difference of 32 J/mol/K in the magnitudes of ΔSa_immob (Table 4). The potential importance of ΔSa_immob in modelling vector competence is also shown in Fig. 6a through its affecting both the peak temperature and the maximum transmission efficiency for CHIVK strain 2010–1909. Fig. 6b shows how increasing the number of GP/Cr contacts counteracts the effect of decreasing ΔSa_immob in terms of vector competence as expected from Fig. 3. However information on the magnitude of ΔSa_immob at present is limited and contradictory. Thus while Carnerio et al. (2002) report a huge ΔSa_virus of −16,062 J/mol/K which is consistent with a large ΔSa_immob in agreement with the theory for immobilization of the high molecular weight virion with its large moments of inertia, the models in Fig. 3a, d and e show that large numbers (n = 12 to 18) of GP/Cr contacts are required to overcome a ΔSa_immob of −750 J/mol/K. This is irrespective of whether the GP/Cr binding is enthalpy driven, entropy-assisted or entropy-driven. There is no evidence that such large numbers of contacts occur during virus binding, most being estimated to be between 5 (Xiong et al., 2013) and 7 (Carneiro et al., 2002) contacts. Since using a ΔSa_immob of 0 J/mol/K gives more realistic results in terms of the number of GP/Cr contacts required, it is tempting to conclude that the contribution from ΔSa_immob to ΔSa_virus is small even for the results of Carneiro et al. (2002). For this reason, the arbovirus vector competence models in Figs. 4 and 6 are based on low values for ΔSa_immob (Table 4). Xiong et al. (2013) appear to assume ΔSa_immob is negligible by not considering it. Thus setting ΔSa_immob = 0 J/mol/K in Equation 17 expresses Ka as the reciprocal of Kd_receptor raised to the power of n as assumed by Xiong et al. (2013).

ΔSa_immob is a feature of the whole virus and is not related to either viral GP or host Cr and in this respect could reflect differences in the morphology and aggregation state of the virus particle (Gale, 2017, Gale, 2018). It is interesting to note that only the CHIKV 2010–1909 strain showed decreased transmission efficiency at 28 °C compared to 20 °C for five of the six Ae. albopictus mosquito populations (Zouache et al., 2014). In contrast, CHIKV 06–021 strain showed higher transmission efficiency at 28 °C compared to 20 °C for five of the six mosquito populations. This suggests the decreased transmission at higher temperature is a feature of the 2010–1909 strain and consistent which a property such as ΔSa_immob which is virus specific and independent of the host cell. However, ΔSa_immob alone was not sufficient to change vector competence from decreasing with temperature to increasing with temperature (Fig. 6), suggesting that the kinetic component is also important as is now discussed.

From Fig. 7 it can be seen that the temperature dependence of DENV transmission in Ae. albopictus is more influenced by the effect of temperature on the kinetic component than by the effect of temperature on C.VT. Indeed C.VT for DENV is relatively unaffected particularly at the lower temperatures because the 1000 midgut cells are almost fully saturated with virus with Ka > 1015 M−1 (Table 4). Thus DENV C.VT only starts to fall significantly at temperatures above 25 °C as Ka decreases well below 1014 M−1 (Fig. 7). In contrast for the CHIKV model, Ka is less than 1014 M−1 even at 10 °C (Table 4) and the midgut cells are not saturated. Thus CHIKV C.VT experiences the full effect of increasing temperature, but DENV does not in this model. This is solely due to the difference of 32 J/mol/K for ΔSa_immob which decreases the Ka for CHIKV to below 1014 M−1 such that Ka is 46-fold lower for CHIKV than for DENV at all temperatures (Table 4). In terms of overall impact of temperature on vector competence, the strongly-affected kinetic component due to the large value of EA for DENV dominates over the less affected C.VT and vector competence increases overall with temperature for DENV and for WNV (Fig. 4) which also has a large EA (Table 4). In contrast for CHIKV, the effect of temperature on the kinetic component is much weaker (Fig. 7) due to the smaller EA (Table 4) and there is a peak in the vector competence (Fig. 6) due to the opposing effects of temperature on C.VT and pcompleteT. Thus for CHIKV 2010–1909 in Brazilian Ae. albopictus, the strong effect of temperature on reducing C.VT wins at higher temperatures over the slowly increasing pcompleteT giving a peak followed by a decrease in transmission efficiency at temperatures above 18 °C (Fig. 6). The transmission probabilities for increasing competence for DENV and WNV in Fig. 4 would also peak, but at higher temperatures around 36 °C-38 °C as the vector/virus enzymes are denatured at 41 °C (Scholtessik and Rott, 1969) although some WNV strains vary in their temperature sensitivity (Kinney et al., 2006). This peak is not addressed here and arthropod temperatures are unlikely to exceed 35 °C particularly in northern Europe.

The thermodynamic scenario for GP/CR binding markedly affects vector competence at low temperature

The results here show that information on whether GP/Cr binding is enthalpy-driven, entropy-assisted, or entropy-driven is important because it affects the temperature sensitivity of Ka (Fig. 3) and therefore C.VT (Fig. 5) and markedly influences the effect of temperature on arbovirus transmission efficiency particularly at lower temperature (Fig. 4c, d). Thus the increased virus/host cell binding affinity at lower temperature of the enthalpy-driven GP/Cr contacts (Fig. 3a) greatly promotes transmission efficiency at lower temperatures giving a good fit to the data of Vogels et al. (2016) for WNV transmission by Cx. pipiens in Fig. 4c, d. In effect higher C.VT at lower temperatures with the enthalpy-driven GP/Cr binding (Fig. 5) compensates for the lower replication kinetic rates (Fig. 7) to some degree such that significant transmission can still occur at lower temperatures. For the entropy-assisted scenario, virus/host cell binding affinity (Fig. 3d) and C.VT (Fig. 5) are little affected by temperature such that the effect of temperature on vector transmission efficiency (Fig. 4c) is predominantly determined by the kinetic component (Equation 2). For entropy-driven GP/Cr binding, the decrease in C.VT with decreasing temperature (Fig. 5) together with the decreasing kinetic rate (Fig. 7), act synergistically to greatly diminish transmission efficiency and hence vector competence particularly at lower temperature compared to enthalpy-driven GP/Cr binding (Fig. 4d). Although not shown in the results here, decreasing the value of EA and increasing pcomplete283 in the kinetic component (Table 4) shifts the curve (dotted line in Fig. 4c) for entropy-driven GP/Cr binding back to the left giving a reasonable fit to the WNV data points of Vogels et al. (2016). However at 10 °C, the transmission efficiency predicted is only a third of that of the enthalpy-driven model (solid line in Fig. 4d). Thus there is a complex interplay between the thermodynamic scenario used for GP/Cr binding and the kinetic component that influences transmission efficiency at low temperatures.

Minimum temperature thresholds

The minimum temperature threshold is an important consideration for assessing the risks of BTV transmission by Culicoides midges in northern Europe. Carpenter et al. (2011) and Wittmann et al. (2002) specify minimum temperature thresholds of 11 to 13 °C for different orbiviruses in Culicoides midge species. Indeed the predicted EIP is infinitely long at a temperature of 10.4 °C according to the model of Mullens et al. (2004) such that transmission is zero at this temperature. The model here does not have a low temperature threshold because pcomplete283 is not zero (Table 4) and even with entropy-driven GP/Cr binding, the predicted transmission efficiency is not zero at 10 °C (dotted line Fig. 4d). At 10 °C the model for DENV in Ae. albopictus with n = 4 (or 5) enthalpy-driven GP/Cr contacts predicts a ptransmissionT of 0.0008 (solid line Fig. 4b), such that 1 in 1250 arthropods would transmit at this low temperature. Even with n = 3 GP/Cr contacts, the model predicts 1 in 15,300 arthropods would transmit at 10 °C (crosses, Fig. 4b). This absence of an absolute threshold predicted by the thermodynamic model here is entirely consistent with the Maxwell-Boltzmann distribution for velocities of gas molecules which shows that a proportion of molecules have high kinetic energies even at lower temperatures, and also with the fact that at low temperatures adult arthropods live for longer, giving more time for virus replication to complete. This however does not take into account biomolecular mechanisms for a threshold which could include cold lability of proteins (i.e. some proteins denature at low temperatures) and also molecular phase transitions e.g. in phospholipid bilayers at lower temperatures (Gale and Watts, 1991). Thus lower lateral diffusion rates of Cr molecules within the plane of the membrane due to reduced fluidity of the phospholipid bilayer may restrict recruitment of a sufficient number of Cr molecules for high affinity virus binding. For this reason values of ptransmissionT may be much lower than suggested by the model here at lower temperatures. If there is molecular evidence that replication cannot occur below (or above) a certain temperature, then pcompleteT should be set to zero over those temperatures, and ptransmissionT will also equal zero according to Equation 7. In addition there are temperature-specific biological factors to consider including for example the effect of temperature on core inhibition of RdRp in the case of BTV (van Dijk and Huismans, 1982) and host immune responses (Foxman et al., 2015). These biological factors could mean that the change in magnitude of pcompleteT with temperature is not as straight forward as suggested by Equation 2.

The temperature paradox presented here is that according to the Van't Hoff Isochore, the stronger the “molecular fit” between GP and Cr (i.e. the more negative in magnitude ΔHa_receptor is in Equation 15), then the greater the decrease in Ka with temperature and the weaker the binding at higher temperatures (Table 5 and Fig. 2). This is consistent with the data of Houk et al. (1990) for binding of 35S-labelled WEEV to mosquito BBFs. Thus the specific binding of WEEV decreased from 14% at 20 °C (at which binding was optimal) to just 3% at 40 °C for BBFs from susceptible Cx. tarsalis mosquitos. Fei et al. (2015) showed experimentally that the GP/Cr binding affinity for enthalpy-driven binding of AIV HA to SA glycans decreased by 10-fold over a 10 °C increase in temperature. Indeed results in Fig. 3b,c for AIV HA/sialyllactose binding suggest that AIV would not be very infectious at mammalian or avian body temperatures which is clearly not correct. The effect of temperature on Ka is potentially huge as demonstrated for VSV in Fig. 2. Paradoxically, this does not always relate to the observed effect of temperature on published virus infectivity data or even to all published virus binding data. Indeed the opposite is observed in some cases. Thus, Nunes-Correia et al. (1999) reported that the binding affinity of AIV H1N1 to canine kidney cells increased by over 10-fold from 4 °C to 37 °C. In the case of filovirus infectivity, Miller et al. (2016) demonstrated that at 72 h, titres of Ebola virus Makona and Marburg virus Angola were 1000-fold and 100-fold higher, respectively, at 41 °C than at 37 °C in Rousettus aegyptiacus cell lines (although they were both similar at the 96 h end point). This is not consistent with a negative ΔHa_receptor for GP binding to Cr and indeed suggests a positive ΔHa_receptor as in entropy-driven binding (Fig. 3e, f and Fig. 5). Of course temperature kinetic effects on rates of viral macromolecule replication (Scholtessik and Rott, 1969) would also be important here (Fig. 7), but without efficient cell binding, infection would not initiate in the first place.

The work here identifies three thermodynamic/mechanistic reasons why strong virus binding and/or infectivity could still occur at high temperatures. First, increasing the number of GP/Cr contacts made on virus/cell binding by just one or two at higher temperature compensates the decrease in Ka due to temperature for enthalpy-driven GP/Cr binding (Fig. 3a). The impact on vector competence is shown in Fig. 4a with transmission efficiency increasing greatly at the higher temperatures as n is increased from 3 to 5. This mechanism requires n to increase with temperature as has been reported for HIV binding through increasing host membrane fluidity allowing effective recruitment of more Cr molecules to bind the virus (Harada et al., 2004) and is consistent with the observed increase in specific binding of WEEV to mosquito BBFs from 3% to 14% on increasing the temperature from 4 °C to 20 °C (Houk et al., 1990). The second is that GP/Cr binding is entropy-driven for some viruses as shown for HRV3 such that C.VT increases with temperature (Fig. 5). Alternatively the temperature sensitivity of binding can be eliminated through entropy-assisted binding as ΔHa_virus →0 kJ/mol (Fig. 3d and Fig. 5). The third (which applies to infectivity and not to binding affinity per se) is that the Ka values are >1015 M−1 over the biologically relevant temperature range (Fig. 2) such that changes in Ka, although large in magnitude, have no observed effect on C.VT and hence on infectivity. An example of how temperature-induced changes in Ka outside the biologically significant range of Ka have relatively little effect (over the low temperature range at least) on infectivity is apparent in the plot of C.VT as a function of temperature in Fig. 7 for DENV in Ae. albopictus.

In Fig. 2 it is shown how changes in ΔSa_virus of ±500 J/mol/K are sufficient to select the temperature window over which Ka for VSV binding is either within or outside the biologically significant range of 106 to 1015 M−1. It is suggested here that through the and terms in Equation 10, the IDRs in the GP and Cr molecules (Fong et al., 2009, Dolan et al., 2015) could control ΔSa_receptor and ΔSa_virus and hence dictate the temperature range over which Ka is biologically significant. Thus, based on the estimated change in conformational entropy (ΔΔSconf) of ∼ 6.1 J/mol/K per amino acid residue in an IDR of a protein (Rajasekaran et al., 2016), “adjusting” ΔSa_virus by ±500 J/mol/K would require a change in the order/disorder status of 82 amino acids in the VSV GP/Cr system. This is across seven VSV G proteins (Carneiro et al., 2002) and hence seven Cr molecules and would equate to ∼12 amino acid residues per GP/Cr pair. IDRs in proteins may be up to 23 to 30 amino acid residues in length (Fong et al., 2009, Liu and Huang, 2014). Thus a change in the order/disorder status of 12 amino acid residues per GP/Cr pair would be quite feasible and across seven GP/Cr pairs involved in binding would be sufficient to change ΔSa_virus by ±500 J/mol/K. VSV is an arbovirus and its life cycle involves horizontal transmission from insects to mammals (Za'arate and Novella, 2004). Using the ΔSa_virus value of −16,062 J/mol/K from Carneiro et al. (2002) predicts that VSV binds to phospholipid bilayers very strongly at arthropod temperatures (25 °C) but not at mammalian temperatures (37 °C) (Fig. 2). Strong binding at 37 °C is achieved with ΔSa_virus = −15,562 J/mol/K. Although speculative, it is suggested here that the order/disorder status of the IDRs of the VSV GP/mammalian cell Cr pair could shift ΔSa_virus by +500 J/mol/K such that VSV binds strongly at mammalian body temperature. Although outside the scope of the work here, some mechanism is needed to facilitate arbovirus binding to mammalian host cells which in the case of avian species for WNV could be at temperatures as high as 44 °C (Kinney et al., 2006) at which Ka would be greatly diminished compared to 25 °C, depending on the thermodynamics of GP/Cr binding.

Conclusion

A model for the effect of temperature on the transmission efficiency of an arbovirus by an arthropod vector is developed based on a simulated arthropod midgut and a thermodynamic/kinetic framework that naturally accommodates temperature. The model can explain both increasing and decreasing arbovirus transmission efficiency with temperature by allowing for the opposing effects of temperature on the affinity of virus binding (as defined by Ka) and on the rate of viral replication and spread to the salivary glands. The model is based on four thermodynamic/kinetic parameters, ΔHa_virus, ΔSa_virus, EA and pcomplete283, some of which may be measured experimentally in biochemical “test tube” experiments, and it is anticipated, may ultimately be related to gene sequence data for the arbovirus and its arthropod vector. The thermodynamic nature of GP/Cr binding (i.e. enthalpy-driven, entropy-assisted or entropy-driven) influences the sensitivity of Ka, and hence vector competence, to temperature. GP/Cr binding for many viruses is enthalpy-driven and there is evidence for entropy-driven binding for at least one virus. Specifically enthalpy-driven GP/Cr binding promotes transmission efficiency at low temperature while entropy-driven GP/Cr binding diminishes transmission efficiency at low temperatures although there is no minimum temperature threshold according to the model. Enthalpy-driven binding, however, presents a paradox in that the stronger the binding at lower temperatures, then the weaker the binding at the higher temperatures at which virus replication rates are optimal. This may not necessarily result in lower infectivity at higher temperatures. Thus, according to the arthropod midgut model developed here, changes in Ka above 1015 M−1 have little effect on the number of midgut cells with bound virus hence masking the underlying effect of temperature on Ka. Furthermore higher temperatures may allow formation of more GP/Cr contacts which have a dominant effect (over the decrease in Ka with temperature for enthalpy-driven binding) on increasing Ka. According to the model increasing the number of enthalpy-driven GP/Cr contacts from three to five greatly increases transmission efficiency at the higher temperatures. Although entropy-assisted binding has not yet been documented for virus binding, it would minimize the effect of temperature on Ka. The potential role of the entropy of intrinsically disordered regions of GP and Cr in controlling the temperature sensitivity window of virus binding should be investigated in relation to binding at arthropod temperatures and the higher temperatures of the avian or mammalian host reservoirs. The thermodynamic approach developed here identifies the entropy of immobilization of the whole virus particle on binding to the host cell (ΔSa_immob) as a key parameter, the “repulsive” effects of which can be offset by increasing the number of GP/Cr contacts made. According to the model, even small differences in ΔSa_immob strongly influence the effect of temperature on the number of midgut cells with a bound virus and hence infectivity. It is anticipated that thermodynamic and kinetic data obtained from biochemical experiment and the interpretation of gene sequences could be applied to modelling the effect of temperature on arthropod vector competence through this model.