# Using a dynamical model to study the impact of a toxoid vaccine on the evolution of a bacterium: The example of diphtheria

Florian Lecorvaisier, Dominique Pontier, Benoît Soubeyrand, David Fouchet

► To cite this version:

Florian Lecorvaisier, Dominique Pontier, Benoît Soubeyrand, David Fouchet. Using a dynamical model to study the impact of a toxoid vaccine on the evolution of a bacterium: The example of diphtheria. Ecological Modelling, 2024, 487, pp.110569. 10.1016/j.ecolmodel.2023.110569. hal-04330503

HAL Id: hal-04330503

https://hal.science/hal-04330503v1

Submitted on 8 Dec 2023

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

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

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

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

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

# Using a dynamical model to study the impact of a toxoid vaccine on the evolution of a bacterium: The example of diphtheria

Florian Lecorvaisier<sup>a, *</sup>, Dominique Pontier<sup>a</sup>, Benoît Soubeyrand<sup>b</sup>, David Fouchet<sup>a</sup>

## ARTICLE INFO

**Keywords:**

Virulence evolution

Diphtheria

Imperfect vaccination

Competition

## ABSTRACT

Vaccines are one of humankind's greatest weapons against infectious diseases. However, their usefulness is often questioned and the public tends to distrust vaccines. A mathematical model published in the early 2000s predicts the selection of more virulent strains of pathogens when populations are protected with imperfect vaccines, i.e., vaccines which reduce but do not entirely block pathogen transmission, such as toxoid vaccines. In this study, we built a disease-specific competition model to analyze the evolution of diphtheria's virulence under the pressure of a toxoid vaccine. Our results show that i) high vaccine coverage favors the emergence and increase prevalence of avirulent (or less virulent) strains of Corynebacterium diphtheriae (the etiologic agent of diphtheria) and ii) that competition between strains is crucial in the eradication of toxigenic strains when toxoid vaccines are used. We conclude that the use of toxoid vaccines could lead to disease eradication if the interaction between strains is taken into account. Our results could extend to biologically similar systems such as pertussis.

## 1. Introduction

Vaccination, one of the great achievements of the 20th century, was initially thought of as a means of protection against the disease or its severe symptoms – for vaccinated individuals. Nevertheless, its action goes beyond the individual to extend the population. When a pathogenic agent is introduced within a population, we define the basic reproduction number of the infection, denoted $R_0$, as the number of secondary cases induced by one infectious individual within this population. Theoretically, knowing this epidemic threshold parameter is an important step to estimate the vaccine coverage of an infectious disease to reach for this $R_0$ to fall below one (Fine, 1993), as an epidemic likely fades out when $R_0$ becomes less than unity.

On May 8, 1980, the World Health Organization (WHO) declared the eradication of smallpox worldwide. Eradication became possible thanks to a high vaccine coverage leading to what we call herd immunity. This term encapsulates the capacity of a population to protect itself from a pathogen thanks to a collective immunity level high enough for the pathogen to be unable to spread among the population (Fine, 1993). As variola virus could only infect humans and transmission could no longer be sustained in the vaccinated populations, smallpox could only die out.

This kind of population-level effect of vaccination can only be seen if

a high enough proportion of the population is vaccinated with a vaccine blocking the transmission of the pathogen it targets. It is classically accepted that the critical proportion of vaccine coverage to reach for a disease to go extinct is given with the formula $p_c = 1 - 1/R_0$ (Fine, 1993). In practice, this is not this simple, as populations are highly heterogeneous, and the behavior of individuals toward vaccination may have complex and sometimes hardly predictable effects on the outcome of an epidemic (Lu et al., 2019; Zou et al., 2023). This phenomenon has been exemplified in the recent COVID-19 pandemic, where the influence of anti-vaccination groups, notably on social networks, may have increased the pandemic length and burden (e.g., Pullan and Dey, 2021).

Nevertheless, with both individual- and population-level benefits, vaccination seemed to be the panacea to every infectious disease. Yet, vaccines can also have deleterious effects. Some vaccines are known to have severe side effects, such as for example whole-cell pertussis vaccines (Edwards and Decker, 2012). Mass vaccination also had unexpected consequences, as for example it has been linked to an increase in the mean age of first infection for measles (Li et al., 2017) and pertussis (Lavine et al., 2011; Rohani et al., 2010), two infections that disproportionately affect children. It has been shown, among adolescents and adults, who are now more frequently infected by childhood diseases like measles and pertussis, that age increases the risk of severe symptoms and

<footnote><sup>a</sup> UMR CNRS 5558 Laboratoire de Biométrie et Biologie Evolutive (LBBE), Université Claude Bernard Lyon 1, Université de Lyon, 69622, Villeurbanne, France</footnote>

<footnote><sup>b</sup> Blossom Vaccinology, 13 rue d'Algérie, 69001, Lyon, France</footnote>

<footnote>* Corresponding author.
E-mail address: florian.lecorvaisier@univ-lyon1.fr (F. Lecorvaisier).</footnote>
hospitalization (e.g., De Serres et al., 2000 for pertussis; Giladi et al., 1987; Yasunaga et al., 2010 for measles). Vaccination, by targeting specific strains in a population, also shapes the evolution of pathogens by altering selection pressures on them, thus favoring bacterial or viral strains evading the vaccine-induced immunity (Gandon and Day, 2008; McLeod and Gandon, 2022).

Another possible effect of mass vaccination on the evolution of pathogens is the selection of strains of different virulence. In the early 2000s, Gandon and coauthors, using a mathematical model, have theorized that vaccines which block the onset of symptoms but not the transmission nor the growth of pathogens – so called “imperfect” or “leaky” vaccines – could lead to the selection of more virulent strains because they remove the cost of producing virulence factors in terms of host mortality without removing the benefit or producing them (Gandon et al., 2003, 2001; but see also McLeod and Gandon, 2022).

The pioneer work by Gandon et al. (2001) brought to light the
possible problems linked to vaccination programs intending to vaccinate
a maximum of individuals without questioning their possible long-term
(i.e., evolutionary) consequences and the need for theoretical works in
the establishment of vaccination policies (Bull and Antia, 2022; Walter
and Lion, 2021). Nevertheless, it has been criticized for its over-
simplification of biological mechanisms (Soubeyrand and Plotkin,
2002). Notably, while the model proposed by Gandon et al. (2001) de-
scribes how vaccines inducing antitoxin immunity should select for
pathogen strains with higher toxin production (i.e., higher virulence), it
fails to account for the metabolic cost associated with such production
(Soubeyrand and Plotkin, 2002). Gandon and coauthors responded to
Soubeyrand and Plotkin by incorporating a cost to toxin production in
their model and concluded that when the cost of toxin production is high
a high vaccine efficacy should select for lower toxin production (Gan-
don et al., 2002). While this conclusion matches empirical observations
for diphtheria, with countries now reporting a high proportion of
non-toxigenic strains over toxigenic strains (e.g., Badell-Ocando et al.,
2022; Benamrouche et al., 2016; Czajka et al., 2018; Dangel et al.,
2018), the model is still general, meaning it cannot fully apprehend the
mechanisms specific to each host-pathogen interaction. More specific
models, which would be built based on the specificity of exemplary
diseases, are needed to better characterize the mechanisms through
which antitoxin vaccines may affect the evolution of their targeted
bacteria.

Thus, the model of Gandon et al. (2001) is fatally limited in face of
the great diversity in host-pathogen interactions and mechanisms at
play. For example, it does not take into account phenomena such as
competition between strains or waning immunity, which could have an
impact on the evolutionary outcome of the use of vaccines. Diphtheria
and pertussis are two human diseases for which the prevalence of
non-toxigenic strains has increased in countries with a high toxoid
vaccine coverage (Guiso et al., 2022; Soubeyrand and Plotkin, 2002).
Diphtheria is an acute and serious disease caused by strains of Coryne-
bacterium species, the most common being Corynebacterium diphtheriae,
a gram-positive bacillus. The disease exists in three forms: invasive,
cutaneous and respiratory, the latter being the most common form of
the disease for indigenous cases in temperate climates (Tiwari and Whar-
ton, 2012). Respiratory diphtheria is caused by the invasion of the upper
respiratory tract by the bacterium, followed by the production and
emission of the diphtheria toxin (DT), a toxin encoded by a bacterio-
phage gene called tox (Pappenheimer and Gill, 1973).

The production of DT leads to the death of the host's pharynx cells
and thus to the formation of a characteristic pseudomembrane. This
pseudomembrane on the pharynx is one of the most common symptoms
of the disease (Dobie, 1979). It can provoke a respiratory distress that
may be fatal, but most of the deaths are imputed to a secondary effect of
the toxin: it can be absorbed by other organs and induce severe com-
plications such as cardiopathy and neuropathy (Tiwari and Wharton,
2012).

There are two ways to fight against the disease: therapy and

prophylaxis (Tiwari and Wharton, 2012). Once the disease has been
clinically confirmed in a patient, therapy can consist in the use of an
antitoxin to reduce morbidity and antibiotics to kill the bacteria. Pro-
phylactic vaccines against diphtheria are composed of a deactivated
form of the DT, called a toxoid (Tiwari and Wharton, 2012). While a
recent study showed that currently used vaccines may contain up to
hundreds of C. diphtheriae proteins (Möller et al., 2019), these vaccines
are designed to stimulate the immunity against this toxin particularly
and not against the whole bacterium, to protect against the deleterious
effects of DT (Scheifele and Ochnio, 2009). DT can thus be considered
the virulence factor of C. diphtheriae. Toxoid vaccines can protect against
the symptomatic forms of the infection (by stimulating the host's de-
fenses against the virulence factor) but not the infection itself. Such
vaccines are thus unable to block the transmission (Fine, 1993). This
particularity is what makes diphtheria an interesting model in the study
of the impact of an imperfect (i.e., non-transmission blocking) vaccine
on the evolution of a pathogen, as it has been predicted that such vaccine
should select for strains of higher virulence (Gandon et al., 2001), but
empirical evidence for diphtheria, in the form of a multiplication of
isolated non-toxigenic strains in different countries with a high vaccine
coverage, weakens this theory (Pappenheimer, 1984; Soubeyrand and
Plotkin, 2002).

In the present study, we build a mathematical model based on the
biology of diphtheria to understand how the vaccine pressure as well as
ecological relationships between a virulent (i.e., toxin-producing) and
an avirulent (i.e., toxin-deficient) strain of C. diphtheriae could drive the
evolution of this bacterium. Our model relies on the hypothesis that the
incidence of the avirulent strain should increase when the proportion of
the vaccinated population increases because it does not pay the cost
associated with the production of a toxin, as it has been theorized by
Soubeyrand and Plotkin (2002) and later by Guiso and coauthors
(2022). Thus, the outcome of the competition between avirulent and
virulent strains in highly vaccinated populations is closely dependent on
the cost of toxin production and may lead to the disappearance of
virulent strains.

2. Methods

2.1. *Model description*

We propose to model the dynamical relationship between two hy-
pothetical strains of *C. diphtheriae* in a host population with different
states of immunity. The first strain is a *tox* positive strain, i.e., able to
produce the DT and therefore it induces damages in non-immunized
individuals. The second strain is a *tox* negative strain, i.e., unable to
produce the toxin. For the sake of simplicity, we consider in this work
that a strain that does not produce the toxin is avirulent (i.e., non-
pathogenic). Here, *tox* positive and *tox* negative strains are referred to
as *Cd*₁ and *Cd*₂, respectively.

The model represents the evolution of the disease within a closed
homogeneous and uniformly mixed population. This disease is trans-
mitted by direct contact with an infected individual. More specifically,
the model is an extension of the standard SIR one, with five compart-
ments that classifies individuals as susceptible (S), infected with the Cd₁
strain (I₁), infected with the Cd₂ strain (I₂), vaccinated (V) and recovered
(R). Fig. 1 represents the flows of individuals among the five epidemi-
ological compartments.

Host reproduction occurs at a constant rate *b*. A *p*₁ fraction of offspring, representing the part of the population vaccinated at birth, enters the V compartment and a 1 − *p*₁ fraction enters the S compartment. All hosts die at a rate *d* representing the natural mortality rate of the population.

Susceptible individuals get vaccinated at a rate $\omega$, which is the product of the frequency of the boosters $f$ and the proportion of the population taking the boosters $p_2$. We allow for vaccinal immunity to wane at a rate $\eta$, moving individuals to the susceptible compartment.
<figure><img src="image_4.png" /><figcaption>Fig. 1. Schematic diagram of the model described in equations (2). Definition of the parameters is given in Table 1 and in the main text. Definition of the compartments is given in the main text. Compartment-to-compartment arrows represent epidemiological transitions between states while void-to-compartment and compartment-to-void arrows represent demographic events (births and deaths, respectively).</figcaption></figure>

We consider three transmission rates: $\beta_1$, $\beta_1'$ and $\beta_2$, representing the transmission rate of $Cd_1$ to a susceptible individual, the transmission rate of $Cd_1$ to a vaccinated individual, and the transmission rate of $Cd_2$ to any kind of host, respectively.

Production of DT confers an ecological advantage to $Cd_1$ in susceptible hosts by increasing the availability of iron that is essential to the bacterial development (Oram and Holmes, 2006). But in vaccinated hosts, DT production becomes disadvantageous for $Cd_1$ because of the cost associated with producing a “useless” toxin (i.e., deactivated by the host’s immune system). This advantage/disadvantage has consequences on the pathogen’s capacity to exploit its host, thus to multiply, and therefore to transmit to a susceptible host. For $Cd_2$, transmission is the same whether the host is susceptible or vaccinated, since this strain does not produce the toxin and thus is not affected by the vaccination status of its host. This is true for a pure toxoid vaccine, and while a recent study showed that other immunogenic proteins can be found in diphtheria vaccines it has not yet been demonstrated that such proteins are present in a high enough concentration to induce protective immunity against them (Möller et al., 2019). We thus consider in this work that diphtheria vaccine is a pure toxoid vaccine.

We define $\beta_2$ as the “basal” transmission rate of C. diphtheriae, i.e., the transmission rate without any cost nor benefits associated with the production of the toxin. For toxin production to be advantageous and thus selected, it must confer a fitness gain, at least in a susceptible population. If the population is immunized against the toxin, then the cost of its production is likely to exceed its benefit. We thus obtain the following relationship between $\beta_1$, $\beta_1'$ and $\beta_2$:

$$ \beta_1 > \beta_2 > \beta_1' \qquad (1) $$

which is the fundamental inequation of our model.

For the sake of simplicity, we consider no superinfection, i.e., an individual infected with one strain cannot be infected with the other strain at the same time. This means that the two strains are in competition for access to hosts: their growth (multiplication) is limited by the growth of the other strain. Moreover, we consider no excess mortality induced by the infection. Virulence is therefore not defined as the mortality induced by the pathogen; we simply consider that $Cd_1$ is virulent and $Cd_2$ is avirulent.

Infected individuals recover after a time period $\gamma$ and enter the R compartment with a newly acquired immunity. Recovered individuals can lose their immunity at a rate $\eta$ and return to the susceptible state (S compartment). In absence of data for diphtheria we assume that both naturally and artificially induced immunity (vaccination) wane at the same rate $\eta$. This may represent an oversimplification since it is presumably wrong for pertussis, a close host-vaccine-pathogen system (Wendelboe et al., 2005).

Nevertheless, immunity acquired through vaccination is expected to differ from that acquired naturally. Vaccination provides immunity against DT whereas natural infection provides immunity against the “whole pathogen” (DT and other bacterial antigens). This means that a non-vaccinated individual infected with $Cd_1$ develops an immunity against DT and other antigens whereas an individual infected with $Cd_2$ develops an immunity only against these other antigens if it has not been vaccinated, or against these antigens and DT if it has been vaccinated. For the sake of simplicity, we consider that a recovered individual cannot be reinfected because we assume that the immunity developed with exposure to these other antigens protects against any new C. diphtheriae infection.

With these assumptions, the epidemiological dynamics can be described by the following ordinary differential equations system:

$$ \frac{dS(t)}{dt} = (1 - p_1)b + \eta(V(t) + R(t)) - \beta_1 \frac{S(t)}{N(t)} I_1(t) - \beta_2 \frac{S(t)}{N(t)} I_2(t) - \omega S(t) - dS(t) \qquad (2a) $$

$$ \frac{dI_1(t)}{dt} = \beta_1 \frac{S(t)}{N(t)} I_1(t) + \beta_1 \frac{V(t)}{N(t)} I_1(t) - \gamma I_1(t) - dI_1(t) \qquad (2b) $$

$$ \frac{dI_2(t)}{dt} = \beta_2 \frac{S(t)}{N(t)} I_2(t) + \beta_2 \frac{V(t)}{N(t)} I_2(t) - \gamma I_2(t) - dI_2(t) \qquad (2c) $$

$$ \frac{dV(t)}{dt} = p_1 b + \omega S(t) - \beta_1 \frac{V(t)}{N(t)} I_1(t) - \beta_2 \frac{V(t)}{N(t)} I_2(t) - \eta V(t) - dV(t) \qquad (2d) $$

$$ \frac{dR(t)}{dt} = \gamma(I_1(t) + I_2(t)) - \eta R(t) - dR(t) \qquad (2e) $$

where $S(t)$, $I_1(t)$, $I_2(t)$, $V(t)$ and $R(t)$ represent the number of individuals in each compartment of S, $I_1$, $I_2$, V and R, respectively, at time t.

## 2.2. Invasion analysis

The basic reproduction number ($R_0$) is defined as the number of secondary cases generated by a single infected individual in an entirely susceptible population at equilibrium. In our model, when focusing on $Cd_1$, we have:

$$ R_0 = \frac{\beta_1}{\gamma + d} \qquad (3) $$

If individuals in the population begin to be vaccinated, the population is no longer “entirely susceptible”, and Eq. (3) cannot describe the invasive capacity of a pathogen in such a context. We thus define new
equations for the reproduction number of each strain.

We consider a hypothetical population at the disease-free equilibrium (DFE). We define here a population at the DFE as either a mixed population (i.e., composed of both susceptible and vaccinated individuals) or as a fully susceptible population. In the simplest case where only one strain spreads, this strain can invade the population and start an epidemic if its reproduction number is superior to one. Let $R_{0,x}$ be the $R_0$ (i.e., the reproduction number) of the $Cd_x$ strain ($x \in (1, 2)$). We consider that an invasion starts with the arrival of a unique infected individual in the population (i.e., $I_x(0) = 1$ for an invasion by $Cd_x$).

In a mixed population, $Cd_1$ can invade if $R_{0,1} > 1$ with:

$$R_{0,1} = \frac{\beta_1 s_{DFE}^*}{\gamma + d} + \frac{\beta_1 v_{DFE}^*}{\gamma + d} \quad (4)$$

where $s_{DFE}^*$ is the proportion of susceptible individuals at the DFE and $v_{DFE}^*$ is the proportion of vaccinated individuals at the DFE. The value of $R_{0,1}$ thus depends on the vaccine coverage.

Reciprocally, $Cd_2$ can invade a population at the DFE if $R_{0,2} > 1$ with:

$$R_{0,2} = \frac{\beta_2}{\gamma + d} \quad (5)$$

which means the value of $R_{0,2}$ is independent of the vaccine coverage. This is expected since vaccination only acts on the transmission of pathogenic strains.

We investigated the conditions promoting the long-term coexistence of the two strains based on the criterion of mutual invasibility, which means two strains can coexist if each can invade a host population (i.e., have a positive growth rate; for a review see Grainger et al., 2019) at the equilibrium in which the other strain would be present alone. Technically, after equilibrium has been established with the first strain, the second strain is introduced in the population with only one infected individual at time $t = 0$. We then search for the conditions under which $dI_1/dt > 0$ (respectively $dI_2/dt > 0$) when the population is at the endemic equilibrium with $Cd_2$ (respectively $Cd_1$) alone. We obtain:

$$\frac{s_1^*}{v_1^*} < \frac{\dot{e}^*}{\varepsilon} < \frac{s_2^*}{v_2^*} \quad (6)$$

with $\dot{e}^* = \beta_2 - \dot{\beta}_1$, $\varepsilon = \beta_1 - \beta_2$ and $s_x^*$ (respectively $v_x^*$) the proportion of susceptible (respectively vaccinated) individuals at the equilibrium where the $Cd_x$ strain is endemic ($x \in (1, 2)$). Knowing (1), $\dot{e}^* > 0$ and $\varepsilon > 0$.

## 2.3. Parametrization

Table 1 summarizes the model parameters. The description of the assigned values can be found in Hammarlund et al. (2016) for duration of immunity, in Truelove et al. (2020) for duration of infection, basic reproduction number, and reduction of $Cd_1$ transmission to vaccinated individuals. The host mortality rate is taken from the Institut National de la Statistique et des Etudes Economiques (INSEE) website (https://www.insee.fr/fr/outil-interactif/5367857/tableau/20_DEM/23_DME, last accessed on January 2nd, 2023). When relevant data are lacking, parameters (see below) have been calculated or chosen arbitrarily. Note that literature is contradictory in the definition of an immunizing antibody level, with threshold varying to authors. In this work, we consider an antibody level of 0.01 UI.mL⁻¹ as protective, as proposed by Amanna et al. (2007) and Hammarlund et al. (2016), and choose a value of $\eta$ adapted to this threshold.

The WHO recommends that adults take booster doses of diphtheria toxoid every ten years (World Health Organisation, 2017). This recommendation is followed in most countries with one exception being France where boosters are recommended at ages 25, 45, 65 and every 10 years beyond that age (Ministère des Solidarités et de la Santé, 2022). The proportion of adults taking the booster is not given by the WHO but

<caption>**Table 1**</caption>

Parameters and their values used in the simulations of the model described in equations (2). d: days.

<table><thead><tr><th>Parameter</th><th>Signification</th><th>Value</th><th>Reference</th></tr></thead><tbody><tr><td>b</td><td>Birth rate</td><td>1/3 d<sup>-1</sup></td><td>Arbitrary</td></tr><tr><td>d</td><td>Death rate</td><td>1/(83 × 365) d<sup>-1</sup></td><td>INSEE</td></tr><tr><td>p<sub>1</sub></td><td>Children vaccine coverage</td><td>Varies</td><td>None</td></tr><tr><td>p<sub>2</sub></td><td>Adults vaccine coverage</td><td>Varies</td><td>See Figs. 2 and 3</td></tr><tr><td>f</td><td>Booster frequency</td><td>1/(10 × 365) d<sup>-1</sup></td><td>World Health Organization, 2017</td></tr><tr><td>ω</td><td>Booster rate</td><td>fp<sub>2</sub></td><td>See text</td></tr><tr><td>η</td><td>Immunity loss</td><td>1/(42 × 365) d<sup>-1</sup></td><td>Hammarlund et al., 2016</td></tr><tr><td>γ</td><td>Recovering rate</td><td>1/(18.5) d<sup>-1</sup></td><td>Truelove et al., 2020</td></tr><tr><td>R<sub>0</sub></td><td>Basic reproduction number</td><td>2.6</td><td>Truelove et al., 2020</td></tr><tr><td>r<sub>1</sub></td><td>Transmission reduction of Cd<sub>1</sub></td><td>0.60</td><td>Truelove et al., 2020</td></tr><tr><td>r<sub>2</sub></td><td>Transmission reduction of Cd<sub>2</sub></td><td>0.42</td><td>See Section 2.3</td></tr><tr><td>β<sub>1</sub></td><td>Transmission rate of Cd<sub>1</sub> to S individuals</td><td>R<sub>0</sub>(γ + d) d<sup>-1</sup></td><td>(7)</td></tr><tr><td>β<sub>1</sub></td><td>Transmission rate of Cd<sub>1</sub> to V individuals</td><td>(1 - r<sub>1</sub>)β<sub>1</sub> d<sup>-1</sup></td><td>(8)</td></tr><tr><td>β<sub>2</sub></td><td>Transmission rate of Cd<sub>2</sub></td><td>(1 - r<sub>2</sub>)β<sub>1</sub> d<sup>-1</sup></td><td>(9)</td></tr></tbody></table>

it can be estimated from multiple studies in different countries: France (Baratin et al., 2012), Mexico (Cruz-Hervert et al., 2013), the United States (Williams et al., 2017) and Germany (Poethko-Müller and Schmitz, 2013). For a review of the status of adult vaccination in Europe, see Kanitz et al. (2012).

Usually, the value of the $R_0$ for a disease is calculated (or rather, should be calculated) at the beginning of its spread in a completely susceptible population (Delamater et al., 2019). In the case of diphtheria, $R_0$ is usually calculated during diphtheria outbreaks in a partially immune population, and thus cannot be rigorously called $R_0$ (e.g., Matsuyama et al., 2018). Nevertheless, using statistical analysis tools, it is possible to infer the value of $R_0$ as it would be in a susceptible population (Truelove et al., 2020). Once $R_0$ is obtained, it is possible to calculate the value of the transmission rate of the strains. Since only toxigenic strains can be monitored (i.e., $Cd_1$ strains in our model), any $R_0$ from field studies can only describe the transmission of such a strain. The transmission rate of $Cd_1$ in a susceptible population, $\beta_1$, can therefore be obtained from Eq. (3):

$$\beta_1 = \frac{R_0}{\gamma + d} \quad (7)$$

Then the transmission rate of $Cd_1$ in a vaccinated population, $\dot{\beta}_1$, can be obtained using the following equation:

$$\dot{\beta}_1 = (1 - r_1)\beta_1 \quad (8)$$

where $r_1$ represents the decrease in the transmission of $Cd_1$ due to vaccination (i.e., the vaccine efficacy).

As said earlier, $R_0$ cannot be obtained for a nontoxigenic strain, i.e., a $Cd_2$ strain in our model, as such strain is not monitored. It is therefore not possible to calculate the transmission rate of $Cd_2$ from Eq. (3). Instead, we calculate the value of the transmission rate of $Cd_2$, $\beta_2$, as a reduction of the transmission rate of $Cd_1$ in a susceptible population:

$$\beta_2 = (1 - r_2)\beta_1 \quad (9)$$

where $r_2$ represents the decrease in the transmission of $Cd_2$ due to the lack of toxin production in this strain. We set $r_1 > r_2$ in agreement with
(1). Because to our knowledge no information about the transmission dynamics of non-toxigenic strains is available yet, we assigned an arbitrary value for $r_2$. Since the production of the toxin is essential for the iron uptake in an iron-depleted growth medium such as the respiratory tract, we estimate that non-production of the toxin represents a high handicap for the strain and thus that its reduction in transmission can be fixed at a high value (Murphy, 1996; Parveen et al., 2019).

The numerical values of the parameters and their sources are summarized in Table 1.

### 3. Results

The mathematical model presented in (2) was implemented in a C script. Simulations from the model were run using the R language (R Core Team, 2022) and the *deSolve* package (Soetaert et al., 2010). Simulations were run long enough to reach the near-equilibrium state ($10^6$ days, so roughly 2740 years).

#### 3.1. Prevalence of Cd₁

We first investigate the propagation of the pathogenic strain $Cd_1$ in the host population, since it is an important goal from a public health point of view.

We run simulations to assess the effect of vaccine coverage and of changes in seven key parameters values on the prevalence of $Cd_1$ in the population, with $Cd_1$ prevalence estimated as the number of $Cd_1$-infected individuals per 10,000 individuals at the equilibrium (Fig. 2). We started with a population at the equilibrium, either with $Cd_2$ ($Cd_2$-endemic equilibrium or CEE, dashed lines in Fig. 2) or with no bacterium in the population (disease-free equilibrium or DFE, solid lines in Fig. 2). Note that for low values of vaccine coverage the dashed and solid lines

are superimposed.

Vaccine coverage values ranged from 0 to 1. The other six parameters could be parameters directly taken from the model equations (e.g., $r_2$, Fig. 2E), or relevant novel parameters adapted from model parameters (e.g., the proportion of individuals vaccinated at birth that would later take the adult booster, $p_2/p_1$, Fig. 2F). For each parameter, we consider between three and five possible values. For example, we explore four booster frequencies to analyze how public health decisions could impact disease spread (Fig. 2B). When the value of one parameter is changed, the values of the other parameters stay the same (see Table 1).

For any parameter value, the higher the vaccine coverage, the lower the prevalence of $Cd_1$ at equilibrium (Fig. 2A–F). We also observe that the best values (from a public health perspective) for immunity-related parameters, i.e., immunity duration, booster frequency, and proportion of individuals vaccinated at birth taking the adult boosters, lead to a lower $Cd_1$ burden (Fig. 2A, B, F). At the same time, when $Cd_1$ has a better transmission capacity, i.e., a higher $R_0$, its equilibrium prevalence is higher (Fig. 2D). Obviously, the effect of the transmission reduction of $Cd_2$, $r_2$, on the prevalence of $Cd_1$ can only be seen if $Cd_2$ is spreading in the population. In this case, we observe that the lower $r_2$ is, the lower the vaccine coverage needed for $Cd_1$ prevalence to decrease (Fig. 2E).

We can also see from Fig. 2 that when $Cd_2$ circulates in the population, the prevalence of $Cd_1$ at equilibrium drops stiffly above certain values of vaccine coverage, and this drop can occur for a more or less important vaccine coverage depending on the parameters values.

Following this observation of a sharp drop in the “absolute” prevalence of $Cd_1$, we wanted to investigate how the relative prevalence of $Cd_1$ changes under the same conditions as those described above (Fig. 3). By relative prevalence, we mean the proportion of infected individuals carrying $Cd_1$, or in other words the $I_1/(I_1+I_2)$ value at equilibrium. Since it would be useless to observe the relative prevalence of $Cd_1$ if

<figure><img src="image_5.png" /><figcaption>**Fig. 2.** $Cd_1$ prevalence at the system equilibrium as vaccine coverage is moved, for different values of (A) immunity duration ($\frac{1}{\eta} \times 365$ years), (B) time between boosters ($\frac{1}{F} \times 365$ years), (C) basic reproduction number ($R_0$), (D) transmission reduction of $Cd_1$ ($r_1$), (E) transmission reduction of $Cd_2$ ($r_2$), and (F) proportion of individuals vaccinated at birth who take the vaccine booster ($\frac{p_2}{p_1}$). For plots (A)–(E) we considered that $p_2$ changed synchronously with $p_1$ (i.e., $\frac{p_2}{p_1} = 1$). Solid lines indicate simulations where $Cd_2$ is endemic while dashed lines indicate simulations where $Cd_2$ is absent from the population.</figcaption></figure>
<figure><img src="image_6.png" /><figcaption>**Fig. 3.** Cd₁ relative prevalence at the endemic equilibrium as vaccine coverage is moved, for different values of (A) immunity duration ($\frac{1}{\eta \times 365}$ years), (B) time between boosters ($\frac{1}{f_1 \times 365}$ years), (C) basic reproduction number ($R_0$), (D) transmission reduction of Cd₁ ($r_1$), (E) transmission reduction of Cd₂ ($r_2$), and (F) proportion of individuals vaccinated at birth who take the vaccine booster ($\frac{p_2}{p_1}$). For plots (A)–(E) we considered that $p_2$ changed synchronously with $p_1$ (i.e., $\frac{p_2}{p_1} = 1$).</figcaption></figure>

there was no Cd₂ in the population, simulations were run starting with a population at the CEE.

We find similar effects of parameters to those observed previously. This means that if a parameter value makes the absolute prevalence fall faster (respectively slower) along the vaccine coverage gradient, it also makes the relative prevalence fall faster (respectively slower). What is of particular interest is that for most vaccine coverage values, the infected population is either infected with Cd₁ or Cd₂, and coexistence of both strains appears to occur only for a short range of vaccine coverage values (Fig. 3A–F).

## 3.2. Mutual invasion criterion

For both strains to coexist, they must be able to mutually invade each other. The criterion for such co-invasion to occur (and reciprocally, for exclusion to occur) is given by Eq. (6).

Parameter $\dot{\epsilon}$ represents the difference in fitness between Cd₂ and Cd₁ within a population of vaccinated hosts. Based on inequality (1), $\dot{\epsilon}$ is positive, showing the advantage of Cd₂ among vaccinated hosts. In the same way, parameter $\epsilon$ represents the difference of fitness between Cd₁ and Cd₂ within a susceptible population. Again, based on (1), $\epsilon$ is positive and shows the advantage of Cd₁ among susceptible hosts. We can separate the terms of the inequation to better understand the conditions under which both strains can coexist. We obtain from one side:

$$ \frac{\dot{\epsilon}}{\epsilon} < \frac{s_2^*}{v_2} \Leftrightarrow \epsilon s_2^* > \dot{\epsilon} v_2^* \qquad (10) $$

and from the other side:

$$ \frac{\dot{\epsilon}}{\epsilon} > \frac{s_1^*}{v_1^*} \Leftrightarrow \dot{\epsilon} v_1^* > \epsilon s_1^* \qquad (11) $$

Inequation (10) implies that the strain Cd₁ can invade a population where Cd₂ is endemic if the advantage of Cd₁ among susceptible individuals relatively to the proportion of susceptible individuals in the population is higher than the advantage of Cd₂ among vaccinated individuals relatively to the proportion of vaccinated individuals in the population.

Reciprocally, inequation (11) implies that strain Cd₂ can invade a population where Cd₁ is endemic if the advantage of Cd₂ among vaccinated individuals relatively to the proportion of vaccinated individuals in the population is greater than the advantage of Cd₁ among susceptible individuals relatively to the proportion of susceptible individuals in the population.

Invasion here implies competition: if a strain cannot invade a population where the other one is endemic, it means that the latter strain is competitively excluding the former one.

## 3.3. Vaccine efficacy

Since $r_1$ represents the loss of transmission for Cd₁ when the targeted host is vaccinated, it can be interpreted as the vaccine efficacy against transmission of the bacterium. A value of $r_1 = 0$ means that the vaccine fails to prevent the transmission of the pathogen, while a value of $r_1 = 1$ means the vaccine perfectly blocks the transmission (a perfect vaccine).

We ran simulations to investigate the simultaneous effect of both vaccine coverage and vaccine efficacy on the reproduction number of Cd₁ (Fig. 4). The simulations started with two types of population: a population at DFE (Fig. 4A) and a population at CEE (Fig. 4B). A single individual infected with Cd₁ was then included in the population and simulations were re-run until a near-equilibrium state was reached again. In both cases, there are values of vaccine efficacy and vaccine coverage for which $R_0 < 1$, which is the condition under which the pathogen does not transmit sufficiently and goes extinct (herd
<figure><img src="image_7.png" /><figcaption>Fig. 4. Effect of the vaccine coverage and vaccine efficacy on the reproduction number of Cd₁ (A) in a population without Cd₂ or (B) with Cd₂ at the endemic equilibrium. For (A) and (B) the white dashed line indicates R₀ = 1. (C) Combination of vaccine coverage and vaccine efficacy for which R₀ = 1 in a population without Cd₂ (solid line) or in a population with Cd₂ at the endemic equilibrium (dashed line, r₂ = 0.42). The ribbon surrounding the dashed line indicates possible positions of the dashed line for r₂ values moving between r₂ = 0.34 (lower boundary) and r₂ = 0.50 (upper boundary).</figcaption></figure>

immunity). This threshold is reached for lower values of vaccine coverage and vaccine efficacy when Cd₂ is circulating in the population than when there is no competitive strain (Fig. 4C).

We ran more simulations with a range of values to investigate how the transmission capacity of Cd₂ could affect this threshold. When the transmission rate of Cd₂ was high (i.e., the value of r₂ was low), lower vaccine coverage and efficacy were needed to reach the extinction threshold (Fig. 4C).

## 4. Discussion

In this study, we developed a disease-specific model to investigate the evolutionary direction followed by C. diphtheriae strains in a context of mass vaccination. We considered the cost associated with toxin production as well as the effect of a toxoid vaccine against toxin-producing strains and concluded that the use of such vaccine should select for bacteria strains of lower virulence. Using a two-strains model, we evaluated the role of different parameters on the selection of each strain (one toxin-deficient, one toxin-producing) and the conditions under which coexistence is possible. We found that when vaccine coverage is high both for children (primary vaccination) and adults (booster vaccination) the toxin-deficient strain is selected. A shorter booster interval in adults also contributes to the selection of toxin-deficient strains. These three parameters, primary coverage, booster coverage, and booster intervals, vary widely among countries, and are among the easiest to influence. While vaccination against diphtheria at an early age is mandatory in most countries, vaccination coverage may be low in some regions, during some periods or at certain ages (Abdi et al., 2021; Lee and Choi, 2017; Marchal et al., 2021; Murhekar et al., 2009; Sein et al., 2016), reducing selective pressure against toxigenic strains. Booster vaccination is often recommended but not mandatory, and the frequency of boosters may differ between countries. In France, it is recommended that adults receive boosters at 25, 45, and 65 years of age, and then every 10 years (Ministère des Solidarités et de la Santé, 2022), but the percentage of individuals receiving the booster is low and

decreases sharply with age (Marchal et al., 2021).

Our findings match clinical data for diphtheria for example in Algeria (Benamrouche et al., 2016), France (Badell-Ocando et al., 2022), Germany (Dangel et al., 2018) and Poland (Czajka et al., 2018; Zasada, 2013; Zasada and Rzeczkowska, 2019), where the use of toxoid vaccines has been associated with an increase in the prevalence of tox- strains. Countries with a good monitoring system on the incidence of C. diphtheriae, like France, still report finding tox+ strains circulating in the population at low prevalence (Badell-Ocando et al., 2022). In addition to the phenomenon described above that vaccination can drop, these countries probably did not reach their equilibrium state since vaccination was introduced less than a century ago. In our simulations, equilibrium was long to reach because of the inclusion of demographic parameters and immunity decay (data not shown), two biological phenomena acting on a time scale much longer than the epidemiological scale. This could explain the persistence of tox+ strains in these countries. Nevertheless, we argue that in the end toxoid vaccines should select for toxin-deficient strains, because toxin production represents a potentially high metabolic cost, which is not compensated by any benefit if the host is resistant to the toxin's action.

Other parameters may be more difficult to act upon, as they depend on the biology of both the pathogen and the host. As expected, the higher the reproductive capacity of a strain, the better it persists in the population. Thus, the greater the relative advantage of a toxigenic strain over a non-toxigenic strain, the more difficult it will be for the latter to invade the population. The transmission of a pathogen cannot be directly controlled, but public health measures can help reduce its R₀, as was demonstrated during the recent COVID-19 pandemic (Talic et al., 2021). Combined with contact restriction, reducing the infectious period (e.g., with antibiotics) could also reduce the transmission of pathogenic strains (Truelove et al., 2020). These measures, if individual-based, could reduce the transmission of symptomatic strains (i.e., pathogenic strains) while having little to no effect on the transmission of non-pathogenic strains.

Our results also illustrate how competition between strains for access
to hosts can contribute to the elimination of pathogenic strains, as non-toxigenic strains reduce the pool of the susceptible population that could be used by toxigenic strains to maintain in the population. In our simulations, we found that, over a wide range of model parameters, the toxigenic strain could only be eliminated if the nontoxigenic strain was also present in the population. This is because the vaccine does not block either the infection by the toxigenic bacterium or its transmission. Instead, it reduces its transmission by removing the benefit of producing the toxin, leaving the bacterium to produce a costly toxin without any benefit, which impacts the transmission rate of the bacterium. In a scenario where the cost of producing the toxin is high enough that the transmission rate decreases sufficiently and the bacterium can no longer spread in the population, vaccination could result in the extinction of the toxigenic strain. Data are currently lacking on the extent to which vaccine reduces transmission of diphtheria, but studies on the similar pertussis system show that toxoid vaccines (which also contains non-toxin antigens) are not sufficient to stop transmission of the pathogen in a baboon model (Warfel et al., 2014).

This study focuses on one pathogen, *C. diphtheriae*, the etiologic agent of diphtheria, but our results may explain the occurrence of non-toxigenic strains in other host-pathogen systems such as pertussis. Most developed countries abandoned conventional whole-cell vaccines (WCVs) in favor of acellular vaccines (ACVs) in the last three decades. These new generation vaccines are similar to diphtheria vaccines in that they are not composed of the whole pathogen, but of certain antigens selected for their role in triggering disease. One of these antigens, present in all ACVs, is pertussis toxin, which is thought to be responsible for most of pertussis symptoms (Carbonetti, 2010; Gregg and Merkel, 2019; Scanlon et al., 2019a, 2019b). Interestingly, at least two countries using ACVs have reported the detection of toxin-deficient strains in recent years (Bouchez et al., 2021; Williams et al., 2016). It is possible that these toxin-deficient strains, rather than being rare, are frequent but difficult to detect, due to the lack of symptoms in the host they infect. Thus, the use of ACVs (i.e., vaccines with a toxoid component) could potentially result in a similar mechanism of competition between toxigenic and non-toxigenic strains, a hypothesis that has not been substantiated so far. Note that for pertussis the massive use of vaccines has been linked to the emergence of strains producing more toxin *in vitro*, but this finding has never been confirmed clinically (Mooi et al., 2009).

One general and important conclusion of our work is that there may be some advantages in having avirulent strains (or at least strains of reduced virulence) of a bacterium circulating in a host population. Our results may support a new and ecological way of vaccine use. Instead of targeting the eradication of pathogens, the aim of some immunization programs could be to reduce the disease to nil by favoring the selection and spread out within the human population of non-pathogenic (i.e., commensal) strains.

This new way of considering the role of vaccination can also be illustrated in rabbits, where the impact of RHDV (Rabbit Haemorrhagic Disease Virus), a highly virulent virus that has decimated rabbit populations worldwide, is reduced in populations with a high level of pre-existing circulation of non-pathogenic strains of the virus (Forrester et al., 2007; Strive et al., 2013; White et al., 2002). The same question arises for COVID-19. Indeed, new variants of importance are described every few months (Ciotti et al., 2022) with characteristics (in terms of transmission and pathogenicity) that are not fully characterized until they invade many countries. This leaves our populations vulnerable to undesirable biological events that could make a new variant both more transmissible and more virulent. Efficient circulation of the pre-existing variant, combined with high vaccination coverage of the most vulnerable part of the population, could be an interesting strategy to reduce the impact of such an unpredictable event.

Finally, it is important to keep in mind that our model has some limitations. As with many SIR-type models, we focus on the biology of the pathogen but not of the host. Diphtheria, like pertussis, can be spread by asymptomatic individuals carrying the pathogen without even

knowing they are infected (Tiwari and Wharton, 2012). Our model does not distinguish susceptible individuals based on their risk of developing a more or less severe form of the disease, but it would be interesting to investigate how asymptomatic individuals might shape the evolutionary dynamics of the disease. Symptomatic individuals are more likely to receive therapeutic treatments (antitoxin and antibiotics) than asymptomatic individuals, thus creating another pressure against virulent strains, which are more likely to be toxigenic. We also do not distinguish individuals based on their “willingness to be vaccinated”. In our model, individuals receive booster vaccination randomly, whereas one would expect the choice to receive a booster to be determined by i) the vaccination status during childhood and ii) whether or not they received a booster in previous years. A recent work also showed that despite being theoretically a toxoid vaccine, diphtheria vaccines also contain other proteins of *C. diphtheriae*, some of which showing immunogenic properties (Möller et al., 2019). This aspect could be incorporated into the model. Finally, our model could be completed by taking into account migratory flux, with new arrivals of *tox+* individuals in the population, participating in the dynamical of maintenance of toxigenic strains in countries with high vaccine coverage as it has already been described in Europe (European Centre for Disease Prevention and Control, 2022).

## CRediT authorship contribution statement

**Florian Lecorvaisier:** Conceptualization, Methodology, Software, Formal analysis, Writing – original draft, Visualization. **Dominique Pontier:** Conceptualization, Supervision, Writing – review & editing. **Benoît Soubeyrand:** Conceptualization, Writing – review & editing. **David Fouchet:** Methodology, Conceptualization, Supervision, Writing – review & editing.

## Declaration of Competing Interest

The authors declare the following financial interests/personal relationships which may be considered as potential competing interests: Benoît Soubeyrand reports a relationship with Blossom Vaccinology that includes: board membership. Co-author BS previously worked for Sanofi.

## Data availability

No data was used for the research described in the article.

## Acknowledgments

We thank the Laboratoire d'Excellence (LabEx) ECOFECT 'Eco-Evolutionary Dynamics of Infectious Diseases' (ANR-11-LABX-0048) of Université de Lyon. We thank the Laboratoire de Biométrie et Biologie Evolutive (LBBE). This work was performed using the computing facilities of the CC LBBE/PRABI. We also thank Nicole Guiso for her expertise on the biology of diphtheria and for her advices in the revision of the manuscript. We also thank the three anonymous reviewers for their helpful comments that led to the improvement of the quality of this manuscript.

## References

* Abdi, I., Gidding, H., Leong, R.N., Moore, H.C., Seale, H., Menzies, R., 2021. Vaccine coverage in children born to migrant mothers in Australia: a population-based cohort study. Vaccine 39, 984–993. https://doi.org/10.1016/j.vaccine.2020.12.058.

* Amanna, I.J., Carlson, N.E., Slifka, M.K., 2007. Duration of humoral immunity to common viral and vaccine antigens. N. Engl. J. Med. 357, 1903–1915. https://doi.org/10.1056/NEJMoa066092.

* Badell-Ocando, E., Toubiana, J., Brisse, S., 2022. Rapport Annuel D'activité Du Centre national De Référence des Corynèbactéries du Complexe Diphtheriae – Bilan des Activités 2017-2021. Institut Pasteur, Paris, France.
Baratin, D., Del Signore, C., Thierry, J., Caulin, E., Vanhems, P., 2012. Evaluation of adult dTPaP vaccination coverage in France: experience in Lyon city. 2010–2011. BMC Public Health 12, 940. https://doi.org/10.1186/1471-2458-12-940.

Benamrouche, N., Hasnaoui, S., Badell, E., Guettou, B., Lazri, M., Guiso, N., Rahal, K., 2016. Microbiological and molecular characterization of Corynebacterium diphtheriae isolated in Algeria between 1992 and 2015. Clin. Microbiol. Infect. 22 https://doi.org/10.1016/j.cmi.2016.08.013, 1005e1-1005e7.

Bouchez, V., Guillot, S., Landier, A., Armatys, N., Matczak, S., the French pertussis microbiology study group, Toubiana, J., Brisse, S., 2021. Evolution of *Bordetella pertussis* outbreak over a 23-year period in France, 1996 to 2018. Eurosurveillance 26, 2001213. https://doi.org/10.2807/1560-7917.ES.2021.26.37.2001213.

Bull, J.J., Antia, R., 2022. Which 'imperfect vaccines' encourage the evolution of higher virulence? Evol. Med. Public Health 10, 202–213. https://doi.org/10.1093/emph/eoa015.

Carbonetti, N.H., 2010. Pertussis toxin and adenylate cyclase toxin: key virulence factors of *Bordetella pertussis* and cell biology tools. Future Microbiol. 5, 455–469. https://doi.org/10.2217/fmb.09.0133.

Ciotti, M., Ciccozzi, M., Pieri, M., Bernardini, S., 2022. The COVID-19 pandemic: viral variants and vaccine efficacy. Crit. Rev. Clin. Lab. Sci. 59, 66–75. https://doi.org/10.1080/10408363.2021.1979462.

Cruz-Hervert, L.P., Ferreira-Guerrero, E., Díaz-Ortega, J.L., Trejo-Valdivia, B., Téllez-Rojo, M.M., Mongua-Rodríguez, N., Hernández-Serrato, M.I., Montoya-Rodríguez, A., García-García, L., 2013. Cobertura de vacunación en adultos y adultos mayores en México. Salud Pública México 55, S300-S306. https://doi.org/10.21149/spm.v55s2.5128.

Czajka, U., Wiatrzyk, A., Mosiej, E., Formińska, K., Zasada, A.A., 2018. Changes in MLST profiles and biotypes of Corynebacterium diphtheriae isolates from the diphtheria outbreak period to the period of invasive infections caused by nontoxigenic strains in Poland (1950–2016). BMC Infect. Dis. 18, 121. https://doi.org/10.1186/s12879-018-3020-1.

Dangel, A., Berger, A., Konrad, R., Bischoff, H., Sing, A., 2018. Geographically diverse clusters of nontoxigenic *Corynebacterium diphtheriae* infection, Germany, 2016–2017. Emerg. Infect. Dis. 24, 1239–1245. https://doi.org/10.3201/eid2407.172026.

Delamater, P.L., Street, E.J., Leslie, T.F., Yang, Y.T., Jacobsen, K.H., 2019. Complexity of the basic reproduction number (R). Emerg. Infect. Dis. 25, 1–4. https://doi.org/10.3201/eid2501.171901.

De Serres, G., Shadmanini, R., Duval, B., Boulianne, N., Déry, P., Douville Fradet, M., Rochette, L., Halperin, S.A., 2000. Morbidity of pertussis in adolescents and adults. J. Infect. Dis. 182, 174–179. https://doi.org/10.1086/315648.

Dobie, R.A., 1979. Clinical features of diphtheria in the respiratory tract. JAMA 242, 2197. https://doi.org/10.1001/jama.1979.03300200027017.

Edwards, K.M., Decker, M.D., 2012. Pertussis vaccines. Vaccines. Elsevier, Amsterdam, pp. 447–492.

European Centre for Disease Prevention and Control, 2022. Increase of reported diphtheria cases due to *Corynebacterium diphtheriae* among migrants in Europe.

Fine, P.E.M., 1993. Herd immunity: history, theory, practice. Epidemiol. Rev. 15, 265–302. https://doi.org/10.1093/oxfordjournals.epirev.a036120.

Forrester, N.L., Trout, R.C., Gould, E.A., 2007. Benign circulation of rabbit haemorrhagic disease virus on Lambay Island, Eire. Virology 358, 18–22. https://doi.org/10.1016/j.viroi.2006.09.011.

Gandon, S., Day, T., 2008. Evidences of parasite evolution after vaccination. Vaccine 26, C4–C7. https://doi.org/10.1016/j.vaccine.2008.02.007.

Gandon, S., Mackinnon, M., Nee, S., Read, A., 2003. Imperfect vaccination: some epidemiological and evolutionary consequences. Proc. R. Soc. London 270, 1129–1136. https://doi.org/10.1098/rspb.2003.2370.

Gandon, S., Mackinnon, M.J., Nee, S., Read, A.F., 2002. Antitoxin vaccines and pathogen virulence. Nature 417. https://doi.org/10.1038/417610a, 610–610.

Gandon, S., Mackinnon, M.J., Nee, S., Read, A.F., 2001. Imperfect vaccines and the evolution of pathogen virulence. Nature 414, 751–756. https://doi.org/10.1038/414751a.

Gilardi, M., Schulman, A., Kedem, R., Danon, Y.L., 1987. Measles in adults: a prospective study of 291 consecutive cases. BMJ 295. https://doi.org/10.1136/bmj.295.6609.1314, 1314–1314.

Grainger, T.N., Levine, J.M., Gilbert, B., 2019. The invasion criterion: a common currency for ecological research. Trends Ecol. Evol. (Amst.) 34, 925–935. https://doi.org/10.1016/j.tree.2019.05.007.

Gregg, K.A., Merkel, T.J., 2019. Pertussis toxin: a key component in pertussis vaccines? Toxins (Basel) 11, 557. https://doi.org/10.3390/toxins11100557.

Guiso, N., Soubeyrand, B., Macina, D., 2022. Can vaccines control bacterial virulence and pathogenicity? *Bordetella pertussis*: the advantage of fitness over virulence. Evol. Med. Public Health 10, 363–370. https://doi.org/10.1093/emph/eoa028.

Hammarlund, E., Thomas, A., Poore, E.A., Amanna, I.J., Rynko, A.E., Mori, M., Chen, Z., Slifka, M.K., 2016. Durability of vaccine-induced immunity against tetanus and diphtheria toxins: a cross-sectional analysis. Clin. Infect. Dis. 62, 1111–1118. https://doi.org/10.1093/cid/ciw066.

Kanitz, E.E., Wu, L.A., Giambi, C., Strikas, R.A., Levy-Bruhl, D., Stefanoff, P., Mereckiene, J., Appelgren, E., D'Ancona, F., 2012. Variation in adult vaccination policies across Europe: an overview from VENICE network on vaccine recommendations, funding and coverage. Vaccine 30, 5222–5228. https://doi.org/10.1016/j.vaccine.2012.06.012.

Lavine, J.S., King, A.A., Bjornstad, O.N., 2011. Natural immune protection in pertussis dynamics and the potential for long-term vaccine failure. Proc. Natl. Acad. Sci. 108, 7259–7264. https://doi.org/10.1073/pnas.1014394108.

Lee, H.-J., Choi, J.-H., 2017. Tetanus-diphtheria-acellular pertussis vaccination for adults: an update. Clin. Exp. Vaccine Res. 6, 22. https://doi.org/10.7774/cevr.2017.6.1.22.

Li, S., Ma, C., Hao, L., Su, Q., An, Z., Ma, F., Xie, S., Xu, A., Zhang, Y., Ding, Z., Li, H., Cairns, L., Wang, H., Luo, H., Wang, N., Li, L., Ferrari, M.J., 2017. Demographic transition and the dynamics of measles in six provinces in China: a modeling study. PLoS Med. 14, e1002255 https://doi.org/10.1371/journal.pmed.1002255.

Lu, Y., Geng, Y., Gan, W., Shi, L., 2019. Impacts of conformist on vaccination campaign in complex networks. Physica A: Stat. Mech. Appl. 526, 121124 https://doi.org/10.1016/j.physa.2019.121124.

Marchal, C., Belhassen, M., Guiso, N., Jacoud, F., Van Ganse, E., Le Pannerer, M., Cohen, R., Verdier, R., Uhart, M., 2021. Vaccination coverage rates for diphtheria, tetanus, poliomyelitis and pertussis booster vaccination in Europe between 2013 and 2017: learnings from an analysis of National Health System Real-World Data. Vaccine 39, 505–511. https://doi.org/10.1016/j.vaccine.2020.12.021.

Matsuyama, R., Akhmetzhanov, A.R., Endo, A., Lee, H., Yamaguchi, T., Tsuzuki, S., Nishiura, H., 2018. Uncertainty and sensitivity analysis of the basic reproduction number of diphtheria: a case study of a Rohingya refugee camp in Bangladesh. November–December 2017 PeerJ 6, e4583. https://doi.org/10.7717/peerj.4583.

McLeod, D.V., Gandon, S., 2022. Effects of epistasis and recombination between vaccine-escape and virulence alleles on the dynamics of pathogen adaptation. Nat. Ecol. Evol. 6, 786–793. https://doi.org/10.1038/s41559-022-01709-y.

Ministère des Solidarités et de la Santé, 2022. Calendrier Des Vaccinations Et Recommandations Vinales 2022.

Möller, J., Kraner, M., Sonnewald, U., Sangal, V., Tittlbach, H., Winkler, J., Winkler, T., H. Melnikov, V., Lang, R., Sing, A., Mattos-Guaraldi, A.L., Burkovski, A., 2019. Proteomics of diphtheria toxoid vaccines reveals multiple proteins that are immunogenic and may contribute to protection of humans against Corynebacterium diphtheriae. Vaccine 37, 3061–3070. https://doi.org/10.1016/j.vaccine.2019.04.059.

Mooi, F.R., van Loo, I.H.M., van Gent, M., He, Q., Bart, M.J., Heuvelman, K.J., de贵西, G.F., Diavatopoulos, D., Teunis, P., Nagelkerke, N., Mertsola, J., 2009. *Bordetella pertussis* strains with increased toxin production associated with pertussis resurgence. Emerg. Infect. Dis. 15, 1206–1213. https://doi.org/10.3201/eid1508.081511.

Murhekar, M.V., Bitragunta, S., Hutin, Y., Ckakravarty, A., Sharma, H.J., Gupte, M.D., 2009. Immunization coverage and immunity to diphtheria and tetanus among children in Hyderabad, India. J. Infect. 58, 191–196. https://doi.org/10.1016/j.jinf.2008.12.010.

Murphy, J.R., 1996. Corynebacterium diphtheriae. In: Baron, S. (Ed.), Medical Microbiology. University of Texas Medical Branch at Galveston, Galveston, Texas.

Oram, D.M., Holmes, R.K., 2006. Diphtheria toxin. The Comprehensive Sourcebook of Bacterial Protein Toxins. Academic Press, Paris, pp. 245–256.

Pappenheimer, A.M., 1984. Diptheria. In: Germanier, R. (Ed.), Bacterial Vaccines. Academic Press, Orlando, pp. 1–36.

Pappenheimer, A.M., Gill, D.M., 1973. Diptheria: recent studies have clarified the molecular mechanisms involved in its pathogenesis. Science 182, 353–358. https://doi.org/10.1126/science.182.4110.353.

Parveen, S., Bishai, W.R., Murphy, J.R., 2019. *Corynebacterium diphtheriae*: diphtheria toxin, the tox operon, and its regulation by Fe2⁺ activation of apo-DtxR. Microbiol. Spectr. 7 https://doi.org/10.1128/microbiolspec.GPP3-0063-2019, 7.4.4.

Poetko-Müller, C., Schmitz, R., 2013. Vaccination coverage in German adults. Bundesgesundheitsbl 56, 845–857. https://doi.org/10.1007/s00103-013-1693-6.

Pullan, S., Dey, M., 2021. Vaccine hesitancy and anti-vaccination in the time of COVID-19: a Google trends analysis. Vaccine 39, 1877–1881. https://doi.org/10.1016/j.vaccine.2021.03.019.

R Core Team, 2022. R: A Language and Environment for Statistical Computing. R Foundation for Statistical Computing, Vienna, Austria.

Rohani, P., Zhong, X., King, A.A., 2010. Contact network structure explains the changing epidemiology of pertussis. Science 330, 982–985. https://doi.org/10.1126/science.1194134.

Scanlon, K., Skerry, C., Carbonetti, N., 2019a. Role of major toxin virulence factors in pertussis infection and disease pathogenesis. In: Fedele, G., Ausiello, C.M. (Eds.), Pertussis Infection and Vaccines, Advances in Experimental Medicine and Biology. Springer International Publishing, Cham, pp. 35–51. https://doi.org/10.1007/5584_2019.403.

Scanlon, K., Skerry, C., Carbonetti, N., 2019b. Association of pertussis toxin with severe pertussis disease. Toxins (Basel) 11, 373. https://doi.org/10.3390/toxins11070373.

Scheifele, D.W., Ochnio, J.J., 2009. Module 2: diphtheria. The Immunological Basis for Immunization Series. Genève.

Sein, C., Tiwari, T., Macneil, A., Wannemuehler, K., Soulaphy, C., Souliphone, P., Ryburn, R., Ramirez Gonzalez, A., Watkins, M., Goodson, J.L., 2016. Diphtheria outbreak in Lao People's Democratic Republic, 2012–2013. Vaccine 34, 4321–4326. https://doi.org/10.1016/j.vaccine.2016.06.074.

Soetaert, K., Petzoldt, T., Woodrow Setzer, R., 2010. Solving differential equations in R: package deSolve. J. Stat. Softw. 33, 1–25. https://doi.org/10.18637/jss.v033.i09.

Soubeyrand, B., Plotkin, S.A., 2002. Antitoxin vaccines and pathogen virulence. Nature 417, 609–610. https://doi.org/10.1038/417609b.

Strive, T., Elsworth, P., Liu, J., Wright, J.D., Kovaliski, J., Capucci, L., 2013. The non-pathogenic Australian rabbit calicivirus RCV-A1 provides temporal and partial cross protection to lethal rabbit haemorrhagic disease virus infection which is not dependent on antibody titres. Vet. Res. 44, 51. https://doi.org/10.1186/1297-9716-44-51.

Talic, S., Shah, S., Wild, H., Gasevic, D., Maharaj, A., Ademi, Z., Li, X., Xu, W., Mesa-Eguiagaray, I., Rostron, J., Theodoratou, E., Zhang, X., Motee, A., Liew, D., Ilic, D., 2021. Effectiveness of public health measures in reducing the incidence of covid-19, SARS-CoV-2 transmission, and covid-19 mortality: systematic review and meta-analysis. BMJ, e068302. https://doi.org/10.1136/bmj-2021-068302.
Tiwari, T.S.P., Wharton, M., 2012. Diphtheria toxoid. Vaccines. Elsevier, Amsterdam, pp. 153–166.

Truelove, S.A., Keegan, L.T., Moss, W.J., Chaisson, L.H., Macher, E., Azman, A.S., Lessler, J., 2020. Clinical and epidemiological aspects of diphtheria: a systematic review and pooled analysis. Clin. Infect. Dis. 71, 89–97. https://doi.org/10.1093/cid/ciz808.

Walter, A., Lion, S., 2021. Epidemiological and evolutionary consequences of periodicity in treatment coverage. Proc. R. Soc. B 288, 20203007. https://doi.org/10.1098/rspb.2020.3007.

Warfel, J.M., Zimmerman, L.I., Merkel, T.J., 2014. Acellular pertussis vaccines protect against disease but fail to prevent infection and transmission in a nonhuman primate model. Proc. Natl. Acad. Sci. 111, 787–792. https://doi.org/10.1073/pnas.1314688110.

Wendelboe, A.M., Van Rie, A., Salmaso, S., Englund, J.A., 2005. Duration of immunity against pertussis after natural infection or vaccination. Pediatr. Infect. Dis. J. 24, S58–S61. https://doi.org/10.1097/01.inf.0000160914.59160.41.

White, P.J., Norman, R.A., Hudson, P.J., 2002. Epidemiological consequences of a pathogen having both virulent and avirulent modes of transmission: the case of rabbit haemorrhagic disease virus. Epidemiol. Infect. 129, 665–677. https://doi.org/10.1017/s095026880200777x.

Williams, M.M., Sen, K., Weigand, M.R., Skoff, T.H., Cunningham, V.A., Halse, T.A., Tondella, M.L., CDC Pertussis Working Group, 2016. Bordetella pertussis strain

lacking pertactin and pertussis toxin. Emerg. Infect. Dis. 22, 319–322. https://doi.org/10.3201/eid2202.151332.

Williams, W.W., Lu, P.-J., O'Halloran, A., Kim, D.K., Grohskopf, L.A., Pilishvili, T., Skoff, T.H., Nelson, N.P., Harpaz, R., Markowitz, L.E., Rodriguez-Lainz, A., Fiebelkorn, A.P., 2017. Surveillance of vaccination coverage among adult populations — United States, 2015. MMWR Surveill. Summ. 66, 1–28. https://doi.org/10.15585/mmwr.ss6611a1.

World Health Organisation, 2017. Diphtheria vaccine: WHO position paper. Wkly. Epidemiol. Rec. 92, 417–435.

Yasunaga, H., Shi, Y., Takeuchi, M., Horiguchi, H., Hashimoto, H., Matsuda, S., Ohe, K., 2010. Measles-related hospitalizations and complications in Japan, 2007–2008. Intern. Med. 49, 1965–1970. https://doi.org/10.2169/internalmedicine.49.3843.

Zasada, A.A., 2013. Nontoxigenic highly pathogenic clone of Corynebacterium diphtheriae, Poland, 2004–2012. Emerg. Infect. Dis. 19 https://doi.org/10.3201/eid1911.130297.

Zasada, A.A., Rzeczkowska, M., 2019. Nontoxigenic Corynebacterium diphtheriae Infections. Europe. Emerg. Infect. Dis. 25, 1437–1438. https://doi.org/10.3201/eid2507.180995.

Zou, R., Duan, X., Han, Z., Lu, Y., Ma, K., 2023. What information sources can prevent the epidemic: local information or kin information? Chaos, Solit. Fractals 168, 113104. https://doi.org/10.1016/j.chaos.2023.113104.