Susceptibility to organophosphates pesticides and the development of infectious-contagious respiratory diseases

Susceptibility to organophosphates pesticides and the development of infectious-contagious respiratory diseases

Journal Pre-proof

Susceptibility to organophosphates pesticides and the development of infectious-contagious respiratory diseases ´ ´ ˜ J.P. Gutierrez-Jara, F.D. Cordova-Lepe, M.T. Munoz-Quezada, G. Chowell PII: DOI: Reference:

S0022-5193(19)30502-8 https://doi.org/10.1016/j.jtbi.2019.110133 YJTBI 110133

To appear in:

Journal of Theoretical Biology

Received date: Revised date: Accepted date:

11 April 2019 29 November 2019 19 December 2019

´ ´ ˜ Please cite this article as: J.P. Gutierrez-Jara, F.D. Cordova-Lepe, M.T. Munoz-Quezada, G. Chowell, Susceptibility to organophosphates pesticides and the development of infectious-contagious respiratory diseases, Journal of Theoretical Biology (2019), doi: https://doi.org/10.1016/j.jtbi.2019.110133

This is a PDF file of an article that has undergone enhancements after acceptance, such as the addition of a cover page and metadata, and formatting for readability, but it is not yet the definitive version of record. This version will undergo additional copyediting, typesetting and review before it is published in its final form, but we are providing this version to give early visibility of the article. Please note that, during the production process, errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain. © 2019 Published by Elsevier Ltd.

Highlights • An intoxication-diseases dynamic model with genetic susceptibility levels is studied • R0 value is expressed as a weighted average by levels susceptibility to the toxic • The genotypic distribution to the toxic impacts the level epidemic of the disease

1

Susceptibility to organophosphates pesticides and the development of infectious-contagious respiratory diseases J.P. Guti´errez-Jaraa,∗, F.D. C´ ordova-Lepeb , M.T. Mu˜ noz-Quezadac , d G.Chowell a Doctorado

en Modelamiento Matem´ atico Aplicado, Universidad Cat´ olica del Maule, Avenida San Miguel 3605, Talca, Chile. b Facultad de Ciencias B´ asicas, Universidad Cat´ olica del Maule, Avenida San Miguel 3605, Talca, Chile. c Facultad de Ciencias de la Salud, Universidad Cat´ olica del Maule, Avenida San Miguel 3605, Talca, Chile. d School of Public Health, Georgia State University, Decatur 140 Georgia, Atlanta, USA.

Abstract In this paper we develop an SIRS compartmental model to investigate the dynamic interplay between pesticide intoxication and the spread of infectiouscontagious respiratory diseases. We are particularly interested in investigating three levels of genetic susceptibility to the toxic. The genotypic distribution, corresponding to the Q and R alleles, is proposed and parameterized according to ethnic variation using real population data from published studies. We then use the mathematical models to illustrate the impact of this distribution on the spread of hypothetical influenza-like respiratory disease in a population exposed to the organophosphate pesticide. In this context, we show how an initial basic reproductive number below the epidemic threshold of 1.0 could be enhanced to support epidemic outbreaks in agricultural populations that employ chlorpyrifos pesticides. We further illustrate our modeling framework to study the effect of ethnic group variation in Singapore (Malay, Indian and Chinese) using genetic distribution data from published studies. ∗ Corresponding

author Email addresses: [email protected] (J.P. Guti´ errez-Jara), [email protected] (F.D. C´ ordova-Lepe), [email protected] (M.T. Mu˜ noz-Quezada), [email protected] (G.Chowell)

Preprint submitted to Journal of Theoretical Biology

December 20, 2019

Keywords: Pesticides, Intoxications, Mathematical model SIRS, Genetic susceptibility 2010 MSC: 92B05 - 92D30.

1. Introduction Pesticides are a broad class of chemical substances that are widely employed in agriculture for controlling, preventing or eradicating pests (FAO, 2003). However, their inappropriate can exert adverse health effects (Bakhsh et al., 2016; Damalas and Koutroubas, 2016; Mu˜ noz-Quezada et al., 2013, 2016). Because of their widely use in agricultural activities (Marcelino et al., 2019), agricultural workers and their families are at an increased risk of pesticide exposure. In particular, organophosphates (OP) exposure in humans inhibits the acetylcholinesterase, leading to accumulation of acetylcholine, which, in turn, affects the functioning of the nervous impulse (Dutta and Bahadur, 2018). The paraoxonases are a family of enzymes in mammals composed of three elements: PON1, PON2 and PON3. The most studied and applied is PON1 (M´ arquez-Llanos, 2015), which has the ability to hydrolyze OP insecticides and their corresponding metabolites (Costa et al., 2004; Herrera and Ticona, 2011). Consequently, a low hydrolysis rate of PON1 is linked to greater susceptibility to OP (Costa et al., 2013; Mackness et al., 2003). Importantly, prior studies support genetic susceptibility to organophosphate pesticide (Costa et al., 2000). Specifically, three levels (high, medium and low) of genetic susceptibility to this toxic have been previously identified (Costa et al., 2000). Moreover, these levels of genetic predisposition are linked to the genetic polymorphism at codon 192 with its alternative forms Q (Glutamine) and R (Arginine), commonly denoted by PON1 192 Q/R, depending on the OP under consideration (Hofmann et al., 2010; Singh et al., 2011). For example, the Q allele is less efficient than the R allele in hydrolyzing Paraoxon or Chlorpyrifos. Further, the Q allele is more efficient at hydrolyzing, Diazoxon, Soman and Sarin. Also, both alleles are similarly efficient at hydrolyzing Phenyl acetate

3

(Z´ un ˜iga-Venegas et al., 2015). In summary, heterogeneous genetic susceptibility in humans highlights the need to understand human variation to pesticide exposure and the underlying genetic variation across different geographic or ethnic groups (Povey, 2010). Meanwhile, the dynamics of infectious disease outbreaks (e.g., the probability of invasion) are heavily influenced by susceptibility variation within a population (Guti´errez-Jara et al., 2019; Hickson and Roberts, 2014). A case in point in this paper is the extent of pesticide intoxication on the dynamics of respiratory diseases (Hoppin et al., 2007; Tual et al., 2013; Valcin et al., 2007). Indeed, prior studies suggest that intoxicated individuals with pesticides are at a higher risk of acquiring bronchitis-like disease (Hoppin et al., 2007;Mostafalou and Abdollahi, 2017; Tual et al., 2013; Valcin et al., 2007; Ye et al., 2017). Because exposure to toxic pollutants leads to increased dilation of the respiratory tract (Amara, 2014; Hern´ andez et al., 2011), this phenomenon could enhance susceptibility to infection with a respiratory pathogen, which we numerically explore using our modeling framework. Because certain subpopulations (e.g., agricultural workers) in many areas of the world are at an increased risk to pesticide intoxication, there is a need to better understand their role in the dynamics of respiratory infectious diseases (and possibly other diseases as well) in the context of genetic susceptibility variation. In this paper we designed an SIRS (Susceptible-Infectious- RemovedSusceptible) compartmental model that includes a subpopulation that is intoxicated with pesticides and whose susceptibility depends on genetic variation (three levels) of the population in question. Our analytic results include an expression of the basic reproductive number and its relationship with the basic reproductive numbers of the corresponding subsystems with a single level of genetic susceptibility to the pesticide. We show that the basic reproductive number can be expressed as a weighted average of the subsystem basic reproductive numbers, where the corresponding weights are the genotypic proportions of the population under study. Our modeling framework is illustrated in the context of a respiratory infec4

tious disease with a basic reproductive number below 1.0 at the baseline level in the absence of toxic exposure. We then illustrate how this infectious disease has the potential to spread in certain populations exposed to the toxic with a realistic genetic makeup based on data published in prior studies (Chia et al., 2009; Hofmann et al., 2009; Holland et al., 2006; L´ opez-Flores et al., 2009; M´arquez-Llanos, 2015; Singh et al., 2011; Z´ un ˜iga-Venegas et al., 2015) while epidemic potential is not enhanced in others. We further illustrate the role of genetic variation in the transmission potential of the disease in a population with ethnic variability.

2. Model 2.1. Definition of variables and parameters The model incorporates three epidemiological states related to the infectious disease, which are: susceptible (S), infected (I) y removed (R). For each epidemiological state, individuals are further classified in two conditions: intoxicated (P) and not intoxicated (D) by pesticides. Finally, genetic susceptibility to the toxic substance will have three levels: low (L), medium (M) and high (H). Therefore, this model incorporates a total of 18 compartments to keep track of the different individual stages. The notation of the compartments is shown in Table 1 Susceptible

Infected

Removed

Not Intoxicated

SDX

IDX

RDX

Intoxicated

SP X

IP X

RP X

Table 1: Notation of compartments of the model for X ∈ {L, M, H}.

Two contagion rates are considered, one for the intoxicated individuals (βP ) and another for the not intoxicated individuals (βD ) where βD < βP , since the model explores this possibility. There is no difference between recovery rates with respect to the disease, so in all cases, it will be denoted by γ. The rate α

5

models the transition from the removed compartment back to the susceptible compartment. Intoxication rates will be represented by µX (X ∈ {L, M, H}) where µL < µM < µH , similarly, detoxification rates are presented by ωX , with ωL > ωM > ωH . Finally, the entry rate according to the population size (N ) will be b(N ) = b which will be distributed according to the genotypic proportions P λX (with λX = 1) and exit rate by d, with d = b. It is important to note that the model is relatively simple as it does not in-

corporate other sources of heterogeneity such as age or gender. Moreover, we do not consider the simultaneous role of other neurotoxicants. On the other hand, all individuals (agricultural workers) enter the population without intoxication and free of disease. We model equal rates of entry (b) and exit (d) of agricultural workers (study population). Note that deaths, resignations and layoffs are part of the exit rate. Finally, a normalized population will be assumed (N = 1). 2.2. Flow chart The dynamics of the model is represented by the following flow chart:

Figure 1: Dynamics of the states in which an individual can be found. X, Y ∈ {L, M, H} and A ∈ {D, P }

Figure 1 shows the dynamics in which an individual (agricultural worker) could be exposed to pesticides. The entries to system are classified according to the classes of genetic susceptibility to the toxic substance with respect to 6

the susceptibility of the disease and without intoxication, which is reflected in the distribution λX for each one of them. Exits to system can occur in any compartment. The contagion rate of the disease differentiates intoxicated and not intoxicated individuals, incorporating the complexity of contagion in the values of β (βD < βP ), while the recovery (γ) is assumed constant as it is not influenced by the genetic susceptibility to the pesticide (Ferrer, 2003). It is important to highlight that the transfer from intoxicated to non-intoxicated and vice versa can occur in all compartments corresponding to the disease (susceptible, infected, removed). In addition, transmission of the disease stems from the interaction between susceptible individuals and all infectious individuals. For example, an individual with moderate genetic susceptibility to the toxic substance that begins to work in the agricultural sector belongs to the SDM state where, he/she can remain in his/her condition, or (i) get intoxicated or (ii) become infected. In the case of (i) it goes to state SP M while in the case of (ii) it goes to state IDM . In this last case (ii) the individual can be infected by another infectious individual belonging to any of the six possible infectious states while the individual maintains its intoxication condition βD . In the case of (ii), there are two possibilities, (ii.1) the individual recovers from the disease so it goes to the state RDM or it becomes intoxicated (ii.2) arriving at the state IP M , which (ii.2) or recovering by going to the state RP M or is detoxified by returning to state IDM . Finally, recovered individuals lose immunity and become susceptible to the disease, whether intoxicated (SP M ) or not (SDM ), depending on whether it was in the RP M or RDM state respectively.

7

2.3. System of differential equations The mathematical model is given by the following system of differential equations:  X  IAY − (µX + d)SDX + ωX SP X + αRDX + λX b  S˙ DX = −βD SDX    A,Y    X    S˙ P X = −βP SP X IAY − (ωX + d)SP X + µX SDX + αRP X     A,Y   X   I˙ = +βD SDX IAY + ωX IP X − (µX + γ + d)IDX DX A,Y   X   ˙P X = +βP SP X  I IAY + µX IDX − (ωX + γ + d)IP X     A,Y      R˙ DX = +γIDX + ωX RP X − (µX + α + d)RDX      ˙ RP X = +γIP X + µX RDX − (ωX + α + d)RP X , (1) where A ∈ {D, P } y X, Y ∈ {L, M, H}.

This model preserves the genetic pool, since the distribution of entries, which is given by λX =

X C,A

CAX ,

C ∈ {S, I, R},

A ∈ {D, P },

X ∈ {L, M, H},

remains constant. In fact, summing all the equations of system (1), the sum of the C˙ AX equals zero, therefore λ0X = 0. 3. Basic reproductive number To derive the expression of the basic reproductive number, R0 (see Appendix A.1), we employed the classical method of the next generation matrix (Diekmann et al., 2000; Van den Driessche and Watmought, 2002). We obtained the following R0 expression for our model: R0 =

∗ ∗ ∗ βD (SDL + SDM + SDH ) + βP (SP∗ L + SP∗ M + SP∗ H ) , γ+d

8

(2)

λX (ωX + d) λX µX and SP∗ X = , X ∈ {L, M, H}, are µX + ωX + d µX + ωX + d the susceptibility components of the infection-free equilibrium. ∗ where SDX =

It is natural to determine the basic reproductive number when the labor population has a single level of susceptibility to the toxic. Accordingly, if we assume that in this population all individuals have a low susceptibility level towards the toxic substance (L), the system (1) is reduced to:                                     

S˙ DL S˙ P L I˙DL

= −βD SDL = −βP SP L =

+βD SDL

X A

X A

X A

I˙P L

=

+βP SP L

X A

IAL − (µL + d)SDL + ωL SP L + αRDL + b IAL − (ωL + d)SP L + µL SDL + αRP L IAL + ωL IP L − (µL + γ + d)IDL

(3)

IAL + µL IDL − (ωL + γ + d)IP L

R˙ DL

=

R˙ P L

+γIDL + ωL RP L − (µL + α + d)RDL

=

+γIP L + µL RDL − (ωL + α + d)RP L .

We note that in this case λL = 1, which explains why this value does not appear in the system. The respective basic reproductive number (RL 0 ) is equal to (see Appendix A.2): RL 0 = where

∗ SDL =

ωL + d µL + ωL + d

∗ βD SDL + βP SP∗ L , γ+d

and SP∗ L =

(4)

µL . µL + ωL + d

Similarly, the basic reproductive number is obtained for cases in which the population only includes individuals with a medium susceptibility level towards H the toxic (RM 0 ) or only high (R0 ), which is given by:

RM 0 =

∗ βD SDM + βP SP∗ M , γ+d

9

(5)

where

∗ SDM =

ωM + d µM + ωM + d

and SP∗ M =

RH 0 = where

∗ SDH =

ωH + d µH + ωH + d

µM , or µM + ωM + d

∗ βD SDH + βP SP∗ H , γ+d

and SP∗ H =

(6)

µH , respectively. µH + ωH + d

Therefore, we can deduce from (4)(5) and (6) that (2) can be written as: M H R0 = λL RL 0 + λM R0 + λH R0 .

(7)

That is, the basic reproductive number is given by a linear combination of the basic reproductive numbers of the subsystems corresponding to the variation in the genetic susceptibility to the pesticide. 3.1. Basic reproductive number in relation to the classic basic reproductive number of a SIRS model For the classic SIRS models (Britton, 2003; van den Driessche, 2017), the respective basic reproductive number (denoted here by Rc0 ) is given by Rc0 = βD /(γ + d), remember that βD is the not intoxicated contagion rate. So, it is interesting to represent the basic reproductive number obtained from our model (R0 ) in terms of the Rc0 , in order to show how Rc0 is affected by the presence of pesticides and levels genetic susceptibility to the toxic substance. For this, we note that expression (7) can be written as follows: R0 = Rc0 (λL TL + λM TM + λH TH ), where TX =

(8)

(βP /βD )µX + ωX + d , for X ∈ {L, M, H}. µX + ωX + d

It is known that βP /βD is always greater than one, so TX > 1 for all X, so

we check the expected, that R0 > Rc0 , and with no pesticides (which is not our

case) is obtained R0 = Rc0 . From this mode, if Rc0 > 1, always R0 > 1. Given a respiratory infectious disease with basic reproductive number less than one in 10

a population not exposed to the toxic (Rc0 < 1), it is interesting to determine the conditions under which the disease can be installed in a population exposed to pesticides (agricultural workers distributed according to levels of genetic susceptibility to a pesticide) at an epidemic level (R0 > 1). Which from (8) you get that if Rc0 < 1, then R0 > 1 if and only if X X

λX TX >

γ+d . βD

(9)

It is also interesting to quantify how the basic reproductive number increases in the presence of pesticides. From (8) and (9) it is deduced that the increase c of the classic basic reproductive number (R ! 0 ) to the basic reproductive number X in our model (R0 ) is of λX TX − 1 . Therefore, the percentage increase

is given by

X

100

X X

!

λX TX − 1 %.

(10)

4. Basic reproductive number behavior in relation to genotypic proportions Based on PON1 biomarker, studies have demonstrated significant geographic variation in the distribution of genetic susceptibility to the pesticide (Hofmann et al., 2009; Holland et al., 2006; L´ opez-Flores et al., 2009; M´ arquez-Llanos, 2015; Singh et al., 2011; Z´ un ˜iga-Venegas et al., 2015) (see table 2). These studies indicate that in Latin America, individuals are more susceptible to Chlorpyrifos, which is one of the OP pesticides most frequently used in the region (L´opez-G´ alvez et al., 2018; Mu˜ noz-Quezada et al., 2017; Z´ un ˜iga-Venegas et al., 2015). Indeed, this is explained by the corresponding distribution of genetic susceptibility to the toxic substance which is mostly centered between the QQ and QR genotypes. In addition, even within a single country it is possible to observe considerable variation according to ethnic groups (Chia et al., 2009) (see Table 2). From (7) we can observe how the distributions of genotypic proportions (λL , λM , λH ) affect the final value of the basic reproductive number. So we 11

Country

QQ (%)

QR (%)

RR (%)

Reference

Chile

39

40

21

(Z´ un ˜iga-Venegas et al., 2015)

Colombia

12

88

0

(M´ arquez-Llanos, 2015)

Mexico

16

47

27

(L´ opez-Flores et al., 2009)

United States

30

50

20

(Hofmann et al., 2009)

India

36

51

13

(Singh et al., 2011)

4

64

32

Singapore : - Indian

43

53

4

- Malay

13

58

29

- Chinese

(Chia et al., 2009)

Table 2: Genotypic percentage distribution, by means of PON1, according to the country or ethnic group.

will proceed to study the dynamics of the basic reproductive number, setting all the values of the parameters, except the λ’s. The baseline parameter values are shown in Table 3. It is worth mentioning that all the rates in our model have units of days−1 . The rates corresponding to the infectious respiratory disease (βD and γ) and the exit rate from the system (d), were calibrated to obtain a baseline basic reproductive number < 1.0 in the absence of pesticide exposure (Rc0 = βD /(γ+d) ≈ 0.8703). Hence, the value assigned to βD expresses that, approximately, for each encounter between a susceptible individual and an infected one, there is a probability of contagion of 1/5. On the other hand, given that 1/γ models the average duration of the infectious period, the value assigned to this rate indicates that an individual (agricultural worker) remains infectious for four days on average. The rates of entry (b) and exit (d) are interpreted based on labor turnover (b = d). In particular, their values indicate that an agricultural worker, on average, remains at work for one year (an estimated value that does not escape reality). The value of βP was set to be slightly higher than the value of βD since that the model explores this possibility. In order to illustrate graphically that R0 can (for some genetic distribution profiles) be greater than one, we assume a high intoxication rate for the agricultural 12

population exposed to pesticides. In this sense, we have assigned the value for µM (the intermediate one), and for the rates µL and µH we assume values 2/3 lower and 2/3 higher, respectively. Regarding the detoxification rates (ω’s), the OP takes between 24 and 48 hours to leave the human body according to published studies (Jones et al., 2002; Krieger et al., 2000; Meuling et al.,2005). Thus, ωX ∈ [0.5, 1], so the value chosen for ωM (the intermediate) is within the indicated range. For the other two parameters, ωL and ωH , we assign the minimum and maximum values. We emphasize that the major focus of this work is to illustrate how the transmission potential of the disease is differentially affected across populations exposed to the toxic according to the variation in genetic susceptibility to the pesticide. Symbol

Meaning

Values

βD

Infection rate of non-intoxicated individuals

0.22

βP

Infection rate of intoxicated individuals

0.32

γ

Recuperation rates

0.25

d

Exit rate

0.0028

b

Entry rate

0.0028

µL

Intoxication rate of individuals low GSP

0.10

µM

Intoxication rate of individuals medium GSP

0.30

µH

Intoxication rate of individuals high GSP

0.50

ωL

Detoxification rate of individuals low GSP

1.00

ωM

Detoxification rate of individuals medium GSP

0.67

ωH

Detoxification rate of individuals high GSP

0.50

Table 3: Parameter definitions and baseline values employed in simulations. GSP corresponds to genetic susceptibility to the pesticide.

To better illustrate the system dynamics, Figure 2 displays the relationship between the fixed values of RX 0 (with X ∈ {L, M, H}) and the location obtained after different values or conditions of the λ’s. These different values (λ’s), establish the division of the triangle, separating it in R0 greater or smaller than one (see Figure 2). It is worth mentioning that for this analysis, the OP Chlor13

pyrifos was considered, so QQ, QR and RR correspond to the high, medium and low genetic susceptibility to the pesticide respectively.

λHRH 0

1

0.5

R0>1 R0<1

0 1

1

0.5 λLRL0

0.5 0

0

λMRM 0

Figure 2: Value of basic reproductive number (greater or smaller than one) with respect to distribution of the susceptible proportions (λ’s). The triangle corresponds to the graph of the plane 1 =

x RL 0

+

y RM 0

+

z RH 0

M H , where given λ’s , x = λL RL 0 , y = λM R0 and z = λH R0 .

According to the data given in Table 2, given our baseline respiratory infectious disease with reproductive number less than one, in a population exposed to pesticides, results suggest that in countries such as Chile (R0 = 1.0035), Colombia (R0 = 1.0013) and India (R0 = 1.0081), epidemic transmission in agricultural workers exposed to pesticides is plausible. On the other hand, in countries such as Mexico (R0 = 0.8818) or the United States (R0 = 0.9976) our results suggest that the disease is still not able to reach epidemic transmission potential (see Figure 3(a)). In fact, within the same country, like Singapore, among the three ethnic entities studied, those corresponding to Chinese (R0 = 0.9677) and Malays (R0 = 0.9770) do not present an epidemic, but the Indians (R0 = 1.0212) do (see Figure 3(b)).

14

Chile Colombia Mexico USA India H

1

0.5

λHR0

λHRH 0

1

Chinese Indian Malay

R0>1

0.5

R0>1

R <1

R0<1

0

0 1

1

0.5 L λLR0

0 1

1

0.5

0.5 0 0

L λLR0

λ RM M 0

(a) Country

0.5 0 0

λ RM M 0

(b) Ethnicity in Singapore

Figure 3: Distribution of (a) countries and (b) ethnic groups, according to their genetic susceptibility to the toxic in relation to the value of the basic reproductive number with respect to the same respiratory infectious disease.

These dynamics emphasize the importance of the population distribution with respect to the genetic susceptibility to the pesticide as a function of geography and ethnicity. Another interesting question stems from the properties of the pathogen whereby different genetic profiles for different human populations affect its propagation. Hence, we investigate the conditions under which R0 is greater or less than one with respect to the relationship between the susceptibleintoxicated state (SP∗ ) and the difference in contagion rates (which we will denote by ). From equation (2) it is obtained that R0 =

βD ∗ βP ∗ S + S , γ+d D γ+d P

(11)

∗ ∗ ∗ ∗ where SD = SDL + SDM + SDH and SP∗ = SP∗ L + SP∗ M + SP∗ H . In addition,

since βP > βD , we can establish that βP = βD + , with  > 0. So (11) can be expressed by R0 = Rc0 +

 S∗ . γ+d P

(12)

Therefore, R0 = 1 if and only if SP∗ =

1 (1 − Rc0 )(γ + d), 

15

(13)

where SP∗ =

X X

λX µX . µX + ωX + d 1

0.8

S*P

0.6 R0 > 1

0.4 0.2 0

R0 < 1 0

0.5

1

1.5

2

Figure 4: Relationship between increased rate of intoxication () infection and the amount of ∗ ), where ϑ ≈ 0.0328. intoxicated population free of infection (SP

From Figure 4 and the context of the values given in Table 3, we can infer that if the difference between transmission rates is less than 0.0328, then the disease will not spread irrespective of the size of the intoxicated population. Moreover, as this value increases () the fraction of the population intoxicated that is free of infection, necessary so that R0 > 1, decreases rapidly. 5. Discussion and conclusion In this paper, we employed mathematical analysis and a relatively simple structured compartmental model to illustrate the potential effects of the genotypic distribution with respect to the Q and R alleles (corresponding to QQ, QR and RR) on the transmission potential of an infectious respiratory illness. This effect further depends on other factors such as geography and type of pesticide. For example, in Latin America, the most commonly used OP is Chlorpyrifos (L´ opez-G´ alvez et al., 2018; Mu˜ noz-Quezada et al., 2017; Z´ un ˜iga-Venegas et al., 2015), so individuals with a QQ genotype exhibit higher susceptibility to the toxic while those with an RR genotype have low susceptibility. Based on their 16

distribution in the population, countries such as Chile and Colombia my be more susceptible to pesticide intoxication while this may not be the case for other countries such as Mexico. We note that the mathematical model introduced in this paper accommodates high, medium and low susceptibility levels to pesticide intoxication, so it can be applied to any genotypic distribution with respect to the Q and R alleles since this classification depends on the pesticide under study. The proposed model could be expanded to incorporate the cost variable in public health, casting an optimization problem, since the use of pesticides and pesticide intoxication pose significant economic impacts (Ram´ırez-Santana et al., 2014; Soares and de Souza Porto, 2012). Our modeling framework is likely to find more applications as it represents a baseline framework that is useful to generate approximate assessments of the risk posed by pesticides on the transmission dynamics of infectious diseases. Hence, our work may have public health implications. A straightforward intervention will prohibit certain pesticides particularly those for which the population exhibits higher susceptibility such as the already mentioned Chlorpyrifos. Moreover, our models could be refined and tailored to specific populations. In particular, our results indicate that the public health policy will depend on various factors including the pesticide type and the distribution of the population according to the corresponding genetic susceptibility to the pesticide. In summary, our modeling framework, analysis and numerical results could help monitor the potential effects of certain pesticides on the spread of respiratory infectious disease. Moreover, because the distribution of the population according to its genetic susceptibility to the toxic substance has an impact on the transmission potential of a respiratory infectious disease, other studies connecting pesticide intoxication and respiratory infectious disease should consider the effects of the underlying genetic distribution of the population.

17

Acknowledgements This research was supported by the Doctorado en Modelamiento Matem´ atico Aplicado (DM2 A-UCM), a program belonging to the Universidad Cat´ olica del Maule, and School of Public Health corresponding to Georgia State University.

Funding This work was funded by a research project of the National Commission on Science and Technology (CONICYT) of Chile, FONDECYT of Initiation #11150784.

Appendix A. To calculate the basic reproductive numbers, techniques and notation given by Van den Driessche and Watmought (2002) will be used: A.1. BRN (R0 ) The total population is defined as follows: N = SDL + SDM + SDH + SP L + SP M + SP H + IDL + IDM + IDH + IP L + IP M + IP H + RDL + RDM + RDH + RP L + RP M + RP H .

In order to find the BRN, the

equation (1) need to be arranged so that the infection appears in the first m = 6 equations.

By rearranging the system is obtained:

(IDL , IDM , IDH , IP L , IP M , IP H , SDL , SDM , SDH , SP L , SP M , SP H , RDL , RDM , RDH , RP L , RP M , RP H )t . As a result the disease free equilibrium is given by: ∗ ∗ ∗ X ∗ = (0, 0, 0, 0, 0, 0, SDL , SDM , SDH , SP∗ L , SP∗ M , SP∗ H , 0, 0, 0, 0, 0, 0), where

∗ SDX =

λX (ωX + d) , µX + ωX + d

SP∗ X =

λX µX , µX + ωX + d

whit X ∈ {L, M, H}.

The uniqueness of the solution can be computed from the determinant of the system matrix.

18

The matrices F and  β S∗  D DL     ∗ βD SDM      ∗  βD SDH   F =    βP SP∗ L       βP SP∗ M     βP SP∗ H



V −1

V −1 evaluated in X ∗ are: ∗ βD SDL

∗ βD SDL

∗ βD SDL

∗ βD SDL

∗ βD SDL

∗ βD SDM

∗ βD SDM

∗ βD SDM

∗ βD SDM

∗ βD SDM

∗ βD SDH

∗ βD SDH

∗ βD SDH

∗ βD SDH

∗ βD SDH

βP SP∗ L

βP SP∗ L

βP SP∗ L

βP SP∗ L

βP SP∗ L

βP SP∗ M

βP SP∗ M

βP SP∗ M

βP SP∗ M

βP SP∗ M

βP SP∗ H

βP SP∗ H

βP SP∗ H

βP SP∗ H

βP SP∗ H

ωL + γ + d  e      0        0   1  =  γ+d  µL   e      0       0



              ,             

0

0

ωL e

0

ωM + γ + d f

0

0

ωM f

0

ωH + γ + d g

0

0

0

0

ωL + γ + d e

0

µM f

0

0

ωM + γ + d f

0

µH g

0

0

where e = µL + ωL + γ + d f = µM + ωM + γ + d y g = µH + ωH + γ + d. The next generation matrix can be shown by F V −1 . It defines R0 = ρ(F V −1 ) 19

0



      0       ωH   g   ,    0       0      ωH + γ + d  g

where ρ is the spectral radius of the matrix: F V −1 =

1 ·F . γ+d

Thus, the basic reproductive number is given by:

R0 =

∗ ∗ ∗ βD (SDL + SDM + SDH ) + βP (SP∗ L + SP∗ M + SP∗ H ) . γ+d

A.2. Basic reproductive number of individuals with low susceptibility level The total population is defined as follows: N = SDL + SP L + IDL + IP L + RDL + RP L . In order to find the basic reproductive number, the equation (3) need to be arranged so that the infection appears in the first m = 2 equations. By rearranging the system is obtained: (IDL , IP L , SDL , SP L , RDL , RP L )t . As a ∗ , SP∗ L , 0, 0), where result the disease free equilibrium is given by: X ∗ = (0, 0, SDL

∗ = SDL

ωL + d , µL + ωL + d

SP∗ L =

µL . µL + ωL + d

The uniqueness of the solution can be computed from the determinant of the system matrix. The matrices F and V −1 evaluated in X ∗ are:  ∗ βD SDL β S∗  D DL  F =  βP SP∗ L βP SP∗ L

V −1 =



   (γ + d)(µL + γL + d)  1



  , 

ωL + γ + d

µL

ωL



  .  ωH + γ + d

The next generation matrix can be shown by F V −1 . It defines R0 = ρ(F V −1 ) where ρ is the spectral radius of the matrix:

20

F V −1

 ∗ βD SDL 1   =  γ+d βP SP∗ L

∗ βD SDL

βP SP∗ L

Thus, the basic reproductive number is given by: RL 0 =



  , 

∗ βD SDL + βP SP∗ L . γ+d

References [1] Amara, A. (2014). Pesticides and asthma: Challenges for epidemiology. Frontiers Public Health, 2. [2] Bakhsh, K., Ahmad, N., Kamran, M. A., and et. al. (2016). Occupational hazards and health cost of women cotton pickers in pakistani punjab. BMC Public Health, 16:961. [3] Britton, N. F. (2003). Essential mathematical biology. Springer undergraduate mathematics series. [4] Chia, S. E., Ali, S. M., Ericc, Y. P. H., Gan, L., and et al. (2009). Distribution of pon1 polymorphismspon1q192r and pon1l55m among chinese, malay and indian males in singapore and possible susceptibility to organophosphate exposure. NeuroToxicology, 30:214–219. [5] Costa, L., Cole, T., Jansen, K., and Furlong, C. (2000). Paraoxonase (pon1) and organophosphate toxicity. Dordrecht, 6:209–220. [6] Costa, L., Cole, T., Vitalone, A., and Furlong, C. (2004). Measurement of paraoxonase (pon1) status as a potential biomarker of susceptibility to organophosphate toxicity. International journal of clinical chemistry, 352:37– 47. [7] Costa, L., Giordano, G., Cole, T., Marsillach, J., and Furlong, C. (2013). Paraoxonase 1 (pon1) as a genetic determinant of susceptibility to organophosphate toxicity. Toxicology, 307:115–122. 21

[8] Damalas, C. and Koutroubas, S. (2016). Farmers exposure to pesticides: Toxicity types and ways of prevention. Toxics, 4:1–10. [9] Diekmann, O., Heesterbeek, H., and Britton, T. (2000). Mathematicals tools for understanding infectious diseases dynamics. Princeton series in theoretical and computational biology. [10] Dutta, S. and Bahadur, M. (2018). Effect of pesticide exposure on the cholinesterase activity of the occupationally exposed tea garden workers of northern part of west bengal, india. Biomarkers, 78:1–16. [11] FAO (2003). International code of conduct on the distribution and use of pesticides. Food and Agriculture Organization of the United Nations. [12] Ferrer, A. (2003). Pesticide poisoning. Anales del Sistema Sanitario de Navarra, 26:155–171. [13] Guti´errez-Jara, J. P., C´ ordova-Lepe, F., and Mu˜ noz-Quezada, M. T. (2019). Dynamics between infectious diseases with two susceptibility conditions: A mathematical model. Mathematical Biosciences, 309:66–77. [14] Hern´ andez, A., Parr´ on, T., and Alarc´ on, R. (2011). Pesticides and asthma. Current Opinion in Allergy and Clinical Immunology, 11:90–96. [15] Herrera, M. Q. and Ticona, S. Q. (2011). [paraoxonase an inducible enzyme: prevents the risk of cardiovascular disease]. SCIENTIFICA, La Paz, 9:46–50. [16] Hickson, R. and Roberts, M. (2014). How population heterogeneity in susceptibility and infectivity influences epidemic dynamics. Journal of Theoretical Biology, 350:119–124. [17] Hofmann, J., Keifer, M., Checkoway, H., Roos, A. D., and et al. (2010). Biomarkers of sensitivity and exposure in washington state pesticide handlers. Adv Exp Med Biol., 660:19–27. [18] Hofmann, J., Keifer, M., Furlong, C., De Roos, A., and et al. (2009). Serum cholinesterase inhibition in relation to paraoxonase-1 (pon1) status 22

among organophosphate-exposed agricultural pesticide handlers. Environmental Health Perspectives, 117:1402–1408. [19] Holland, N., Furlong, C., Bastaki, M., Richter, R., and et al. (2006). Paraoxonase polymorphisms, haplotypes, and enzyme activity in latino mothers and newborns. Environmental Health Perspectives, 114:985–991. [20] Hoppin, J., Valcin, M., Henneberger, P., Kullman, G., Umbach, D., London, S., Alavanja, M., and Sandler, D. (2007). Pesticide use and chronic bronchitis among farmers in the agricultural health study. American Journal of Industrial Medicine, 50:969–979. [21] Jones, S. G. K., Mason, H., and Cocker, J. (2002).

Exposure to the

organophosphate diazinon: data from a human volunteer study with oral and dermal doses. Toxicology Letters, 134:105–113. [22] Krieger, R., Bernard, C., Dinoff, T., Fell, L., Osimitz, T., Ross, J., and Thongsinthusak, T. (2000). Biomonitoring and whole body cotton dosimetry to estimate potential human dermal exposure to semivolatile chemicals. Journal of Exposure Science & Environmental Epidemiologyvolume, 10:50–57. [23] L´ opez-Flores, I., Lacasa˜ na, M., Blanco-Mu˜ noz, J., Aguilar-Gardu˜ no, C., and et al. (2009). Relationship between human paraoxonase-1 activity and pon1 polymorphisms in mexican workers exposed to organophosphate pesticides. Toxicology Letter, 188:84–90. [24] L´ opez-G´ alvez, N., Wagoner, R., Beamer, P., de Zapien, J., and Rosales, C. (2018). Migrant farmworkers exposure to pesticides in sonora, mexico. International Journal of Environmental Research and Public Health, 15:2651. [25] Mackness, B., Durrington, P., Povey, A. C., Thomson, S., Dippnall, M., Mackness, M., Smith, T., and Cherry, N. (2003). Paraoxonase and susceptibility to organophosphorus poisoning in farmers dipping sheep. Pharmacogenetics, 13:81–88.

23

[26] Marcelino, A. F., Wachtel, C. C., and de Castilhos Ghisi, N. (2019). Are our farm workers in danger? genetic damage in farmers exposed to pesticides. Int J Environ Res Public Health., 16:358. [27] M´ arquez-Llanos, D. C. (2015). [characterization of the enzymatic activity and genetic polymorphisms of paraoxonase-1 (pon-1), in workers exposed to organophosphorus pesticides in the municipality of soacha 2014]. Universidad Nacional de Colombia, Tesis de investigaci´ on:–. [28] Meuling, W., Ravensberg, L., Roza, L., and van Hemmen, J. (2005). Dermal absorption of chlorpyrifos in human volunteers. International Archives of Occupational and Environmental Health, 78:44–50. [29] Mostafalou, S. and Abdollahi, M. (2017). Pesticides: an update of human exposure and toxicity. Arch Toxicol., 91:549–599. [30] Mu˜ noz-Quezada, M. T., Iglesias, V., Lucero, B., and et. al. (2013). Neurodevelopmental effects in children associated with exposure to organophosphate pesticides: A systematic review. Neurotoxicology, 39:158–168. [31] Mu˜ noz-Quezada, M. T., Iglesias, V., Lucero, B., and et. al. (2016). [Organophosphate pesticides and neuropsychological and motor effects in the Maule Region, Chile]. Gac Sanit., 30:227–231. [32] Mu˜ noz-Quezada, M. T., Lucero, B., Iglesias, V., and et. al. (2017). Exposure to organophosphate (op) pesticides and health conditions in agricultural and non agricultural workers from maule, chile. International Journal of Environmental Health Research, 27:82–93. [33] Povey, A. C. (2010). Geneenvironmental interactions and organophosphate toxicity. Toxicology, 278:294–304. [34] Ram´ırez-Santana, M., Iglesias-Guerrero, J., Castillo-Riquelme, M., and Scheepers, P. (2014). Assessment of health care and economic costs due to episodes of acute pesticide intoxication in workers of rural areas of the coquimbo region, chile. Value in health issues, 5:35–39. 24

[35] Singh, S., Kumar, V., Thakur, S., Banerjee, B. D., and et al. (2011). Paraoxonase-1 genetic polymorphisms and susceptibility to dna damage in workers occupationally exposed to organophosphate pesticides. Toxicology and Applied Pharmacology, 252:130–137. [36] Soares, W. L. and de Souza Porto, M. F. (2012). Pesticide use and economic impacts on health. Sa´ ude P´ ublica, 46:209–217. [37] Tual, S., Clin, B., Levque-Morlais, N., Raherison, C., Baldi, I., and Lebailly, P. (2013). Agricultural exposures and chronic bronchitis: findings from the agrican (agriculture and cancer) cohort. Annals of Epidemiology, 23:539–545. [38] Valcin, M., Henneberger, P., Kullman, G., Umbach, D., London, S., Alavanja, M., Sandler, D., and Hoppin, J. (2007). Chronic bronchitis among nonsmoking farm women in the agricultural health study. Journal of Occupational and Environmental Medicine, 49:574–583. [39] van den Driessche, P. (2017). Reproduction numbers of infectious disease models. Infectious Disease Modelling, 2:288–303. [40] Van den Driessche, P. and Watmought, J. (2002). Reproduction numbers and subthreshold endemic equilibria for compartmental models of disease transmission. Mathematical Biosciences, 180:29–48. [41] Ye, M., Beach, J., Martin, J., and Senthilselvan, A. (2017). Pesticide exposures and respiratory health in general populations. journal of environmental science, 51:361–367. [42] Z´ un ˜iga-Venegas, L., Aquea, G., Taborda, M., Bernal, G., and et al. (2015). Determination of the genotype and phenotype of serum paraoxonase 1 (pon1) status in a group of agricultural and nonagricultural workers in the coquimbo region, chile. Journal of Toxicology and Environmental Health, Part A: Current Issues, 78:357–368.

25

Author Contributions Section. Juan Pablo Gutiérrez: Conceptualization, Methodology, Software, Validation, Formal analysis, Investigation, Writing – Original Draft, Writing – Review & Editing, Visualization, Supervision, Project administration. Fernando Córdova: Conceptualization, Methodology, Formal analysis, Investigation, Writing – Review & Editing, Supervision. María Teresa Muñoz: Conceptualization, Methodology, Validation, Investigation, Writing – Review & Editing, Visualization, Supervision. Gerardo Chowell: Conceptualization, Methodology, Investigation, Writing – Review & Editing, Supervision.