DISCRETEANDCONTINUOUSWebsite:http://AIMsciences.orgDYNAMICALSYSTEMS{SE - PDF document

DISCRETEANDCONTINUOUSWebsite:http://AIMsciences.orgDYNAMICALSYSTEMS{SERIESBVolume4,Number2,May2004pp.479{495ANEPIDEMIOLOGYMODELTHATINCLUDESALEAKYVACCINEWITHAGENERALWANINGFUNCTIONJulienArinoDepartmentofMathematicsandStatisticsUniversityofVictoria,VictoriaB.C.,CanadaV8W3P4K.L.CookeDepartmentofMathematicsPomonaCollege,Claremont,CA91711-6348USAP.vandenDriesscheDepartmentofMathematicsandStatisticsUniversityofVictoria,VictoriaB.C.,CanadaV8W3P4J.Velasco-HernandezProgramadeMatem¶aticasAplicadasyComputaci¶onInstitutoMexicanodelPetr¶oleo,EjeCentralL¶azaroC¶ardenas152SanBartoloAtepehuacan,D.F.07730,Mexico(CommunicatedbyLindaAllen)Abstract.Vaccinationthatgivespartialprotectionforbothnewbornsandsusceptiblesisincludedinatransmissionmodelforadiseasethatconfersnoimmunity.Ageneralformofthevaccinewaningfunctionisassumed,andtheinterplayofthistogetherwiththevaccinee±cacyandvaccinationratesisdiscussed.Theintegro-di®erentialsystemdescribingthemodelisstudiedforaconstantvaccinewaningrate,inwhichcaseitreducestoanODEsystem,andforaconstantwaningperiod,inwhichcaseitreducestoasystemofdelaydif-ferentialequations.Forsomeparametervalues,themodelisshowntoexhibitabackwardbifurcation,leadingtotheexistenceofsubthresholdendemicequi-libria.Numericalexamplesarepresentedthatdemonstratetheconsequenceofthisbifurcationintermsofepidemiccontrol.Themodelcanalternativelybeinterpretedasoneconsistingoftwosocialgroups,witheducationplayingtheroleofvaccination.1.Introduction.AclassicalSISepidemicmodelforadiseasethatconfersnoimmunityhasonlyadiseasefreeequilibrium(DFE)whenparametersrenderthebasicreproductionnumberR01,andhasoneendemicstableequilibriumforparametersmakingR01.ThusthediseasediesoutifR01.Byde¯nition(seee.g.,AndersonandMay[1]),R0is\theexpectednumberofsecondarycasesproduced,inacompletelysusceptiblepopulation,byatypicalinfectiveindividual".Aprecisede¯nitionofR0foranODEmodelcanbegivenasthespectralradiusofthenextgenerationmatrix,see,e.g.,[4,25].Intermsofstability,R0isathresholdparameter,suchthatifR01,thentheDFEislocallyasymptoticallystable,andunstableifR01.Inthelattercase,thediseasecangotoanendemiclevelandcontrolstrategiesareusuallyimplementedtoeradicatethediseaseoratleast 1991MathematicsSubjectClassi¯cation.92D30,34K18,34K20.Keywordsandphrases.epidemicmodel,vaccination,backwardbifurcation,delaydi®erentialequationsystem.479 480ARINO,COOKE,VANDENDRIESSCHEANDVELASCO-HERNANDEZtoloweritsprevalencetoreasonablelevels.Suchstrategiesincludetreatmenttocureorincreasethelifeexpectancyofinfectedindividuals,andvaccinationasaprophylacticmeasuretopreventinfection.Inmodelingsuchstrategies,theaimistodeterminethenecessaryamountoftreatmentorvaccination(usuallymodeledasnumberofindividualstreatedorvaccinatedperunittime)sothatthediseasediesout.However,treatmentandvaccinationarenotcompletelye±cient.Vaccinesmayhavelowe±cacyandbe\leaky"(i.e.,successfullyvaccinatedindividualsmayhaveonlypartialprotectionfrominfection).Datasupportthefactthatavacciusuallywanes,thusprovidingonlytemporaryprotection(i.e.,afteracertaintimevaccinatedindividualsbecomesusceptibleagain).Di®erentwaysinwhichavaccinecanfailarediscussedbyMcLeanandBlower[20]andreferencestherein.Kribs-ZaletaandVelasco-Hern¶andezintroducedin[18]amodelofinfectiousdis-easetransmissionforadiseasethatconfersnoimmunityinwhichthesusceptiblepopulationisvaccinatedataratebutthevaccinegivesonlypartialprotection.ThisSISwithvaccinationmodelisappropriatefordiseasessuchaspertussis,tu-berculosis[18],andhepatitisB[17].TheyshowedthatundercertainparameterconditionsthisquitesimpleODEmodel[18]admitsabackwardbifurcation.Toshowthedependenceofthebasicreproductionnumberwithvaccinationonthevaccinationrate,theydenoteditR(),withR(0)=R0.Forcertainparametervalues,thereexistsacriticalvalueRcsuchthatforRcR()1ahysteresise®ectmayarisewithmultipleendemicequilibria.Inthiscase,therearethreeequi-libriaforthesystem:astabletrivialone(correspondingtotheDFE),anunstableendemicequilibriumandalargerstableendemicequilibrium.Thisresulthasim-portantconsequencesforthevaccinationstrategy,andtheaimofthecampaignmustbetoreduceR()belowRc,notmerelybelowone.ForvaluesofR()suchthatRcR()1,thesuccessorfailureofthestrategydependsonthenumberofindividualswhoareinitiallyinfected.Therangeofparametervaluesforwhichbackwardbifurcationispossibleinthismodelisnotnegligible[18].AbackwardbifurcationisalsofoundwhenthisODEmodelisextendedtoincludearecoveredclass[2].BackwardbifurcationhasalsobeenobservedinsomeotherODEepidemicmod-els,includingHIV/AIDSmodels[5,7,16]andmultigroupmodels[8,9],inwhichthebackwardbifurcationisconnectedtoasymmetriesbetweenthedi®erentgroups.Greenhalghetal[6]usedanSISImodelwithvaccinationforanimalinfectionswithincompleteimmunity(e.g.,bovinerespiratorysyncytialvirus,pseudorabiesvirusinpigs)andfoundthatabackwardbifurcationcanoccur.ThisphenomenonisalsofoundinasimpleSISmodel,butwithanonconstantcontactrate[24].Othermodelswithvaccinationhavebeenconsideredintheliterature.Forex-ample,HethcoteandYorke[15,Section4.5]formulatedanSIScore/noncoreODEmodelforgonorrheaandexaminedthee®ectivenessoftwovaccinationstrategies.Ifavaccineforgonorrheabecomesavailable,theypredictedthatitislikelytogivonlytemporaryimmunity[15,p.45]andtheyassumedthatsuchavaccinewouldbetotallye®ective.Theirmodelwithvaccinationofindividualschosenatrandomfromthepopulationatriskshowsthatsuchastrategywouldbeverye®ectiveincontrollinggonorrhea.Multigroupmodelsthatincludetotallye®ectivevaccinationofnewbornsandsusceptibleshavealsobeenconsidered(see,e.g.,[11,13]).In[8,9]modelsfortwosocialgroups\normal"and\educated"wereformulated.Theeducatedgroupcanberegardedasavaccinatedgrouphavingalowertransmissionrateforthedisease.Individualscanmovebackandforthbetweenthetwogroups, DELAYINAVACCINATIONMODEL481accordingtosuchchangesaseducationalstatus,publichealthpolicies.TheSISmodelsstudiedin[8]and[9,Equation(6)]aresimilartotheODEmodelthatwanalyzeinSection3.IntheSISmodelwithvaccinationformulatedin[18],thevaccineisassumedtowaneexponentially.InSection2,weassumeamoregeneralformofthewaningfunctioninformulatinganSISmodelwithvaccinationofthepopulationatriskandafractionofthenewborns.Ouraimistodemonstratetheinterplaybetweenvaccinee±cacy,vaccinationratesandvaccinewaningonthedynamicsofadiseasethatconfersnoimmunity.Asnotedin[18],thelimitingcaseofnorecoveryleadstoanSImodelthatismoreappropriateforfataldiseaseswithvaccinationoreducation.InSection3,wespecializetothecaseinwhichthevaccinewanesexponentiallyandinSection4wespecializetothecaseinwhichthevaccinewaningtimeisaconstant.Wetakeparametersrelevantforsomehumandiseasesanduseacombinationofanalyticalandnumericaltechniquestoconsiderrangesofthevaccinationratesforwhichbackwardbifurcationcanoccur.Thisshowsthatsubthresholdendemicequilibriaarepossible,whichmaybeimportantwhenitcomestodesigningvaccinationstrategies.2.Formulationofthemodel.OurmodelhasthetransferdiagramshowninFigure1.TherearethreeclassesS,IandV,correspondingrespectivelytosus-ceptible,infectiveandvaccinatedindividuals,withnumbersineachclassgivenbyS(t),I(t),V(t),respectively.Asnotedintheintroduction,V(t)mayalternativelycorrespondtoaneducatedclass,butwecontinuetorefertoitasvaccinated.Indi-vidualsmovefromoneclasstotheotherastheirstatuswithrespecttothediseaseevolves.Newindividualsarebornwithabirthrated�0,andaswedonotaccountforverticaltransmissionorimmigrationofinfectives,thisin°owdoesnotentertheIclass.Allindividuals,whatevertheirstatus,aresubjecttodeath,whichoccurswiththerated.Sinceitisassumedthatthediseasedoesnotcausedeath,thetotalpopulationNSIVisconstant.Weassumethataproportion[0;1)ofnewbornsarevaccinatedatbirth;thus®dNentertheVclass,andtheremainder(1)dNentertheSclass.Susceptibleindividuals(regardlessofwhethertheyhavebeenpreviouslyvaccinated)arefurthervaccinatedattherate. Figure1.The°owdiagramforthegeneralmodel.Diseasetransmissionisassumedtobeofstandardincidencetype(see[19]forarecentdiscussionoftransmissionterms),sothatthenumberofinfectivesproducedbyrandomcontactsbetweenIinfectiveandSsusceptibleindividualsisgivenby 482ARINO,COOKE,VANDENDRIESSCHEANDVELASCO-HERNANDEZ¯IS=N,where¯�0isthetransmissioncoe±cient,representingthenumberofadequatecontactsperindividualpertimeunit.Diseaseimmunityisinducedbyvaccination,butthosesuccessfullyvaccinatedmaybeonlypartiallyprotectedfrominfection,resultingininfectedindividualscomingfromthevaccinatedclass.Thenumber1isthedegreeofvaccinee±cacy,with[0;1].If0¾1,thenthevaccineisleaky.If=0,thenthevaccineistotallye®ective.If=1,thenthevaccineisuseless,theVandSclassesareidenticalandthemodelreducestoaclassicalSISmodel.Henceforthweconsider0¾1.ThenumberofinfectivesproducedbyrandomcontactbetweenIinfectiveandVvaccinatedindividualsisgivenby¾¯IV=N.Theparametercanbeinterpretedasthefactorbywhichvaccinationreducesdiseasetransmission.Forsimplicity,itisassumedthatifanindividualintheVclassisinfected,thenthatindividualisequallyasinfectiousasanindividualinfectedfromtheSclass.Thevaccinatedindividualscantheneitherdie(atthenaturaldeathrated),orbecomeinfectiveifthevaccineisleaky(i.e.,¾�0),orhavethevaccineprotectionwearo®,thatis,reentertheSclass.Thislastpointisoneofinteresthere.In[18]itissupposedthatthevaccinewanesexponentiallybuthereweassumeamoregeneralwaningfunctionP(t)forthevaccine.WesupposethatthereisafractionP(t)ofthevaccinatedindividualswhoarestillunderprotectionofthevaccinetunitsafterbeingvaccinated.WesupposethatP(t)isanonnegativeandnonincreasingfunctionwithP(0+)=1,andsuchthatR10P(u)duispositiveand¯nite.TwospecialformsofP(t),namelyanegativeexponentialandastepfunctionareconsideredinSections3and4,respectively.Finally,wesupposethatanyinfectiveindividualscanbecured:membersoftheIclasscanreturntothesusceptibleclass(withnoimmunity),anddosoatarate0(therecoveryrate).Thee®ectofvaccinationisassumedtodisappearafteraninfection:thereisnorecoverytotheVclass.Sincethetotalpopulationremainsconstant,itismoreconvenienttousepro-portionsineachclass.Hereafter,weuseI(t)andV(t)todenotetheproportionofinfectiveandvaccinatedindividuals,respectively,withS(t)=1I(t)V(t),theproportionofsusceptibles.LettheinitialsusceptibleandinfectiveproportionsbeS(0)0;I(0)0andletV0(t)betheproportionofindividualswhoareinitiallyinthevaccinatedclassandforwhomthevaccineisstille®ectiveattimet.AnexpressionforV0(t)isobtainedfromthevaccinationclass-agederivation,see(4)below.Withtheaboveassumptions,weobtainthefollowingintegro-di®erentialsystem.dI(t) dt(S(t)+¾V(t))I(t)(d)I(t)(1a)V(t)=V0(t)+t0(ÁS(u)+®d)P(tu)ed(tu)e¾¯RuI()dxdu(1b)Intheintegraltermin(1b),®distheproportionofvaccinatednewborns,ÁS(u)istheproportionofvaccinatedsusceptibles,P(tu)isthefractionoftheproportionvaccinatedstillintheVclasstutimeunitsaftergoingin(i.e.,notreturnedtoS),ed(tu)isthefractionoftheproportionvaccinatednotdeadduetonaturalcauses,ande¾¯RuI()dxisthefractionoftheproportionvaccinatednotgoneto DELAYINAVACCINATIONMODEL483theinfectiveclass.Thustheintegralin(1b)sumstheproportionofthosewhowerevaccinatedattimeuandremainintheVclassattimet.Foreasyreference,theparametersarecollectedbelow.d�0:naturaldeathrate.0:recoveryrate.¯�0:diseaseinfectivity.0:vaccinationrateofsusceptibles.[0;1):fractionofnewbornsvaccinated.011:degreeofvaccinee±cacy.Amethodtoobtain(1b)andanexpressionforV0(t)istoformulatethemodelwithvaccinationclass-age.Considertheequationforv(t;¿),thedensitywithrespecttovaccinationclass-ageoftheproportionofindividualsinvaccinationclass-agewhoarestillvaccinatedattimet,@ @t@ @¿v(t;¿)=(¾¯I(t)+d())v(t;¿)(2)whereV(t)=R10v(t;¿)d¿.Here()isthevaccinewaningratecoe±cient,withtheproportionstillintheVclassatvaccinationclass-agebeingP()=eR0(q)dq,whichisassumedtosatisfythepreviousassumptionsonthegeneralwaningfunction.Thein°owatvaccinationclass-agezeroisv(t;0)=ÁS(t)+®d,andv(0;¿)0isspeci¯ed.Integratingalongcharacteristicsyieldsv(t;¿)=v(t¿;0)eR(¾¯I(q)+d)dqR0(q)dqfor0tv(t;¿)=v(0;¿t)eR0(¾¯I(q)+d)dqR(q)dqfort·1DividingtheintegralforV(t)att,substitutinginthesolutions,andchangingintegrationvariablesgivesV(t)=t0(ÁS(u)+®d)P(tu)eRu(¾¯I()+d)dxdueR0(¾¯I()+d)dx10v(0;u)P(tu) P(u)du(3)wheretheabovede¯nitionsofv(t;0)andP()havebeenused.Thisisequivalentto(1b)withV0(t)=eR0(¾¯I()+d)dx10v(0;u)P(tu) P(u)du(4)TheratioP(tu)=P(u)=eR+uu(q)dqiswellde¯nedfortuu0andboundedaboveby1.SinceV(0)is¯nite,theintegralin(3)and(4)converges.ThusV0(t)isnonnegative,nonincreasingandlimt!1V0(t)=0.Thisvaccinationclass-ageapproach(similartothatusedforaninfectionclass-agein[14])givesanexplicitexpressionforV0(t).Tobeginanalysisofthemodel,de¯nethesubsetDofthenonnegativeorthantbyDf(S;I;V);S0;I0;V0;SIV=1gToensurethatthemodeliswellposed,andthusbiologicallymeaningful,weneedtoverifythatsolutionsremaininD. 484ARINO,COOKE,VANDENDRIESSCHEANDVELASCO-HERNANDEZTheorem2.1.ThesetDispositivelyinvariantunderthe°owof(1)withI(0)0;S(0)0.Proof.From(1a),andtheinitialconditionI(0)0,I(t)=I(0)eR0((S(u)+¾V(u))(d+)du0forallt0.NotethatI(t)canapproachzeroastapproachesin¯nity.Di®erentiating(1b)givesd dtV(t)=d dtV0(t)+ÁS(t)+®d(d¾¯I(t))(V(t)V0(t))+Q(t)(5)wheretosimplifynotation,wedenoteQ(t)=t0(ÁS(u)+®d)dt(P(tu))ed(tu)e¾¯RuI()dxduSupposethatS(t)0untilsome¯nitetimet1,andthatS(t1)=0.Then,using(1a)and(5),d dt(IV)(t1)=d(I(t1)+V(t1))°I(t1)+Q(t1)+®dd dtV0(t1)+(d¾¯I(t1))V0(t1)SinceS(t1)=0,I(t1)+V(t1)=1.FromthehypothesesonP(v),dt(P(tu))0.SoQ(t1)0,sinceS(t)0for0tt1.Therefored dtS(t1)=d dt(IV)(t1)d(1)°I(t1)+Q(t1)+d dtV0(t1)+(d¾¯I(t1))V0(t1)0sinced(1)0and,by(4),d dtV0(t1)·¡(d¾¯I(t1))V0(t1).ThusdS1=dt�0,andsoS(t)0forallt0.ThisinturnimpliesthatQ(t)0forallt0.From(1b),sinceSandIaswellasP(v)arenonnegative,itfollowsdirectlythatV(t)0forallt.Sincethethreevariablesarenonnegative,itthenfollowsfromI(t)+S(t)+V(t)=1thateachvariablestayslessthanorequalto1foralltandthusthesolutionsremaininD. WiththeassumedinitialconditionsinD,itcanbeshownthatthesystemde¯nedby(1a)and(1b)isequivalenttothesystemde¯nedby(1a)and(5).Theequivalenceof(5)with(1b)canbeseenbyusinganintegratingfactortowrite(5)asd dth(V(t)V0(t))eR0(d+¾¯I())dxi=(ÁS(t)+®dQ(t))eR0(d+¾¯I())dxandnotingthattherighthandsideisd dtt0(ÁS(u)+®d)P(tu)eRu0(d+¾¯I())dxduThesystemde¯nedby(1a)and(5)isofstandardform,thereforeresultsofHaleandVerduynLunel[10,p.43]ensurethelocalexistence,uniquenessandcontinuationofsolutionsofmodel(1). DELAYINAVACCINATIONMODEL485De¯nethebasicreproductionnumberwithvaccinationasRvacR0"1(1)(®d)~P 1+~P#(6)inwhichR0 disthebasicreproductionnumberwithnovaccinationand~P=limt!1t0P(v)edvdvistheaveragelengthoftimethatanindividualremainsvaccinatedbeforelosingprotectionordying.Notethat~P1=d.ThenumberRvac,whichdependson~P,istheimportantquantityinthemodelwithvaccination.WhenonlyoneparametervariesinRvac,wesometimesmakethisdependenceexplicit,e.g.,Rvac()indicatesthatisthebifurcationparameterthatvaries.NotethatRvac·R0,andinthecaseofnovaccination,thatis=0,RvacR0.FromthevaluesofSandVattheDFE(givenintheproofofthefollowingtheorem),Rvacisequaltotheproductofthemeaninfectiveperiod1=(d)andthesumofthecontactrateconstantineachofthesusceptibleandvaccinatedclassesmultipliedrespectivelybytheproportioninthatclassattheDFE,namely¯SDFE¾¯VDFE.Theorem2.2.Formodel(1)withageneralwaningfunction,thereisalwaystheDFE.IfR01,thenthisistheonlyequilibrium,thediseasediesout.IfRvac1,theDFEislocallyasymptoticallystable(l.a.s.),ifRvac1itisunstable.Proof.Equation(1a)hasI=0asanequilibriumandusingI=0inequation(1b)givesV(t)=V0(t)+((1V(t))+®d)t0P(tu)ed(tu)duInthelimitV0(t)=0,thereforeV=((1V)+®d)~Past!1,andthediseasefreeequilibriumpointIDFE=0;VDFE(®d)~P 1+~P;SDFE=1VDF(1®d~P) 1+~Palwaysexists.SinceifR01theonlyequilibriumof(1a)isIDFE=0,itfollowsfromabovethattheDFEistheonlyequilibriumofsystem(1)whenR01.SupposenowthatR01,i.e.,that¯d.Thenequation(1a)givesdI dt(d)((S¾V)1)ISinceI0andS¾VSV1,thisinequalityimpliesthatdI=dt0,andsoI(t)0=IDFEast!1,forallinitialconditionsI(0)0.Linearizing(1a)and(1b)abouttheDFEbysettingI(t)=q(t),V(t)=VDFEr(t)givestheequationsdq(t) dt=((SDFE¾VDFE)(d))q(t) 486ARINO,COOKE,VANDENDRIESSCHEANDVELASCO-HERNANDEZr(t)=t0(r(u)+q(u))P(tu)ed(tu)du¾¯t0((1VDFE)+®d)P(tu)ed(tu)tuq(x)dxdu:Lettingq(t)=C1eztandr(t)=C2ezt,theabovesystembecomestriangularandhasanontrivialsolutionifandonlyifz(SDFE¾VDFE)(d)=(d)(Rvac1)(7a)or1=10P(v)e(d+z)vdv(7b)ast!1.TheseequationsgivetheeigenvalueszattheDFE.Letzxiybearootofequation(7b).ThenbytheproofofLemma2in[24],ifx0,theny=0.Butsince0,equation(7b)hasnononnegativerealroot,thusallofitsrootshavenegativerealparts.Henceby(7a)theDFEisl.a.sifRvac1,andunstableifRvac1. 3.CasereducingtoanODEsystem.Ifweassumethatthevaccinewaningrateisaconstantµ�0,i.e.,P(v)=eµv,andV0(t)=V0(0)e(d+)teR0¾¯I()dxfrom(4),then(1a)and(5)givetheODEsystemdI dt(1I(1)V)I(d)I(8a)dV dt(1IV)¾¯IV(d)V®d(8b)whichwithnonewbornvaccination(=0)isthemodelstudiedin[18].FromTheorem2.2,theDFEwithIDFE=0,SDFE+d(1) d++,VDFE+®d d++alwaysexists.AssumethatR01,thenendemicequilibria(positiveIequilibria,denotedbyI)canbeobtainedanalyticallyfromthequadraticequation(I)=AI2BIC=0whereA¾¯B((d))(d¾Á)C=(d)(d)(Rvac1)=¯withRvacR01(1)(®d) dfrom(6).IfIisapositivesolutionof(I)=0,then(8a)impliesthatV(d¯S)=(¾¯).Substitutingthisvalueintotheequationd(IV)=dt=0,showsthatS0.Then(8b)givesV0.ThusallsolutionswithI0lieinDandarebiologicallyfeasible.Notefrom(8a)thatI11=R0.Backwardbifurcationleadingtotwoendemicequilibriaoccursfor¾�0(i.e.,aleakyvaccine)if0(0)=B�0,(0)=C0andB24AC.Onan(Rvac();I)bifurcationdiagram(seeFigure2),thisoccursforRc()Rvac()1,whereRc()isthevalueofRvac()atthesaddlenodebifurcationpointwherethetwovaluesofIcoincide,i.e.,IIcB=(2A).ForRvac()Rc(),thereisno DELAYINAVACCINATIONMODEL487endemicequilibrium(EEP).ForRvac()1,theconstanttermC�0,andthereisauniqueEEP.Notethatifthevaccineistotallye®ective(=0),thereisnoendemicequilibriumforRvac()1,andnobackwardbifurcationoccurs.BystandardplanarODEarguments(seee.g.,[21]),thefollowingbehaviorcanbeshown,forwhichasketchoftheproofisgiven.Theorem3.1.FortheODEsystem(8)withV(0)0,I(0)0,(i)ifRvacRc,thenthediseasediesout,(ii)ifRcRvac1,thentheEEPwithlargerIisl.a.s.,andtheEEPwithsmallerIisunstable,and(iii)ifRvac1,thentheuniqueEEPisgloballyasymptoticallystableinD¡fI0g.Proof.TheJacobianmatrixJofthelinearizedsystematanendemicequilibriumIhastr(J)0anddet(J)=2¾¯2I(IIc).IfRcRvac1,thentwoEEPexist,withdet(J)0forI�Ic,givinglinearstabilityforthelargerIvalue,whereasdet(J)0forIIc,givinginstabilityatthesmallerIvalue.Thisproves(ii).With1=(IV)asamultiplier,theBendixson-Dulaccriterion[21,p.265]rulesoutperiodicsolutionsintheinteriorofD.Theglobalresultsin(i)and(iii)arecompletedbyusingthePoincar¶eBendixsontheorem[21,p.245]. Using[25,Theorem4]thenatureofthebifurcationatRvac=1isdeterminedbysgn(a)wherea¯SDF¾Á2¯VDF(9)atRvac=1.Ifsgn(a)0thenthebifurcationisforward;whereasifsgn(a)0thebifurcationisbackward.Weillustratethebackwardbifurcationwithanumericalbifurcationdiagram(Figure2)usingparametersappropriateforahumandisease,e.g.,pertussis(see,e.g.,[3]),witha3weekaveragediseaseduration(=0:04762)takingthetimeunitasoneday.Averagelifetimeisassumedtobe75years(d=3:6530E05),andtheaveragenumberofadequatecontactsperinfectiveperdayisestimatedat0.4(=0:4).Assumethatmostbabiesarevaccinatedinthe¯rstfewmonthsoflife,andthatthevaccineise®ective,thusthat=0:9(90%ofnewbornsvaccinated)and=0:1(90%protection).Pertussisvaccinebeginstowaneafterabout3years[3,p.378],andtheaveragewaningtimeofthevaccine1=µisassumedtobe5years,giving=5:4794E04.Withtheseparametervalues,thereisbackwardbifurcationforarangeofvaluesgivenby0:02540:1506(i.e.,vaccinationofsusceptiblesonaverageevery1to8weeks).Withtheaboveparametervalues,R0=8:3936andRvac=0:8807for=0:1,whichisintherangeofbackwardbifurcationsincethecriticalvalueRc=0:8669Rvac1,seeFigure2.Inthebackwardbifurcationrange,thevalueofRvac()mustbedecreasedbelowRctoensurethatthediseaseiscontrolled,i.e.,I0;otherwise,ifIisabovetheunstableEEP,thenItendstothestableendemicvalue.NotethatfromFigure2,valuesofthisendemicequilibriummeanthatbetween50%and80%ofthepopulationisinfectiveformostofthebackwardbifurcationrange.Thisbackwardbifurcationpersists(withtheotherparametervalues¯xedasabove)eveninthecaseinwhichthevaccinedoesnotwane(=0,givingRvac()=0:8396andacriticalvalueRc()=0:7599).Varyingoneparameter,butkeep-ingtheotherparameters¯xedasabove,thebackwardbifurcationpersistsforall[0;1),butnotforthecaseinwhichonlynewbornsarevaccinated(=0).If 488ARINO,COOKE,VANDENDRIESSCHEANDVELASCO-HERNANDEZthereisnorecovery(i.e.,=0),itcanbeseenfrom(9)thatthereisnobackwardbifurcationinthisODEmodel.However,foracorrespondingODEmodelwithdiseaseinducedmortality(i.e.,nonconstantpopulation)andnonewbornvaccina-tion,backwardbifurcationispossiblefor=0(seeFigure4[18]).ThiscasemayberelevantforHIV/AIDS,sincethereiscurrentlygreate®ortto¯ndavaccineforthisdisease,andsuchavaccineisunlikelytogivecompleteprotection. 0.8 0.85 0.9 0.95 1 1.05 1.1 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 Rvac(f)I* Stable EEP Rc(f) Figure2.BackwardbifurcationintheODEvaccinationmodel,obtainedbyvaryingin[0:01;0:2]withotherparametersasgiveninthetext.Hereandinotherbifurcationdiagrams,onlytheen-demicequilibriaareshown(thesolidcurveisthelocallystableEEP,thedashedcurveistheunstableEEP).4.Stepfunction(delay)case.Supposethatthevaccinewaningperiodiscon-stantandequalto,thatisthefunctionP(v)takestheformofastepfunctionona¯niteinterval:P(v)=1ifv[0;!]0otherwiseSinceV0(t)=0fort�!,withS=1IVtheintegralequation(1b)becomes,fort�!V(t)=tt!((1I(u)V(u))+®d)ed(tu)e¾¯RuI()dxdu(10)Di®erentiatingthislastexpression(seeequation(5))givesthedelaydi®erentialequation(DDE)modelasthetwodimensionalsystem,fort�!d dtI(t)=(1I(t)(1)V(t))I(t)(d)I(t)(11a)d dtV(t)=(1I(t)V(t))(1I(t)V(t))ed!e¾¯R!I()dx¾¯IVdV®d1ed!e¾¯R!I()dx(11b) DELAYINAVACCINATIONMODEL489Hereafter,weshifttimebysothattheseequationsholdfort�0.ThewellposednessoftheproblemfollowsfromTheorem2.1andfromthefactthatsolutionsof(1)existandareunique.Foraconstantwaningperiod,thebasireproductionnumberwithvaccinationfrom(6)isRvacR01(1)(®d)(1ed!) d(1ed!)(12)WithIDFE=0,fromTheorem2.2SDFEd®d(1ed!) d(1ed!);VDFE(®d)(1ed!) d(1ed!)(13)Fromnullclines,thereexistsone(ormore)(EEP)i®thereexists0I1suchthatV(I)=(I)(14)where(I)=11=R0I 1(15)for¾1,and(I)=((1I)+®d)(1ed!¾¯!I) (1ed!¾¯!I)+d¾¯I(16)Thisseemshard,ifnotimpossible,tosolveexplicitly.Itishoweverpossibletotackletheproblemnumericallyasdescribedinthefollowingsubsections.4.1.Visualisingandlocatingthebifurcation.Fromthenullclineequations,anEEPexistsi®thereexistsanI(0;1]suchthatequation(14)holds.SowestudythezerosofH(I)=(I)(I)=11=R0I 1((1I)+®d)(1ed!¾¯!I) (1ed!¾¯!I)+d¾¯ITostatetheprobleminaformalway,with¯xedd,letAf®;¯;°;!;Á;¾gbethesetofparametersofinterest,anddenoteH(I;A)=(I)(I)(17)toshowthedependenceontheseparameters.Weproceedasfollows.1.Chooseaparameterai2A.2.Fixallotheraj's(ji).3.Chooseai;min,ai;maxand¢aiforai.4.Forallai;kai;minkai(ksuchthatai;kai;max),computeIsuchthatH(I;ai;k)=0.Step4iscarriedoutusingtheMatLabfzerofunction.IntheDDEcase,analytically¯ndingthepointwhere(17)hasauniquezeroin(0;1]isimpossible(compareH(I)withthequadraticfoundfortheODEcaseinSection3).Itishoweverpossibletoobtainabetterestimatethanbymereobser-vationofthenumericallyobtainedbifurcationdiagram,bymakingthefollowingobservations.ItcanbeshownthatH(0)=Rvac1 (1)R0andthat,for¾1H(1)=1 (1)R0®d(1ed!¾¯!) (1ed!¾¯!)+d¾¯0 490ARINO,COOKE,VANDENDRIESSCHEANDVELASCO-HERNANDEZThereforeforRvac1,thereareseveralpossibilities,whichareillustratedinFigure3.IfRvacRc,thenthereisnoEEP.H(0)andH(1)arestrictlynegative,andnumericalsimulationsseemtoindicatethatHhasnorootsin(0;1](i.e.,thatH0onthisinterval).IfRcRvac1,thenthereareendemicequilibria.Here,sinceH(0)andH(1)arestrictlynegative,theonlypossibilityisthustohaveanevennumberofzerosofH.Numericalsimulationsappeartoindicatethatthenumberofendemicequilibriais2.InbetweenthesetwosituationsRvacRcandthereisoneendemicequilibriumI.Usingthesameprocedureasforthevisualisationofthebifurcation,itispossibletocomputeRcby¯ndingthevalueIsuchthatH(I;A)=0andH0(I;A)=0,foragivenparameterai2A.IfRvac1thenH(0)0andsothereisanoddnumberofendemicequilibria.NumericalsimulationsindicatethatthereisauniqueEEP,seeFigure3.4.2.Numericalsimulations.WeusethesameparametervaluesasintheODEcaseofSection3,exceptthattheconstantwaningtime(thedelay)hastobesubstitutedfor.Wetake=1767,sotheaveragelengthoftime~PthatanindividualstaysvaccinatedisthesameasintheODEmodelofSection3.TheseparametersgiveR0=8:3936andRvac()=0:8808,whichisintherangeofthebackwardbifurcationsince(usingtheabovemethod)Rc()=0:8684.ThebifurcationdiagramisverylikethatdepictedinFigure2.NumericalsimulationsoftheDDEmodelindicatethattherearenoadditionalbifurcations;solutionseithergototheDFEortothe(larger)EEP.Asafunctionofor(ratherthan),thebifurcationhastheshapeshownonFigure2,whereasfortheotherparameters,itappears\reversed"ifplottedasafunctionoftheparameter.Indeed,Rvacisanincreasingfunctionofand,whileitisadecreasingfunctionoftheotherparametersfor¾1,butchangesvery 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 -0.25 -0.2 -0.15 -0.1 -0.05 0 0.05 0.1 IH(I*) Rvac=1.32Rvac=0.88Rvac=0.44 Figure3.ValueofthefunctionH(I),inthreecasescorrespond-ingtothreedi®erentvaluesof(frombottomtotop):0.2,0.4and0.6.OtherparametersasinSection4.2.ThecorrespondingvaluesofRvacareindicatedaslegend.Inallthreecases,R01. DELAYINAVACCINATIONMODEL491littleaschanges.Figure4(a)showsthebifurcationfortheseparametervaluesasafunctionof.Thesituationisclearlydi®erentfromthatofFigure2,sinceinFigure4(a)everyvalueofgivesatleastoneendemicequilibrium.Itshouldbenotedthatthisbehaviorisquiteinterestingintermsofepidemiccontrol.LetmbethevalueofdeterminedbysolvingRvac()=1withRvacgivenby(12).Ifallotherparametersare¯xedasabove,andforsmallwaningtime,0!!m=457:032,givingRvac()1,theonlystableequilibriumisalargeendemicone.Thisisofcourseahighlyundesirablestateintermsofepidemiccontrol.Thenincreasing(i.e.,increasingthewaningtime)pastmallowstheDFEtobecomelocallystable,anditisfoundnumericallythatsolutionsstartingwithI(0)belowtheunstableendemicequilibriumtendtotheDFE.Indeed,considerFigure4(b),whichshowsthebehaviorofI(t)asafunctionoftime.Thisisobtainedbyrunningnumericalintegrationsofsystem(11)usingthepackagedde23[23]withI(t)=fort[!;0],varyingfrom0to1bystepsof0.02.ThevalueoftheunstableEEP(0.088)isshownasadashedline,whichseemstoseparateanendemicfromthediseasefreeasymptoticstate.Increasingmoreisine±cientintermsofdiseasecontrol,since(seeFigure4(a))increasingbeyond1000daysdoesnotfurtherraisethevalueoftheunstableendemicequilibriumandthusdoesnotallowlargerinitialvaluesofinfectivestotendtotheDFE. 0 200 400 600 800 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 tI(t) (b) 0 500 1000 1500 2000 2500 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 I*(a)wwm Figure4.(a)ValuesofIasafunctionofbysolvingH(I;A)=0withai.(b)ValueofI(t)versustime,obtainedbynumericalintegrationofsystem(11)with=1767andinitialdataI(t)=,fort[!;0],varyingfrom0to1bystepsof0.02.Otherparametervaluesasinthetext.Thefollowingtableshows,foragivenparameter,the(rounded)rangesinwhichthedi®erenttypesofbehaviorsareobtained,whenallotherparametershavethe 492ARINO,COOKE,VANDENDRIESSCHEANDVELASCO-HERNANDEZvaluesgiveninSection3(using=1767).TheswitchfromthetwoEEPssit-uationtothesingleEEPsituationisobtainedbysolving,forthevariableunderconsideration,theequationRvac=1.Parameter DFEonly 2EEPs 1EEP ¾ (0;0:085) (0:085;0:114) (0:114;1) (0:143;) (0:025;0:143) [0;0:025) (0;0:362) (0:362;0:454) (0:454;) (0:055;) (0:042;0:055) [0;0:042) impossible (457:032;) (0;457:032) impossible [0,1] impossible5.SpecialcasesoftheDDEmodel.5.1.Case=0(vaccinetotallye®ective).Considerthecasewherethevaccineistotallye®ectivegivingcompleteprotection,i.e.,thatvaccinatedindividualscannevermaketheVtoIdirecttransition.Inthiscase,system(11)reducestodI(t) dt(1I(t)V(t))I(t)(d)I(t)(18a)dV(t) dt(1IV)(1I(t)V(t))ed!dV(t)+®d(1ed!)(18b)HereRvacR01(®d)(1ed!) d(1ed!)andwehavethefollowingresult.Theorem5.1.IfRvac1,thensystem(18)admitsaunique,locallyasymptoti-callystableEEP.Proof.Here(15)is(I)=11=R0Iwhile(16)is(I)=((1I)+®d)(1ed!) (1ed!)+dThereforesolving(14)yieldstheexplicitvalueofI,namelyI11 Rvac(1(1ed!))(whichisbiologicallymeaningfulifRvac1),andinturngivesV(1ed!) dR0(1ed!)LinearizationabouttheEEPyieldsthecharacteristicequationdet(12IV)(d)z¯IÁez!ed!Áez!ed!dz=0UsingtheEEP,the(1,1)entryis¯Iz.Sothecharacteristicequationbecomesz2z[d¯I]+d¯IzÁe!(z+d)=0(19)When=0,thisreducesto(zd)(z¯I)=0andtheEEPisl.a.s.whenthedelayiszero. DELAYINAVACCINATIONMODEL493Supposeziy,y�0.Then(19)isy2iy(d¯I)+d¯IiyÁed!(cos!yisin!y)Takingtheabsolutevalueofeachsidegives(y2d¯I)2y2(d¯I)2y22e2d!Nowsettingy2YgivesY2Y(2(1e2d!)+d22I2+2(d¯I))+d22I2=0Sincefor0eachcoe±cientispositive,Ycannotbepositiveandtherecanbenopureimaginarysolution.Also,z=0isnotasolution,assettingz=0in(19)givesd¯I,whichisstrictlypositive.So,bycontinuity,theEEPisl.a.s.forall0whenitexists,namelyforRvac1. CombiningtheaboveresultwithTheorem2.2showsthatforatotallye®ectivevaccine,thebifurcationisforwardandRvacbehavesasa(local)thresholdasinaclassicalmodel;seetheIntroduction.Multigroupmodelswithatotallye®ectivvaccinealsohavenobackwardbifurcation[11,13,15].5.2.Norecoverycase.Anotherinterestingspecialcaseistheoneinwhichthereisnorecovery(i.e.,=0),correspondingtovaccinationforadiseasewithnocure.Forexample,wereavaccineavailableforHIV/AIDS,thismodelcouldroughly¯t.IntheDDEmodel,consideringthebifurcationdiagramofFigure5(inwhichthedirectionofthebifurcationisreversed,sinceaswaspointedoutinSection4.1,Rvacisadecreasingfunctionof),wecanobservethatfortheparametervaluesused,thereonlyexistsonestableEEPfor=0.Numericalsimulationswithotherparametervaluesseemtoindicatethatfor=0thebifurcationisalwaysforward(asobservedintheODEmodelofSection3). 0 0.01 0.02 0.03 0.04 0.05 0.06 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 gI* Figure5.BifurcationdiagramoftheDDEmodel,asafunctionof:valueofIasafunctionof,otherparametersasinSec-tion4.1.6.Concludingremarks.ThequalitativebehaviorofthemodelappearstoberobustwithrespecttothenatureofthegeneralwaningfunctionP(t).Ourre-sults,basedonnumericalevidence,showthatP(t)forthetwocasesexaminedherepreservesthephenomenonofmultistability.Backwardbifurcationcanoccurforarangeofwaningperiods(bothforexponentialwaningandconstantwaningtime), 494ARINO,COOKE,VANDENDRIESSCHEANDVELASCO-HERNANDEZforaleakyvaccine(onethatdoesnotgivetotalprotection),andforarangeofsusceptiblevaccinationrates.Ourmodelandtheparameterschosenarerelevanttothetransmissionofper-tussis[3].However,amoredetailedmodelofpertussisshouldincludeagestructureandmoreclasses,forexample,classescorrespondingtoindividualswithinfectioacquiredimmunity[12].Somediseasesmayneedotherclassestomakethemodelmorerealistic;forexample,achroniccarrierstageneedstobeincorporatedinamodelforfelinecalicivirus[17,22].However,thebackwardbifurcationfoundhereshouldbeexpectedtopersistinsimilarmodelswithmoreclasses.Thepresenceofbackwardbifurcationhasconsequencesforepidemiccontrol,sincethepresenceofahysteresisloopandaseparatrixbetweenastableDFEandastableEEPmeansthattheendemicstatepersistsforalargerrangeofRvacandalsothattheoutcomeisinitialvaluedependent.Toachievediseasecontrol(I0forallinitialvalues),thevalueofRvac(),whichdependsonthevaccinationpolicy,mustbeloweredtolessthanRc().Foragivenvaccine,thevaccinationrateofsusceptibles()caninprinciplebecontrolled.Thissimplemodelindicatesthatitisimportantthatparametervaluesbeaccuratelyestimatedbeforeavaccinationstrategyisestablishedforaparticulardisease.Anotherinterpretation,similartothatof[8,9],canbegivenforourmodel.Con-siderapopulationinwhichthereexisttwosocialgroupsofindividualssusceptibltoagivenpathogen.ThusVrepresentsindividualswithdi®erentsusceptibilitytotheparticularpathogen(duetoeducation,behavioralchanges,environmentalconditions,biologicalcharacteristics,etc.)Thenisarateofchangingbehav-ior,isameasureofthedi®erenceinsusceptibilitytoinfectionbroughtaboutbythischange,andisthelengthoftimeduringwhichthischangeofbehavioroccurs.Ifthechangeofbehaviordecreasestheriskofcontagion,asassumedin[8,9]thenVindividualsarelesslikely(¾1)tocontractthedisease,andthepossibilityofbackwardbifurcationexists.Thisshouldbetakenintoconsiderationwhendesigningeducationandotherpublichealthpolicies.Acknowledgments.J.A.andP.vdD.acknowledgepartialsupportfromNSERCandMITACS.WethankC.ConnellMcCluskeyforhelpfuldiscussionsontheODEcase,andtherefereesforconstructivecomments.REFERENCES[1]R.M.AndersonandR.M.May.InfectiousDiseasesofHumans.OxfordUniversityPress,1991.[2]J.Arino,C.C.McCluskey,andP.vandenDriessche.Globalresultsforanepidemicmodelwithvaccinationthatexhibitsbackwardbifurcation.Toappear,SIAMJ.Appl.Math.,2003.[3]J.Chin,editor.ControlofCommunicableDiseasesManual.17thEdition.AmericanPublicHealthAssociation,2000.[4]O.DiekmannandJ.A.P.Heesterbeek.MathematicalEpidemiologyofInfectiousDiseases:ModelBuilding,AnalysisandInterpretation.Wiley,2000.[5]J.Dusho®,W.Huang,andC.Castillo-Chavez.Backwardsbifurcationsandcatastropheinsimplemodelsoffataldiseases.J.Math.Biol.,36:227{248,1998.[6]D.Greenhalgh,O.Diekmann,andM.C.M.deJong.Subcriticalendemicsteadystatesinmathematicalmodelsforanimalinfectionswithincompleteimmunity.Math.Biosci.,165:1{25,2000.[7]D.Greenhalgh,M.T.Doyle,andF.Lewis.AmathematicaltreatmentofAIDSandcondomuse.IMAJ.Math.Appl.Med.Biol.,18:225{262,2001.[8]K.P.HadelerandC.Castillo-Chavez.Acoregroupmodelfordiseasetransmission.Math.Biosci.,128:41{55,1995. DELAYINAVACCINATIONMODEL495[9]K.P.HadelerandP.vandenDriessche.Backwardbifurcationinepidemiccontrol.Math.Biosci.,146:15{35,1997.[10]J.HaleandS.VerduynLunel.IntroductiontoFunctionalDi®erentialEquations.Springer-Verlag,1993.[11]H.W.Hethcote.Animmunizationmodelforaheterogeneouspopulation.Math.Biosci.,14:338{349,1978.[12]H.W.Hethcote.Oscillationsinanendemicmodelforpertussis.CanadianAppl.Math.Quart.,6:61{88,1998.[13]H.W.HethcoteandH.Thieme.Stabilityoftheendemicequilibriuminepidemicmodelswithsubpopulations.Math.Biosci.,75:205{227,1985.[14]H.W.HethcoteandP.vandenDriessche.TwoSISepidemiologicmodelswithdelays.J.Math.Biol.,40:3{26,2000.[15]H.W.HethcoteandJ.A.Yorke.GonorrheaTransmissionDynamicsandControl,volume56ofLectureNotesinBiomathematics.Springer-Verlag,1984.[16]W.Huang,K.L.Cooke,andC.Castillo-Chavez.Stabilityandbifurcationforamultiple-groupmodelforthedynamicsofHIV/AIDStransmission.SIAMJ.Appl.Math.,52:835{854,1992.[17]C.Kribs-ZaletaandM.Martcheva.Vaccinationstrategiesandbackwardbifurcationinanage-since-infectionstructuredmodel.Math.Biosci.,177&178:317{332,2002.[18]C.Kribs-ZaletaandJ.Velasco-Hern¶andez.Asimplevaccinationmodelwithmultipleendemicstates.Math.Biosci.,164:183{201,2000.[19]H.McCallum,N.Barlow,andJ.Hone.Howshouldpathogentransmissionbemodelled?TrendsEcol.Evol.,16:295{300,2001.[20]A.R.McLeanandS.M.Blower.ModellingHIVvaccination.TrendsinMicrobiology,3:458{463,1995.[21]L.Perko.Di®erentialEquationsandDynamicalSystems.Springer-Verlag,1996.[22]B.Reade,R.G.Bowers,M.Begon,andR.Gaskell.Amodelofdiseaseandvaccinationforinfectionswithacuteandchronicphases.J.Theoret.Biol.,190:355{367,1998.[23]L.F.ShampineandS.Thompson.SolvingDDEsinMATLAB.Appl.Numer.Math.,37:441{458,2001.[24]P.vandenDriesscheandJ.Watmough.AsimpleSISepidemicmodelwithabackwardbifurcation.J.Math.Biol.,40:525{540,2000.[25]P.vandenDriesscheandJ.Watmough.Reproductionnumbersandsub-thresholdendemicequilibriaforcompartmentalmodelsofdiseasetransmission.Math.Biosci.,180:29{48,2002.ReceivedDecember2002;revisedAugust2003.E-mailaddress:arino@math.mcmaster.ca;KLC04747@pomona.eduE-mailaddress:pvdd@math.uvic.ca;velascoj@imp.mx