Laureen Ribassin-Majed

Rachid Lounes

Stephan Clemençon

Deterministic modeling for transmission of Human Papillomavirus 6/11: impact of vaccination

Keywords: equilibrium, stability, reproduction number, vaccination

This paper is devoted to assess the impact of quadrivalent Human Papillomavirus (HPV) vaccine on prevalence of non-oncogenic HPV 6/11 types in French males and females. For this purpose, a non-linear dynamic model of heterosexual transmission for HPV 6/11 types infection is developed, which accounts for immunity due to vaccination in particular. The vaccinated reproduction number R v is derived using the approach described by Diekmann (2010) called the Next Generation Operator approach. The model proposed is analyzed, with regard to existence and uniqueness of the solution, steady-state stability. Precisely, the stability of the model is investigated depending on the sign of R v -1. Prevalence data are used to fit a numerical HPV model, so as to assess infection rates. Our approach suggests that 10 years after introducting vaccination, the prevalence of HPV 6/11 types in females will be halved and that in males will be reduced by one quarter, assuming a sustained vaccine coverage of 30% among females. Using the formula we derived for the vaccinated reproduction number, we show that the non-oncogenic HPV 6/11 types would be eradicated if vaccine coverage in females is kept above 12%.

Human Papillomavirus, deterministic epidemic model,

Introduction

Human Papillomavirus (HPV) is the most common sexually transmitted infection. At least 70 per cent of sexually active people acquire HPV infection at some point in their lives [START_REF] Syrjänen | Prevalence, incidence, and estimated life-time risk of cervical human papillomavirus infections in a nonselected Finnish female population[END_REF]. Nearly one hundred HPV genotypes have been identified, among which there are low risk genotypes, causing benign anogenital lesions, and high risk genotypes, which induce pre-cancerous lesions in the cervix. Epidemiological studies on HPV infections establish the role of these viruses as the primary cause of cervical cancer [START_REF] Muñoz | Human papillomavirus and cancer: the epidemiological evidence[END_REF]. These infections are also the cause of anogenital cancers, head and neck cancers, anogenital warts and recurrent respiratory papillomatosis among women and men. While HPV 16/18 are incriminated in 70% of cervical cancer, HPV 6/11 are the primary cause of almost (90%) all genital warts and of most respiratory recurrent papillomatosis cases (Gissman et al., 1983). While available epidemiological data in France indicate a prevalence of 1% in the general population for genital warts, a prevalence of 10% has been observed in young individuals aged 15 to 25 [START_REF] Monsonego | Traité des infections et pathologies génitales à papillomavirus[END_REF]. Efficacy of curative treatment is limited in presence of a high recurrence rate for genital warts. Two prophylactic vaccines against HPV infections are available. The bivalent vaccine protects individuals from oncogenic HPV 16/18 types. The quadrivalent vaccine protects individuals from oncogenic HPV types 16/18 and non-oncogenic HPV types 6/11. The purpose of this paper is to assess the impact of the quadrivalent HPV vaccine on the prevalence of non-oncogenic HPV 6/11 types in French individuals.

INTRODUCTION

2

Mathematical epidemic modeling provides useful tools for analyzing the spread and control of infectious diseases [START_REF] Hethcote | The Mathematics of Infectious Diseases[END_REF][START_REF] Brauer | Mathematical Epidemiology[END_REF]. In particular, it can be used in order to assess the impact of vaccination. Numerous mathematical models have been introduced in the literature to study epidemics of communicable diseases such as measles, influenza, rubeola, and chickenpox. References are much too numerous to be listed exhaustively, see Nuno et al.(2005) or [START_REF] Feng | A model for tuberculosis with exogenous reinfection[END_REF] for instance for recent accounts of (deterministic) epidemic models.

Since anti-HPV vaccines have been developed, several deterministic or hybrid models have been developed to assess the potential impact of vaccination on HPV prevalence and linked diseases (see e.g. [START_REF] Hughes | The theoretical population-level impact of a prophylactic human papilloma virus vaccine[END_REF]; [START_REF] Barnabas | Epidemiology of HPV 16 and cervical cancer in Finland and the potential impact of vaccination: mathematical modelling analyses[END_REF]; [START_REF] Elbasha | Model for assessing human papillomavirus vaccination strategies[END_REF]; [START_REF] Taira | Evaluating human papillomavirus vaccination programs[END_REF]; Ribassin-Majed et al. (2012a,b). However, the models documented in the literature are based on numerical simulations and offer very limited analytical results. It should be noticed in addition that the complexity of these models does not allow to study local and global stability for the equilibrium points. To the best of our knowledge, only Elbasha et al. studied local and global stability for equilibrium solutions in a simple Susceptible-Infected-Retired (SIR) model [START_REF] Elbasha | Impact of Prophylactic Vaccination against Human Papillomavirus Infection, Contemporary Mathematics 410, Mathematical Studies on Human Disease Dynamics; Emerging Paradigms and Challenges[END_REF][START_REF] Elbasha | Global stability of Equilibria in a Two-sex HPV Vaccination Model[END_REF].

In certain models developed for communicable diseases, the system of differential equations that describes the evolution of the epidemics is not globally asymptotically stable, multiple equilibrium points coexisting [START_REF] Feng | A model for tuberculosis with exogenous reinfection[END_REF][START_REF] Nuño | Dynamics of two-strain influenza with isolation and partial cross-immunity[END_REF]. When considering vaccination, the vaccinated reproduction number R v , in addition to the basic reproduction number R 0 , is estimated to assess the spread of the disease of interest. The basic reproduction number R 0 is a threshold quantity establishing whether an epidemic is likely to spread out or not. It is defined as the expected number of secondary cases of HPV caused by an infected individual during the entire period of infectiousness, in a completely susceptible population [START_REF] Dietz | Transmission and control of arbovirus diseases[END_REF][START_REF] Diekmann | On the definition and the computation of the basic reproduction ratio R 0 in models for infectious diseases in heterogeneous populations[END_REF]. The vaccinated reproduction number represents the threshold in presence of vaccination. These values may determine the likeliest scenario for the evolution of the epidemic disease under study, depending on whether they are below or above the critical value 1. However, bringing the vaccinated reproduction number below 1 may not be sufficient to eradicate endemicity of the disease when multiple locally stable equilibrium solutions coexist [START_REF] Kribs-Zaleta | A simple vaccination model with multiple endemic states[END_REF]. Hence, in dynamic models, investigating the asymptotic behavior is an important issue.

The main objective of this paper is to assess the impact of quadrivalent HPV vaccine on HPV 6/11 prevalences in French females and males. A deterministic model is presented using a system of ordinary differential equations to describe the heterosexual transmission of the virus in the French population, considering the real vaccine coverage observed in France. In addition, a mathematical analysis of the model is provided. Deterministic epidemic modeling for sexually transmitted infections, taking into account sex and vaccination both at the same time, are not well documented in the literature. Human Papillomavirus is the only sexually transmitted infection which can be avoided through vaccination. To the best of our knowledge, only [START_REF] Elbasha | Impact of Prophylactic Vaccination against Human Papillomavirus Infection, Contemporary Mathematics 410, Mathematical Studies on Human Disease Dynamics; Emerging Paradigms and Challenges[END_REF][START_REF] Elbasha | Global stability of Equilibria in a Two-sex HPV Vaccination Model[END_REF] proposed a Susceptible-Infected-Retired (SIR) deterministic model for HPV transmission. This paper is thus the first to investigate specifically the impact of HPV vaccination on non-oncogenic HPV types, through a mathematical/numerical analysis of the dynamic model proposed.

The paper is structured as follows. In section 2, a two-sex model of HPV infection transmission in the sexually active population is introduced. We developed a deterministic model based on a system of ordinary differential equations for HPV 6/11 transmission considering vaccination, to evaluate the effects of vaccination campaigns in France. The basic reproduction number R 0 and vaccinated reproduction number R v are derived from an analytic formula in section 3. The asymptotic behavior of the HPV models (with and without vaccination) is studied in section 3 and the stability of equilibrium solutions of the developed models is illustrated through simulations in section 4. Subsequently, numerical values for the infection rates are assessed by means of a fitting procedure based on prevalence data, and the impact of quadrivalent vaccine on HPV 6/11 prevalence in France for both sexes is quantified in section 5, considering the current vaccine coverage in France. The experimental results are finally discussed in section 6. Technical details are deferred to the Appendix section.

Human Papillomavirus model

In this paper we develop a deterministic model that describes the transmission mechanism of HPV 6/11 infection in a heterosexually active population. Classically, the sexually active population is divided into compartments. Non-oncogenic HPV types do not induce specific natural immunity [START_REF] Monsonego | Traité des infections et pathologies génitales à papillomavirus[END_REF], therefore we use a Susceptible-Infected-Susceptible (SIS) structure. In addition, vaccination is taken into account, see Figure 1.

Figure 1: Flow diagram.

X S X I σ f δ 1 (1 -ϕ f )Λ ϕ f Λ µ f µ f V S V I (1 -τ )σ f δ 1 µ f µ f Y S Y I σ m δ 1 (1 -ϕ m )Λ ϕ m Λ µ m µ m W S W I (1 -τ )σ m δ 1 µ m µ m
Deterministic models are widely used to simulate the spread of sexually transmitted infectious diseases. Such models can be described using a system of ordinary differential equations, with females and males in different compartments. The transmission term is non-linear as it reflects the interaction between infectious and susceptible individuals. The use of a deterministic model allows us here to take into account herd immunity (indirect benefits of vaccination, see Garnett (2005)) in a simple manner, which corresponds to a decrease of HPV infections in non-vaccinated subpopulations due to vaccination coverage of other individuals. When non-vaccinated (respectively vaccinated) women enter in the female sexually active population (of size N f ), they move into the susceptible compartment X S (respectively V S ) at the constant rate [(1ϕ f )Λ] (respectively [ϕ f Λ]) and they leave all compartments at rate µ f . We assume that the exit of the sexually active population balances the entrance in the sexually active population so that the population size N in the model remains constant. In a similar fashion, non-vaccinated (resp. vaccinated) men enter into the male sexually active population (of size N m ) into the susceptible compartment Y S (resp. W S ) at constant rate [(1ϕ m )Λ] (respectively [ϕ m Λ]) and leave all compartments at rate µ. Susceptible individuals are infected with HPV at a per capita rate λ m or λ f (annual rates), depending on their sex. The force of infection depends on infection rates (σ m for men and σ f for women) and on HPV infection prevalence in the opposite sex as well. Then, they move into infected compartments: X I (resp.

V I ) for women, Y I (resp. W I ) for men in non-vaccinated population (resp. in vaccinated population), see the description of parameters on Table 1. The following assumptions shall also be required:

• Vaccinated people can be infected. The degree of protection of vaccine is τ , the relative risk of a vaccinated person experiencing a breakthrough infection is 1 -τ .

• Vaccinated infected individuals are as much infectious as non-vaccinated persons.

• Vaccine immunity does not wane during all sexually active life.

• Women and men who clear HPV infection at rate δ leave infected compartments and go back to susceptible compartments.

Demographic and biological parameters are strictly positive. The notation d dt (•) is used for derivative. The non-linear system of ordinary differential equations that represents this compartmental model is:

dX S dt = (1 -ϕ f )Λ -λ f X S + δX I -µX S dX I dt = λ f X S -(δ + µ)X I dY S dt = (1 -ϕ m )Λ -λ m Y S + δY I -µY S dY I dt = λ m Y S -(δ + µ)Y I dV S dt = ϕ f Λ -(1 -τ )λ f V S + δV I -µV S (2.1) dV I dt = (1 -τ )λ f V S -(δ + µ)V I dW S dt = ϕ m Λ -(1 -τ )λ m W S + δW I -µW S dW I dt = (1 -τ )λ m W S -(δ + µ)W I
We highlight the fact that the system above is nonlinear, the forces of infection depend on infection rates and the prevalences of HPV infection in the opposite sex:

λ f = σ f (Y I + W I ) N m , λ m = σ m (X I + V I ) N f .
The population in the model is assumed to remain constant, that is:

N f = X S + X I + V S + V I , N m = Y S + Y I + W S + W I , N = N f + N m .
And we assume:

N f = N m . Thus N ′ = 2Λ -µN.
Since at equilibrium N * = 2 Λ µ , we only need to analyze the asymptotically autonomous limiting system where N is replaced by its equilibrium value. We consider the system only in the region

D = (X S , X I , Y S , Y I , V S , V I , W S , W I ) ∈ R 8 + , X S + X I + V S + V I = Λ µ = N f and Y S + Y I + W S + W I = Λ µ = N m .
It can be verified that D is positively invariant for this system, which has a unique solution in D. The model is epidemiogically and mathematically well posed. In section 3, the equilibria of the model without vaccination are analyzed and a closed analytical form for the basic reproduction number is given.

Analysis of equilibria and reproduction numbers

In this section, we consider the model without vaccination as a first go. There are two possible equilibria: the Disease Free Equilibrium (DFE in abbreviated form) and the Endemic Equilibrium. In order to analyze the stability of these equilibria, the basic reproduction number R 0 is computed.

The model without vaccination

In the abscence of vaccination, ϕ m = 0 and ϕ f = 0 as well as V S = V I = W S = W I = 0. The system of ordinary differential equations is as follows:

dX S dt = Λ - σ f Y I N m X S + δX I -µX S dX I dt = σ f Y I N m X S -(δ + µ)X I dY S dt = Λ - σ m X I N f Y S + δY I -µY S (3.1) dY I dt = σ m X I N f Y S -(δ + µ)Y I
The equilibria of this model are obtained by setting the right hand sides of the model equations to zero.

The system (3.1) has two equilibria, one at:

P 0 = (X * S , X * I , Y * S , Y * I )= ( Λ µ , 0 , Λ µ , 0
) which is the DFE, and

P 1 = (X * * S , X * * I , Y * * S , Y * * I )
the endemic equilibrium, where:

X * * S = dN (δ + µ) 2σ f , X * * I = Λ µ - dN (δ + µ) 2σ f Y * * S = Λ µ (1 -1/d) + N (δ + µ) 2σ f Y * * I = Λ dµ - N (δ + µ) 2σ f with d = σ m σ f + σ f (δ + µ) σ m σ f + σ m (δ + µ) .
The existence of the DFE is established. Following the Next Generation Matrix approach (NGM) (see [START_REF] Diekmann | The construction of next-generation matrices for compartmental epidemic models[END_REF]; van den Driessche & Watmough ( 2002)), the basic reproduction number R 0 is computed in order to analyze local and global stability of the DFE depending on R 0 values. Here we use the additional notations: ẋ for the temporal derivative of the vector x and T for the transpose. The system (3.2) is defined with the first two components corresponding to compartments of infected individuals and the last two components corresponding to susceptible compartments:

ẋ = ( ẊI , ẎI , ẊS , ẎS ) T = (0, 0, 0, 0) T . (3.2)
Following the Next Generation Matrix Method (see e.g. van den Driessche & Watmough ( 2002)) we break up ẋ into F -V differentiating new infections from all other changes in population, it gives:

F =       σ f Y I Nm X S σ m X I N f Y S 0 0       and V =        (δ + µ)X I (δ + µ)Y I -Λ + λ f X S -δX I + µX S -Λ + λ m Y S -δY I + µY S        .
Then, the Jacobian matrices of F and V are evaluated at the disease-free equilibrium (DFE). Using the relation:

N = 2 Λ µ dF(P 0 ) = F 0 0 0 and dV(P 0 ) =   V 0 W µI 2   , with F =   0 σ f σ m 0   and V =   δ + µ 0 0 δ + µ   .
Here, W is a square matrix and I 2 is the unit matrix of size 2 × 2. So, the next generation matrix for the model [START_REF] Diekmann | On the definition and the computation of the basic reproduction ratio R 0 in models for infectious diseases in heterogeneous populations[END_REF][START_REF] Van Den Driessche | Reproduction numbers and sub-threshold endemic equilibria for compartmental models of disease transmission[END_REF] is given by:

F V -1 = 1 (δ + µ) 0 σ f σ m 0 .
The quantity R 0 is equal to the spectral radius of F V -1 , thus:

R 0 = R 0,f R 0,m , where R 0,f = σ f (δ + µ) and R 0,m = σ m (δ + µ) .
The basic reproduction number depends on parameters which describe the dynamic of infection (clearance rate, male and female infection rates) and on the retirement rate. The basic reproduction number can be less than unity if the infection rate is lower than the sum of the clearance rate and of the rate of exit of the sexually active population for each sex. This corresponds to diseases which clear quickly with low infection rates. Please note that R 0 is the geometric mean of two values. In a one-sex model: R 0,f = R 0,m , we find R 0 = σ/(δ + µ) which is a classic expression of R 0 in simple Susceptible-Infected-Susceptible model [START_REF] Hethcote | The Mathematics of Infectious Diseases[END_REF]. The quantity R 0,f is the number of secondary infections generated by one infected woman in a population of susceptible men during her infectious period. Then, each infected man can infect a mean number of R 0,m susceptible women during his infectious period. We do not know how many secondary infections are generated by one infected individual in a population of susceptibles. Theoretical computation based on the HPV model may thus provide the means of estimating the basic reproduction number. In section 5, the value of the reproduction number for HPV epidemic is estimated using estimates of parameters involved in the expression of R 0 . In addition, the basic reproduction number is a threshold value for the global dynamics of the model. We obtain the global stability of either disease-free or endemic steady state in terms of the basic reproduction number. Now, the following Theorems focus on the stability of the Disease Free Equilibrium. Detailed proofs can be found in the Appendix.

Theorem 3.1 if R 0 < 1 then the DFE is locally asymptotically stable.

To prove this result, we prove that all eigenvalues of the jacobian matrix evaluated at the DFE have strictly negative real parts.

Theorem 3.2 The DFE is globally asymptotically stable if and only if R 0 ≤ 1.

This Theorem is proved using a suitable Lyapunov function.

We proved the local and global stability of the DFE when R 0 ≤ 1 in the model without vaccination. Now, the endemic equilibrium of the model without vaccination is considered. Using the expression of R 0 , we can rewrite P 1 as:

X * * S = Λ µ R 2 0 + R 0,f R 2 0 + R 0,m 1 R 0,f X * * I = Λ µ R 2 0 -1 R 2 0 + R 0,m Y * * S = Λ µ 1 + R 0,f R 2 0 + R 0,f Y * * I = Λ µ R 2 0 -1 R 2 0 + R 0,f
. X * * I > 0 and Y * * I > 0, and then P 1 is feasible in D, if and only if R 0 > 1 (P 1 = P 0 when R 0 = 1). Previously, we proved that the model without vaccination reaches a steady-state which is a Disease Free Equilibrium when R 0 ≤ 1. The asymptotic behavior of the model has to be studied when R 0 > 1. The following theorem provides the condition for local stability of the endemic equilibrium:

Theorem 3.3 The endemic equilibrium is locally asymptotically stable if and only if R 0 > 1.

To prove this result, we use the same method that for the DFE. We prove that eigenvalues of the jacobian matrix of the linearized system have strictly negative real parts when R 0 > 1. Also the endemic equilibrium is probably globally asymptotically stable for R 0 > 1, but the Lyapunov function used in Theorem 3.2 does not work. Numerical calculations suggest asymptotic stability (see numerical simulations in section 4). Now, let us consider the model including vaccination compartments.

The model with vaccination

Here we focus on the model including vaccination. As we did in the previous sections, we want to study stability of DFE and endemic equilibrium. In a first step, we compute the vaccinated reproduction number. The DFE of (2.1) is given by:

Q 0 = (1 -ϕ f ) Λ µ , 0, (1 -ϕ m ) Λ µ , 0, ϕ f Λ µ , 0, ϕ m Λ µ , 0 .
The Next Generation approach uses only equations of infected persons. We define ẋ = ( ẊI , VI , ẎI , ẆI ) T .

We break up ẏ into F 1 -V 1 and compute the Jacobian matrices of F 1 and V 1 , linearized around the DFE Q 0 . The matrices F 1 and V 1 are defined by:

F 1 = dF 1 (Q 0 ) =        0 0 σ f (1 -ϕ f ) σ f (1 -ϕ f ) 0 0 (1 -τ )σ f ϕ f (1 -τ )σ f ϕ f σ m (1 -ϕ m ) σ m (1 -ϕ m ) 0 0 (1 -τ )σ m ϕ m (1 -τ )σ m ϕ m 0 0        , V 1 = dV 1 (Q 0 ) =        (δ + µ) 0 0 0 0 (δ + µ) 0 0 0 0 (δ + µ) 0 0 0 0 (δ + µ)        .
The next generation matrix is defined as

F 1 V -1 1 with F 1 V -1 1 =        0 0 σ f (δ+µ) (1 -ϕ f ) σ f (δ+µ) (1 -ϕ f ) 0 0 (1 -τ ) σ f (δ+µ) ϕ f (1 -τ ) σ f (δ+µ) ϕ f σm (δ+µ) (1 -ϕ m ) σm (δ+µ) (1 -ϕ m ) 0 0 (1 -τ ) σm (δ+µ) ϕ m (1 -τ ) σm (δ+µ) ϕ m 0 0       
.

An estimate of the vaccinated reproduction number R v which is the spectral radius of

F 1 V -1 1 is: R v = R f (ϕ f )R m (ϕ m ) = R 0 [(1 -ϕ m ) + (1 -τ )ϕ m ][(1 -ϕ f ) + (1 -τ )ϕ f ], with R f (ϕ f ) = R 0,f [(1 -ϕ f ) + (1 -τ )ϕ f ] and R m (ϕ m ) = R 0,m [(1 -ϕ m ) + (1 -τ )ϕ m ].
This result is close to the expression for the vaccinated reproduction number using a S-I-R model for HPV found in [START_REF] Elbasha | Impact of Prophylactic Vaccination against Human Papillomavirus Infection, Contemporary Mathematics 410, Mathematical Studies on Human Disease Dynamics; Emerging Paradigms and Challenges[END_REF][START_REF] Elbasha | Global stability of Equilibria in a Two-sex HPV Vaccination Model[END_REF]. The vaccinated reproduction number depends on:

• The basic reproduction number R 0 ,

• Male and female vaccine coverage (ϕ m and ϕ f ),

• Efficacy of vaccine τ .

Please note that the terms inside brackets are less than one, so R v < R 0 . The term under the square root shows how much the vaccination reduces R 0 . The vaccinated reproduction number is a threshold quantity taking into account vaccination, its expression may assist in the identification of important parameters on which we can act. Bringing the vaccinated reproduction number under unity will lead to the eradication of the virus if the DFE is globally asymptotically stable.

The following theorems provide properties of the system 2.1.

Theorem 3.4 The DFE is locally asymptotically stable if and only if R v < 1.

When R v < 1 the DFE is probably also globally asympotically stable. Numerical calculations suggest global asymptotic stability (see simulations in section 4). When R v > 1, we prove the existence and uniqueness of endemic equilibrium of the vaccination model:

Theorem 3.5 If R v > 1 the endemic equilibrium exists and is unique. If R v < 1 there isn't any endemic equilibrium.

To prove these results (Theorem 3.4 and Theorem 3.5), we use the same method as [START_REF] Elbasha | Impact of Prophylactic Vaccination against Human Papillomavirus Infection, Contemporary Mathematics 410, Mathematical Studies on Human Disease Dynamics; Emerging Paradigms and Challenges[END_REF], the HPV model developed therein is quite similar to our model, the difference is the use of a SIR structure which is more convenient for oncogenic HPV types. In our model, a SIS structure is developed corresponding to non-oncogenic HPV types. Proofs are provided in the Appendix. Finally, we have proved that if the endemic equilibrium exists, it must be unique. Furthermore, if R v > 1 the endemic equilibrium is probably globally asymptotically stable. Numerical calculations suggest this asymptotic behaviour (see simulations in section 4). We proved that if R v < 1, the vaccination model has a disease free equilibrium which is locally asymptotically stable and the endemic equilibrium does not exist; whereas if R v > 1 the endemic equilibrium exists and is unique.

Simulations

In this section, we studied the global stability of the DFE in the model with vaccination (Q 0 ) and of endemic equilibrium in both models (2.1) and (3.1) (with and without vaccination) Q 1 and P 1 using numerical simulations.

Figure 2: Some trajectories in the phase plane portrait for the model with vaccination and considering several sets of initial male and female prevalences when R v < 1

Stability of DFE and of the endemic equilibrium in the model with vaccination

In the model with vaccination, we did not find any Lyapunov function to show that the DFE Q 0 is globally asymptotically stable, thus we conducted analyses using simulations in order to study the asymptotic behavior of the model (2.1) (with vaccination) when the vaccinated reproduction number is less than one (R v < 1). We used Scilab-5.1.1 software. In section 3.2, we found an expression for Q 0 , considering several initial conditions (according to parameters which satisfied R v < 1). We simulated a set of 10,000 combinations for initial prevalences with the conditions N f = 500, 000 and N m = 500, 000. We considered several combinations of parameters such as the female and male infection rates, the clearance rate and the vaccination rate that satisfied the condition: R v < 1. Thereafter, we compared the size of compartments at t=100 years to the expression of Q 0 . Results suggest the stability of the DFE Q 0 when R v < 1 (Figures 2 and3).

Figure 3: Some trajectories for the model with vaccination and considering several sets of initial male and female prevalences when R v < 1 (to test that Q 0 is globally asymptotically stable).

We used same method to study the asymptotic behavior of the model with vaccination when R v > 1 to study the hypothesis that Q 1 is globally asymptotically stable. The same set of 10,000 combinations for initial prevalences was used. Several sets of parameters (infection rates, clearance rates, vaccine coverage, efficacy of vaccine, retirement rate) were chosen to verify the condition: R v > 1. Results suggest that the endemic equilibrium is globally asymptotically stable when R v > 1 (Figure 4).

Figure 4: Some trajectories in the phase plane portrait for the model with vaccination and considering several sets of parameters values when R v > 1. The system reaches a steady-state which is an endemic equilibrium

Based on simulations, we the vaccination model to reach a steady-state depending on R v values : 5).

Q 0 when R v < 1 and Q 1 when R v > 1 (Figure
Figure 5: Some trajectories in the phase plane portrait for the model with vaccination considering the same initial female and male prevalences and several initial values for parameters (female and male infection rates, female and male vaccine coverage). The model with vaccination reaches a steady-state equilibrium which is: the

DFE Q * 0 when R v = 0.71, the DFE Q * * 0 when R v = 0.44, the endemic equilibrium Q * 1 when R v = 1.21, the endemic equilibrium Q * * 1 when R v = 4.58

Global stability of the endemic equilibria P 1

We considered the model without vaccination. As described in the previous section, we simulated several set of initial conditions according to initial prevalences and parameters which verified the condition: R 0 > 1. Results are shown on Figure 6 and suggest that P 1 is probably also globally asymptotically stable when R 0 > 1. ) is used to fit the data of HPV 6/11 prevalences in male and female in order to assess the impact of HPV vaccination on HPV 6/11 prevalence.

Calibration

Table 1 describes the values used for the parameters and inputs of the model. We set the population size in the model to N=1,000,000 individuals, equally divided into females (N f = 500, 000) and males (N m = 500, 000). We used the following values for the size of the 4 compartments at t=0 year (initial values):

• X S (0) = 470, 000

• X I (0) = 30, 000 corresponding to a female HPV 6/11 prevalence of 6% [START_REF] Ralston | Type-specific prevalence and persistence of human papillomavirus in women in the United States who are referred for typing as a component of cervical cancer screening[END_REF] • Y S (0) = 487, 000

• Y I (0) = 13, 000 corresponding to a male HPV 6/11 prevalence of 2.6% [START_REF] Nielson | Multiple-type human papillomavirus infection in male anogenital sites: prevalence and associated factors[END_REF] Male and female infection rates are unknown and have to be assessed in the calibration step. Parameters which have to be assessed are not stochastic variables, thus they do not have a probability distribution describing the probabilities of different values occurring. The method which is used is as follows: a set of 250,000 pairs of infection rates were tested, per capita annual infection rates were selected in the interval [0,5]. Each value for the per capita annual male infection rate in the interval [0,5] was tested with 500 values of per capita annual female infection rates in the interval [0,5]. This procedure enables us to test all combinations of infection rates in the space [0, 5] 2 . The combinations of male and female infection rates which reproduced endemic prevalence of HPV 6/11 (before introduction of vaccination) in males and females within a precision of 10% are shown on Figure 7. On this figure, we can see that it is not possible to give a confidence interval for the estimation of the infection rates. Nonetheless, for some combinations of male and female infection rates which reach the target defined in calibration, several female infection rates could be linked with one male infection rate.

Male Infection Rates

Remark 5.1 (on the calibration results): Male infection rates which match with HPV 6/11 prevalences were lower than female infection rates; the explanation is that male HPV 6/11 prevalences used in the fitting procedure are lower than female HPV 6/11 prevalences. In order to fit HPV prevalence, we observed that the product of male and female infection rates was stable; its range was 1.86 to 1.89. This product appears in the expression of the basic reproduction number in section 3. Female and male infection rates which were used in the analysis above are respectively 2.63 and 0.71.

Using estimates for infections rates and expression of R 0 , we found an estimation for the basic reproduction number: R 0 = 1.04. This value is close to unity corresponding to low prevalences used in the fitting procedure. In the analytical formula for the basic reproduction number, there are two kinds of parameters: the clearance rate δ and the retirement rate µ are estimated using data from the litterature, whereas female and male infection rates are assessed in the fitting procedure, these estimates depend on epidemiological data which are used to fit the model. This shows the importance in choosing data, they should be as representative as possible of the epidemiological reality of the infection or disease. HPV is a sexually transmitted infection, thus we need to consider male and female prevalences corresponding to the same sexually active population or the same country in order to model heterosexual transmission for HPV. There is very few published data on male HPV 6/11 prevalences, this is why we used data from the USA [START_REF] Ralston | Type-specific prevalence and persistence of human papillomavirus in women in the United States who are referred for typing as a component of cervical cancer screening[END_REF][START_REF] Nielson | Multiple-type human papillomavirus infection in male anogenital sites: prevalence and associated factors[END_REF]. [START_REF] Ralston | Type-specific prevalence and persistence of human papillomavirus in women in the United States who are referred for typing as a component of cervical cancer screening[END_REF] in females and 2.6% in males [START_REF] Nielson | Multiple-type human papillomavirus infection in male anogenital sites: prevalence and associated factors[END_REF].

Impact of vaccination

Vaccine coverage is assumed to be constant in time in each scenario. Efficacy of vaccine is assumed to be of 90% (Future I/II study group, 2010). We considered 2 scenarios of vaccination: in scenario 1, we assume that 30% of women (who enter annually in the sexually active population) received the 3 doses of vaccine and are protected, this scenario corresponds to the vaccine coverage observed in France at the launch of the vaccination campaign in 2007 [START_REF] Fagot | HPV vaccination in France: Uptake, costs and issues for the National Health Insurance[END_REF]. In scenario 2, we considered that 10% of women received 3 doses of vaccine, corresponding to the decrease in vaccine coverage observed in France a few years after HPV vaccine introduction [START_REF] Fagot | HPV vaccination in France: Uptake, costs and issues for the National Health Insurance[END_REF]. Considering the first scenario (30% of women who enter annually in the model are protected by the vaccine), a 50% reduction of female HPV 6/11 prevalence is expected 10 years after vaccine introduction, and the male HPV 6/11 prevalence is reduced by one quarter. Thirty years after vaccine introduction, HPV 6/11 prevalences in males and females are expected to be very low (respectively 0.005% and 0.31%). Considering scenario 2, a lower vaccine coverage leads to a slower decrease in HPV prevalence (Table 2 and Figure 8). According to R 0 estimates and expression of R v , the minimum female vaccine coverage which is necessary to eradicate HPV 6/11 (i.e. to have R v < 1) varies between 8.2% and 11%. Thus, if current vaccine coverage in France is maintained, HPV 6/11 will be eradicated.

Figure 8: HPV 6/11 Prevalence in males and females. At t=0 (year), introduction of vaccination considering scenario 1 (30% of women who enter annually in the model are vaccinated) and scenario 2 (10% of women who enter annually in the model are vaccinated). Diamond and triangle represent respectively female and male prevalences in scenario 1, + and X represent respectively female and male prevalences in scenario 2

Summary and Discussion

The developed model in our paper may help to appreciate impact of prophylactic vaccination against HPV on HPV 6/11 prevalence. In this paper, we developed a deterministic model of heterosexually HPV 6/11 transmission. A mathematical analysis of the model was carried out. We identify the vaccinated reproduction number R v as a threshold quantity for the stability of equilibria: if R v > 1, the endemic equilibrium exists and is globally asymptotically stable whereas if R v < 1, the infection-free equilibrium exists and is globally asymptotically stable and HPV will be eliminated. The yielded vaccinated reproduction number from our model represents a reliable parameter establishing whether an epidemic can spread or die out. Its expression depends on other parameters such as vaccine coverage for each sex. Thus, such parameters represent consistent tools which can be used by Public Health policy-makers to improve policies that aim to control HPV epidemic. As an example, level of vaccine coverage can be targeted in order to yield a vaccinated reproduction number that precludes the spread of epidemic. However, driving out the basic reproduction number below one is not always sufficient to eradicate the disease [START_REF] Zhang | Backward bifurcation and global dynamics of an SIS epidemic model with general incidence rate and treatment[END_REF]. It is necessary to study the asymptotic behavior of the model. Our analyses found that our model was asymptotically stable and supported the consistency of the derived parameters. In a second part, we fitted the model on HPV 6/11 prevalences in order to estimate male and female infection rates. We studied the impact of French vaccine strategies on the prevalence of HPV 6/11. The basic reproduction number for these non-oncogenic HPV types was estimated at 1.04. As this value is close to unity, a low proportion of individuals have to be vaccinated in order to eradicate HPV 6/11 in the population. We estimated that if more than 12% of women were vaccinated against HPV using quadrivalent vaccine, the non-oncogenic HPV 6/11 types would be eradicated. Non-oncogenic HPV types are responsible for 90% of genital warts, and the lesions caused by non-oncogenic HPV-6/11 types appear between 2 and 8 months after HPV infection. Therefore, we can expect a dramatic reduction in the occurrence of genital warts caused by HPV 6/11 types. Yielded results are based on assumptions made in the modeling process. Vaccine coverage was assumed to be constant in the time in each scenario. French data did not specify vaccine coverage for each available vaccine; we therefore used global vaccine coverage in French females as a parameter in our model. As only one vaccine prevent HPV 6/11 infections, we were likely to overestimate vaccine coverage against HPV 6/11. French recommendations on HPV vaccination changed in 2010. Previously, quadrivalent vaccine was recommended; while since 2010 French guidelines no longer preferentially support any vaccine. In addition, infection rates were assessed in the calibration step and depend on HPV 6/11 prevalence in males and females. In the absence of reliable French data for HPV 6/11 prevalence in males and females, we used prevalence from studies completed in the USA [START_REF] Ralston | Type-specific prevalence and persistence of human papillomavirus in women in the United States who are referred for typing as a component of cervical cancer screening[END_REF][START_REF] Nielson | Multiple-type human papillomavirus infection in male anogenital sites: prevalence and associated factors[END_REF] which represent consistent data to model the heterosexual transmission of this infection. In Australia a high vaccine coverage (70%) is reached due to a school-based program for HPV vaccination. They observe a dramatic reduction of the number of genital warts (approximately a reduction of 90%) in young women (under 21 years) and in heterosexual young men a few years after initiation of vaccination [START_REF] Read | The near disappearance of genital warts in young women 4 years after commencing a national human papillomavirus (HPV) vaccination programme[END_REF]. This observed effect of vaccination matches with results of our modeling. In France, lower vaccine coverage was observed at the launch of the vaccination campaign and the vaccine coverage is currently decreasing [START_REF] Fagot | HPV vaccination in France: Uptake, costs and issues for the National Health Insurance[END_REF]. Thus, the reduction in HPV 6/11 and in genital warts prevalence will be expected to slow down. Our modeling suggests that in the horizon of 10 years after introduction of vaccination, HPV 6/11 prevalence in females will be halved and HPV 6/11 prevalence in males will be reduced by one quarter, assuming a sustained vaccine coverage of 30% among females (annual rate). In this case, HPV 6/11 will be eradicated within 50 years after vaccine introduction. However, these reductions will be slower if vaccine coverage is decreasing. In our modeling, we did not assess directly the impact of vaccination on the occurrence of Recurrent Respiratory Papillomatosis (RRP). RRP is caused by HPV 6/11 and is observed usually in young children but also among adults. RRP in children is due to HPV 6/11 infection in the respiratory area due to mother-to-child transmission of HPV 6/11 during delivery ( D 'Souza et al., 2011). Thus, if female prevalence for HPV 6/11 is decreasing due to vaccination, a reduction of RRP can be expected. However, prevalence of RRP is very low compared to genital warts ( D 'Souza et al., 2011). To study specifically the impact of quadrivalent vaccine on RRP, a model including vertical transmission could be developed in future research. In order to study the asymptotic behavior of the model, we did not include age structure nor sexual behavior in the model. Thus a simplest model was developed whose asymptotic behavior could be assessed. Nonetheless, the majority of HPV dynamic models take into account age or sexual behavior, in future research the vaccinated reproduction number could be assessed in more complex models including age structure. To conclude, vaccination against HPV, using quadrivalent vaccine, represents a strong tool to prevent HPV 6/11 infections. A dramatic decrease of genital warts occurrence may be expected in France, especially among young individuals.

Appendix: Technical proofs

Local and global stability of Disease Free Equilibrium section: Theorem (3.1): if R 0 < 1 then the DFE is locally asymptotically stable.

Proof. The Jacobian matrix of the system without vaccination ( ẋ = (0, 0, 0, 0) T (3.2)) is evaluated at the DFE P 0 :

J (P 0 ) =     -(δ + µ) σ f 0 0 σ m -(δ + µ) 0 0 δ -σ f -µ 0 -σ m -δ 0 -µ     .
We define:

A 1 = -(δ + µ) σ f σ m -(δ + µ) .
Tr(A 1 ) < 0 and det(A 1 ) > 0 if R 0 < 1 . Thus, if R 0 < 1, all eigenvalues of the jacobian matrix linearized system around the DFE have strictly negative real parts, hence by the Routh-Hurwitz Criterion (Gantmacher, 1959) the DFE, P 0 , is locally asymptotically stable if R 0 < 1.

If R 0 > 1, one eigenvalue has positive real part and the DFE is locally unstable.

Theorem (3.2):

The DFE is globally asymptotically stable if and only if R 0 ≤ 1.

Proof. Consider the Lyapunov function in D:

V = X I + R 0,f Y I .
The derivative of V along the solution of (3) is given by

V ′ = X ′ I + R 0,f Y ′ I = R 0,f σ m Y S N f -(δ + µ) X I + σ f X S N m -R 0,f (δ + µ) Y I . Using Y S ≤ Y * S (because Y S ≤ Y S + Y I ≤ Y * S + Y * I ), X S ≤ X * S and N f = N m = Λ µ : V ′ ≤ (R 0,f σ m -(δ + µ)) X I + (σ f -R 0,f (δ + µ)) Y I V ′ ≤ (δ + µ)(R 2 0 -1)X I . If R 0 ≤ 1 then V ′ ≤ 0. We denote A = R 0,f σ m Y S N f -(δ + µ) and B = σ f X S
Nm -R 0,f (δ + µ) . We can rewrite V ′ = AX I +BY I . We can prove that A ≤ 0 and B ≤ 0. If R 0 < 1, then A < 0 and B < 0, the equality V ′ = 0 holds only at the DFE (when X I = 0 and Y I = 0). If R 0 = 1, V ′ = 0 if and only if:

X I = 0 and Y I = 0 or A = 0 and Y I = 0 or B = 0 and X I = 0 or A = 0 and B = 0

The four cases considered lead to P 0 , thus, V ′ = 0 only in P 0 . The Lasalle-Liapunov theory [START_REF] Hale | Ordinary Differential Equations[END_REF] implies that all paths in D approach the largest positively invariant subset of the set E where V'=0. Here, we have proved in that the only positively invariant subset is {P 0 } so P 0 is globally asymptotically stable for R 0 ≤ 1.

Theorem (3.3):

The endemic equilibrium is locally asymptotically stable if and only if R 0 > 1.

Proof. The Jacobian matrix of the system (3.1) is evaluated at the endemic equilibrium P 1 :

J (P 1 ) =              -(δ + µ) σ f Nm X * * S σ f Nm Y * * I 0 σm N f Y * * S -(δ + µ) 0 σm N f X * * I δ - σ f Nm X * * S - σ f Nm Y * * I -µ 0 -σm N f Y * * S δ 0 -σm N f X * * I -µ             
The charasteristic polynomial is:

p(x) = (µ + x) 2 (δ + µ + x) 2 + σ m X * * I N f (δ + µ + x) - σ m σ f N f N m X * * S Y * * S + σ f Y * * I N m δ + µ + x + σ m X * * I N f .
Thus -µ is a double eigenvalue of this matrix. The two other eigenvalues are the roots of the following polynomial:

q(x) = x 2 + 2(δ + µ) + σ m N f X * * I + σ f N m Y * * I x + (δ + µ) 2 + (δ + µ) σ m N f X * * I + σ f N m Y * * I + σ f σ m N f N m (X * * I Y * * I -X * * S Y * * S ) ∆ = 2(δ + µ) + σ m N f X * * I + σ f N m Y * * I 2 -4 (δ + µ) 2 + (δ + µ) σ m N f X * * I + σ f N m Y * * I + σ f σ m N f N m (X * * I Y * * I -X * * S Y * * S ) .
The discriminant is positive:

∆ = σ m N f X * * I - σ f N m Y * * I 2 >0 + 4 σ m σ f N f N m X * * S Y * * S >0
.

Therefore the 2 solutions of q, x 1 and x 2 are eigenvalues of J (P 1 ).

x 1 = -[2(δ + µ) + σm N f X * * I + σ f Nm Y * * I ] - √ ∆ 2 , x 2 = -[2(δ + µ) + σm N f X * * I + σ f Nm Y * * I ] + √ ∆ 2 .
x 1 is negative, we have to study the sign of x 2 . We prove that x 2 is negative if and only if R 0 > 1 in Lemma 6.1. Therefore, the endemic equilibrium is locally asymptotically stable if and only if R 0 > 1.

Lemma 6.1 x 2 < 0 if and only if R 0 > 1.

Proof. We denote the charasteristic polynomial p

(x) = x 2 + bx + c with b = 2(δ+µ)+ σ m N f X * * I + σ f N m Y * * I andc = (δ+µ) 2 +(δ+µ) σ m N f X * * I + σ f N m Y * * I + σ f σ m N f N m (X * * I Y * * I -X * * S Y * * S )
So,

x 2 = -b + √ b 2 -4c 2 = b 2 -1 + 1 - 4c b 2 x 2 < 0 ⇔ 1 - 4c b 2 < 1. Let us prove that c > 0 ⇔ R 0 > 1.
Using values of X * * S , Y * * S , X * * I , Y * * I and N f = N m = Λ µ , we have:

c = (δ + µ) 2 + (δ + µ) 2 R 0,m - R 0,m R 0,f R 2 0 + R 0,f R 2 0 + R 0,m + R 0,f R 2 0 + R 0,m R 2 0 + R 0,f -1 + (δ + µ) 2 R 2 0 1 - 1 R 0,f R 2 0 + R 0,f R 2 0 + R 0,m R 2 0 + R 0,m R 2 0 + R 0,f - 1 R 0,f -(δ + µ) 2 R 2 0 1 R 0,f R 2 0 + R 0,f R 2 0 + R 0,m 1 + 1 R 0,f - R 2 0 + R 0,m R 2 0 + R 0,f
Rearranging and after simplifying

c = (δ + µ) 2 R 2 0 R 0,m + R 0,f + R 2 0 + 1 R 2 0 + R 0,m R 2 0 + R 0,f ) >0 R 2 0 -1
c is the product of a strictly positive number and (R 2 0 -1). We conclude that: c > 0 iff R 0 > 1.

In the model with vaccination, we prove the following results: Theorem 3.4: The DFE is locally asymptotically stable if and only if R v < 1.

Proof. We define the system ( ẊI , VI , ẎI , ẆI , ẊS , VS , ẎS , ẆS ) T = (0, 0, 0, 0, 0, 0, 0, 0) T (6.1)

We compute the Jacobian matrix of the system (6.1) at its DFE Q 0 .

J 1 = A 1 0 A 2 D 1
with A 1 = F 1 -V 1 and D 1 = -µI 4 , I 4 being identity matrix of size 4; A 2 being a square matrix of size 4; Therefore, (-µ) is an eigenvalue of J 1 of order four. We consider the matrix -A 1 :

-

A 1 =             (δ + µ) 0 -σ f (1 -ϕ f ) -σ f (1 -ϕ f ) 0 (δ + µ) -(1 -τ )σ f ϕ f -(1 -τ )σ f ϕ f -σ m (1 -ϕ m ) -σ m (1 -ϕ m ) (δ + µ) 0 -(1 -τ )σ m ϕ m -(1 -τ )σ m ϕ m 0 (δ + µ)            
.

All diagonal elements of -A 1 are non-negative and all off-diagonal entries are non-positive. Therefore, -A 1 is a Z-matrix. Principal minors of -A 1 are:

C 1 = δ + µ C 2 = (δ + µ) 2 C 3 = (δ + µ) 0 -σ f (1 -ϕ f ) 0 (δ + µ) -(1 -τ )σ f ϕ f -σ m (1 -ϕ m ) -σ m (1 -ϕ m ) (δ + µ)
C 4 = det(-A 1 ).

After some calculations, we express more simply C 3 and C 4 to determine their signs, therefore:

C 3 = C 1 C 2 1 -R 2 v + R f (ϕ f ) σ m (δ + µ) (1 -τ )ϕ m C 4 = (δ + µ) 4 (1 -R 2 v )
C 3 and C 4 are positive if and only if R v < 1. If R v < 1, all principal minors of the Z-matrix -A 1 are positive, so -A 1 is a M-matrix [START_REF] Hashimoto | A note on a criterion for M-matrix[END_REF]. The real parts of each eigenvalues of (-A 1 ) are strictly positive, thus real parts of each eigenvalues of (A 1 ) are strictly negative. Then, the DFE Q 0 is locally asymptotically stable if R v < 1. Otherwise, if R v > 1: C 4 < 0, thus det(-A 1 ) < 0, the determinant of (-A 1 ) is equal to the product of 4 eigenvalues. Hence, at least one of the eigenvalues of the matrix (-A 1 ) has strictly negative real part, and at least one of the eigenvalues of the matrix (A 1 ) has strictly positive real part. Therefore, Q 0 is locally asymptotically unstable. We proved that the DFE of the model with vaccination is locally asymptotically stable if and only if R v < 1.

Theorem 3.5: If R v > 1 the endemic equilibrium exists and is unique. If R v < 1 there isn't any endemic equilibrium.

Proof. We solve the system (2.1) in terms of λ f and λ m .

We obtain:

λ f = σ f (1 -ϕ m ) λ m λ m + δ + µ + ϕ m (1 -τ )λ m (1 -τ )λ m + δ + µ λ m = σ m (1 -ϕ f ) λ f λ f + δ + µ + ϕ f (1 -τ )λ f (1 -τ )λ f + δ + µ .
We define 2 level curves:

G f (λ f , λ m ) = -λ f + σ f (1 -ϕ m ) λ m λ m + δ + µ + ϕ m (1 -τ )λ m (1 -τ )λ m + δ + µ , G m (λ f , λ m ) = -λ m + σ m (1 -ϕ f ) λ f λ f + δ + µ + ϕ f (1 -τ )λ f (1 -τ )λ f + δ + µ .
The 2 levels curves go through the origin and limλ k > 0 when λ ′ k go to infinity. To prove the existence and uniqueness of endemic equilibrium, we have to prove that these 2 level curves intersect only once in the first positive quadrant (λ f > 0, λ m > 0), λ f is on the x-axis . Following the method described by [START_REF] Elbasha | Impact of Prophylactic Vaccination against Human Papillomavirus Infection, Contemporary Mathematics 410, Mathematical Studies on Human Disease Dynamics; Emerging Paradigms and Challenges[END_REF], we prove that : if R v > 1 the level curve G m is above G f around the origin, and the slopes of level curves are positive. Furthermore, each level curve intersects only once a ray from the origin. So, if R v > 1, the endemic equilibrium exists. In our modeling, we assumed that vaccination immunity is lifelong, this particular case was studied by [START_REF] Elbasha | Impact of Prophylactic Vaccination against Human Papillomavirus Infection, Contemporary Mathematics 410, Mathematical Studies on Human Disease Dynamics; Emerging Paradigms and Challenges[END_REF], when R v > 1, the both curves are monotically increasing in the positive quadrant, thus the two level curves intersect only once in the positive quadrant, therefore the endemic equilibrium is unique when R v > 1. If R v < 1, following the proof described by [START_REF] Elbasha | Impact of Prophylactic Vaccination against Human Papillomavirus Infection, Contemporary Mathematics 410, Mathematical Studies on Human Disease Dynamics; Emerging Paradigms and Challenges[END_REF], there is no positive endemic equilibrium.

Figure 6 :

 6 Figure 6: Some trajectories in the phase plane portrait for the model without vaccination considering several sets of initial male and female prevalences when R 0 > 1

Figure 7 :

 7 Figure 7: Combinations of annual infection rates for males and females estimated in calibration step.

Table 1 :

 1 Description of variables and parameters.

	Symbol	Description	Estimates	References
	Variables
	Non-vaccinated population
	X S (t)	Susceptible women
	Y S (t)	Susceptible men
	X I (t)	Infected women
	Y I (t)	Infected men
	Vaccinated population
	V S (t)	Susceptible women
	W S (t)	Susceptible men
	V I (t)	Infected women
	W I (t)	Infected men
	λ f	Force of infection for women
	λ m	Force of infection for men
	N f	Number of females	500,000
	N m	Number of males	500,000
	Biological parameters
	σ f	Infection rate for women		calibration
	σ m	Infection rate for men		calibration
	δ	Clearance rate	1.25	Trottier et al, 2008
	Vaccines Parameters
	ϕ f	Female vaccination rate	scenario 1: 30% Fagot et al, 2011
			scenario 2: 10%
	τ	Vaccine degree of protection	90%	Future, 2010
	Demographic parameters
	Λ	Number of individuals in each sex 30,000 who enter annually in the model	Λ = µ * N/2
	µ	retirement rate	6%

  Table2describes the expected reduction in HPV 6/11 prevalence in each scenario of vaccination in an horizon of 10, 20, 30 and 50 years after introduction of vaccination.

	Table 2. Impact of vaccination on HPV 6/11 prevalence
	Number of years after	Female prevalence Male prevalence
	initiation of vaccination
	Scenario 1
	10	3.55%	1.98%
	20	1.52%	0.88%
	30	0.005%	0.31%
	50	4.10 -3 %	2.10 -4 %
	Scenario 2
	10	4.74%	2.54%
	20	3.68%	2.00%
	30	2.83%	1.53%
	50	1.75%	0.95%
	Initial prevalences are 6%


N Izard

E Gaufrès

X Le Roux

S Kazaoui

Y Murakami

D Marris-Morini

E Cassan

S Maruyama

L Vivien

Electroabsorption study of index-defined semiconducting carbon nanotubes A direct probe into carbon nanotube excitonic states

Keywords: nanotube, electroabsorption, modulator

Electroabsorption spectroscopy of well identified index-defined semiconducting carbon nanotubes is reported. The measurement of high definition electroabsorption spectrum allows direct indexation with unique nanotube chirality. Results show that at least for a limited range of diameter, electroabsorption is directly proportional to the exciton binding energy of nanotubes. Electroabsorption is a powerful technique which directly probes into carbon nanotubes excitonic states, and may become an useful tool for in-situ study of excitons in future nanotube based photonic devices such as electroabsorption modulators.

Introduction

Single Wall Carbon Nanotubes (SWNT) are nano-sized unidimensional structures which have been considerably studied in the past ten years for their very peculiar electrical and optical properties. Semiconducting single-wall a e-mail: nicolas.izard@u-psud.fr carbon nanotubes (s-SWNT) with their direct band-gap, present a special interest for nanophotonic and electrooptic applications [1,2]. Electron confinement in one dimensional nanotube leads to strong electron-hole localization with binding energy as high as 0.5 eV, and the formation of the so-called excitonic states [3]. Exciton binding energy strongly depends on the nanotube intrinsic parameters (diameter, chirality) [4], and on external parameters such as strain [5,6] or dielectric constant of the surrounding medium [7,8]. Therefore, the environnement can be used to tailor s-SWNT optical properties; which could strongly affect s-SWNT properties, making them very efficient sensors [9]. As excitonic processes are fundamental, numerous studies buckle down to elucidate them [10][11][12]. Usually, excitons could be classified under two categories, dipoleallowed "bright" excitons, and dipole-forbidden "dark" excitons, the former being probed by linear absorption spectra while the latter are probed by two-photon absorption spectroscopy [3]. In contrast, electroabsorption (EA) spectroscopy should allow access to full nanotube excitonic states, and would be directly related to excitonic binding energy [13,14].

Experimentally, EA spectroscopy of carbon nanotubes is quite challenging due to the mixing of metallic and semiconducting nanotubes, and the mutually overlapped absorptions of different chirality tubes having proximate excitonic energies [15][16][17]. However, recent advances in SWNT separation techniques, either by density gradient ultracentrifugation [18,19] or by polymer-assisted extraction [20,21], allows to obtain well defined semiconducting nanotubes with a limited extent of chiral index, without detectable traces of metallic nanotubes [22]. In this paper, we report the first electroabsorption experiments performed on well identified index-defined semiconducting tubes, which give a direct insight into s-SWNT electronic structure and excitonic states.

Materials and Experimental details

The current study uses two different kinds of s-SWNT, prepared using a polymer-assisted selective extraction technique [22]. The first sample was made to be as close as possible to a single index nanotube distribution. A narrow diameter SWNT distribution synthesized by the CoMo-CAT process (from SouthWest NanoTech. Inc.) was used as a starting material. This raw material was dispersed with poly-9,9-di-n-octyl-fluorenyl-2,7-diyl (PFO) polymer in a toluene / acetic acid solution using sonication which, after an ultracentrifugation step, leads to sample "CoMo-CAT/PFO" made of almost exclusively of (7,5) chirality s-SWNT [23].

The second sample was made such as to have the widest possible range of s-SWNT chiralities, but at the same time having a distribution of well-defined single chirality peaks.

For that purpose, a large diameter SWNT distribution of as-prepared HiPCO SWNT powder (from Unidym Inc.) was used as a starting material. Details of this sample preparation procedure are reported elsewhere [22]. Briefly, the preparation procedure involves a homogeneization step (1 h with a water-bath sonicator and 15 min with a tip sonicator) of the PFO-SWNT mixture, followed by an ultracentrifugation step (2 h, 150.000 g). This method leads to sample "HiPCO/PFO", with a majority of (8,6), (8,7) and (7,6) s-SWNT chiralities, as well as few (7,5) and (9,7) chiralities [24,25].

Both samples were drop casted on a 2x2 mm quartz substrate, with 100 µm spaced interdigitated Cr/Au electrodes prealably deposited, to form thin layers. Those s-SWNT-doped PFO thin layers were subsequently annealed at 180 • C for 15 min, resulting in highly fluorescent layer, with homogeneous s-SWNT density [2,26].

Electromodulation spectroscopy experiments were then made [27,28]. The sample was set inside a vacuum chamber

Results and Discussion

Room temperature electroabsorption spectrum (EA) of

CoMoCAT/PFO at 12 kV/cm is reported in Fig. 1a, while the corresponding absorption spectrum is shown in Fig. 1b.

A striking feature in the EA spectrum is the unique and well defined peak at 1.18 eV. Interestingly, there is a corresponding feature in the absorption spectrum, also at 1.18 eV. The CoMoCAT/PFO absorption spectra is mainly constituted by two main peaks at 1.18 eV and 1.9 eV.

Those peaks correspond to the E 11 and E 22 optical exciton transitions below the ∆ 11 and ∆ 22 continuums. Using well-known assignation rules [29], one can assign the observed excitonic transitions to index-defined nanotubes. In the case of CoMoCAT/PFO sample, peaks at 1.18 eV and 1.9 eV could respectively be attributed to the (7,5)'s E 11 and the (7,5)'s E 22 excitons. The other two small features around 1.24 eV and 1.4 eV are attributed to trace of different nanotubes, probably (6,5) and (6,4) s-SWNT.

Therefore, the observed peak in the EA spectrum could be assigned to the (7,5) SWNT's first exciton.

In order to further investigate the origin of this EA peak, a series of EA spectrum were performed on sam- ple CoMoCAT/PFO, while varying the applied field from 6 to 16.5 kV/cm. Below 6 kV/cm, the EA signal was hardly distinguishable from noise, while above 16.5 kV/cm the voltage amplification system saturates. The (7,5) s-SWNT EA peak was fitted by a single lorenzian, and amplitude ∆α, peak position, and its full width half maximum (FWHM) were extracted. Peak amplitude as a function of the square electrical field F 2 is shown on Fig. 2a, while peak position and FWHM as a function of the electrical field F are shown in Fig. 2b andc, respectively.

F 2 (kV 2 • cm -2 ) 1.0 1.5 2
No significant peak shift or peak broadening could be detected while applying the electrical field from 6 to 16.5 kV/cm. This means that the field applied to the system remains low, and that the excitonic dissociation is negligible. In that case, EA theory [13,14] states that ∆α obey the following relation:

∆α = K 1 • d 2 t E 2 b • F 2 + K 2 (1)
where E b is the excitonic binding energy, d t the SWNT diameter, F the applied electric field, and K 1 and K 2 , constants. As in the current case, there is only one kind of nanotube, E b and d t are both constant, and ∆α should be linear with F 2 , which is what it was experimentally obtained (Fig. 2a). This demonstrate that the observed EA peak at 1.18 eV in Fig. 1a is from excitonic origin, as it was also observed by Kishida et al. [17] Relations between nanotube diameter and EA spectra were also investigated. Previously prepared HiPCO/PFO sample, which contains several well identified s-SWNT chiralities with well separated absorption features, is ideal for that study. Room temperature EA spectrum at 12 kV/cm and the corresponding absorption spectrum are reported in Fig. 3.

EA spectrum is mainly composed of five peaks, with five corresponding peaks in the absorption spectrum. As sponse of each s-SWNT independently from its concentration in the samples. Plot of this normalized EA amplitude as a function of the nanotube diameter is reported in Fig. 4. Exciton binding energy E b could be given by the following analytical expression from Capaz et al. [30]: Exciton binding energy for the five s-SWNT in which we had measured normalized EA amplitude were calculated using Eq. 2 and reported in Fig.

0 0.1 0.2 0.3 Absorption α (a.u.) (a) 0.8 0.9 1.0 1.1 1.2 1.3 Energy (eV)
E b = 1 d t A + B d t + Cζ + Dζ 2 (2

Summary

We report the first electroabsorption experiments performed on well identified index-defined semiconducting nanotubes.

The electroabsorption response is sufficiently defined to allow a peak-to-peak correspondence with the absorption spectrum, leading to an unprecedented indexation of electroabsorption spectrum. Due to its excitonic nature, electroabsorption amplitude is quadratic with the applied electric field. Our results show that at least for a limited range of diameter, electroabsorption is proportional to exciton binding energy, in good agreement with theoretical works.

Electroabsorption is an interesting technique to directly probe excitons in carbon nanotubes. From that aspect, electroabsorption proves to be an invaluable tool for the direct study of excitonic-based effects in future carbon nanotubes-based photonics devices, and in particular, to develop optical modulators based on electroabsorption effects.

N. Izard thanks the Japan Society for the Promotion of Science and CNRS for financial support.

( 10 - 5

 105 torr) to prevent air breakdown under high electric fields. The light source was a high power Halogen/Xenon white lamp, with a Nikkon monochromator. Electrical modulation was assured by a function generator at a frequency f of 333 Hz, with a voltage amplification up to 300 V. Liquid nitrogen-cooled Si and Ge detectors were used for the detection of the second harmonic (2f ) component. No signal was measured at f and 3f . We experimentally obtained the variation of transmission -∆T which could be related to the change in absorption ∆α by the relation: ∆α ∼ -∆T T .

Fig. 1 .

 1 Fig. 1. (color online) Electroabsorption variation (a) and optical absorption (b) spectra of CoMoCAT/PFO extracted sample. Absorption peaks at 1.18 eV and 1.9 eV correspond respectively to E11 and E22 excitonic transitions of (7,5) SWNT. The electroabsorption spectrum was performed at a field of 12 kV/cm.

Fig. 2 .

 2 Fig. 2. (color online) Amplitude ∆α (a), position (b), and full width half maximum (FWHM) (c) as a function of the electric field applied of (7,5) SWNT electroabsorption peak in CoMo-CAT/PFO extracted sample. Errors bars are the standard 2 σ deviation from the lorentzian fit, while the dashed line is the linear regression.

  it was previously demonstrated, it is possible to assign optical transitions to individual s-SWNT E 11 excitons, ranging from(7,5) to (9,7) s-SWNT, with diameter ranging from 0.83 to 1.10 nm. As the absorption peaks of the HiPCO/PFO sample are well separated from each other, the resulting EA spectrum is constituted by the individual contribution from those five s-SWNT, making the assignment of EA peaks straightforward. An asymmetry could be observed on all EA peaks in the HiPCO/PFO sample, as well as to a lesser extent in the CoMoCAT/PFO sample. This effect is particularly clear with the(8,6) EA peak at

Fig. 3 .

 3 Fig. 3. (color online) Electroabsorption (a) and optical absorption (b) spectra of HiPCO/PFO extracted sample. Absorption peaks at 0.923, 0.968, 1.039, 1.094 and 1.19 eV correspond respectively to(9,7),(8,7),(8,6),(7,6) and (7,5) SWNT E11

  ,5) SWNT E11 excitonic transitions. The electroabsorption spectrum was performed at a field of 12 kV/cm.

Fig. 4 .

 4 Fig. 4. (color online) Normalized electroabsorption response ∆α/α of HiPCO/PFO extracted nanotubes (square, left) and excitonic binding energy E 11 b (circle, right) as a function of the nanotube diameter dt.

2 C 2 D

 22 ) where d t is the tube diameter, and ζ = (-1) ν cos 3θ/d t represents the chirality dependence. A, B, C and D are numerical constants, and their values are respectively[30]: A = 0.6724 eV nm B = -4.910 • 10 -2 eV nm = 4.577 • 10 -2 eV nm -= -8.325 • 10 -3 eV nm 3