## Abstract

In Italy, during the course of the past century to the present-day, measles incidence underwent a remarkable decreasing trend that started well before the introduction of the national immunization programme. In this work, we aim at examining to what extent both the demographic transition, characterized by declining mortality and fertility rates over time, and the vaccination programme are responsible for the observed epidemiological pattern. Making use of a non-stationary, age-structured disease transmission model, we show that in the pre-vaccination era, from 1901 to 1982, the decline in birth rates has resulted in a drastic decrease in the effective transmission rate, which in turn has determined a declining trend of measles incidence (from 25.2 to 10.3 infections per 1000 individuals). However, since 1983, vaccination appears to have become the major contributing factor in the decrease of measles incidence, which otherwise would have remained stable as a consequence of the nearly constant birth rates. This led to a remarkable decrease in the effective transmission rate, to a level well below the critical threshold for disease persistence. These findings call for the adoption of epidemiological models, which deviate the age structure from stationary equilibrium solutions, to better understand the biology of infectious diseases and evaluate immunization programmes.

## 1. Introduction

Despite a long history of immunization programmes, measles is still circulating in several regions of Europe, mainly in the form of recurring spatially localized epidemics [1–3]. This poses serious concerns about the possibility of reaching the goal set by the Member States in the WHO Europe Region of eliminating measles by 2015 [4] and it calls for a deeper understanding of measles transmission dynamics.

A lot of work has been done towards understanding some common patterns of measles, such as seasonality and cyclic epidemics [5–11], disease persistence [6,11], transmission potential [12,13], relation between fertility and epidemic size [7,14] and the effects of spatial heterogeneity in the birth rates [11].

The aim of this work is to understand and quantify the roles of the demographic transition and of the implemented routine vaccination programme in shaping the observed pattern of measles incidence over time. Qualitatively, the potential impact of these factors on measles incidence could be inferred theoretically [15]. This study, however, represents the first attempt to decipher and quantify their relative weights.

The proposed investigation is motivated by the analysis of historical incidence of measles in Italy, which shows an almost continuous decreasing trend from the beginning of past century–well before the introduction of the national immunization programme in 1983–to the present day. Specifically, here we aim to answer the following two questions: (i) which are the causes of the observed decreasing trend of incidence before the introduction of measles vaccine in 1983? In fact, while it is well known that changes in birth rate reflect on disease incidence—over the past century, the Italian population underwent dramatic demographic changes—it is not clear whether the observed pattern is fully imputable to demographic transition; (ii) what is the net effect of vaccination on measles incidence? In fact, it is not clear whether the observed decreasing trend of incidence since 1983 is solely ascribable to the introduction of measles vaccine.

To answer these questions, we abandon the classic paradigm of constant population and stable age distribution over time in favour of a non-stationary, age-structured mathematical model of disease transmission, informed with longitudinal demographic data. Calibration and validation of the proposed non-stationary epidemiological model are critical issues. The availability of serological data collected both before (in 1978) and after (in 1997 and 2004) the onset of the vaccination programme [16–18], and the availability of historical data on measles incidence from 1925 to 2009, make Italy particularly suitable for such a study.

## 2. Material and methods

To assess the relative weights of the demographic transition and of the vaccination programme on measles incidence in Italy during the course of past century, a deterministic age-structured model was developed—with the population grouped into 1-year age classes—which is capable of describing both demographic changes in the Italian population as well as measles transmission dynamics from 1901 to 2009.

The model is informed with longitudinal data on age-specific mortality rates, birth rates and migration flows. All these factors are essential to reproduce the demographic structure of the Italian population over time. Specifically, realistic mortality and birth rates are necessary to reconstruct the age structure of the population over time, and thus the pattern that characterized the demographic transition. Migration flows are necessary to obtain realistic predictions of population size over time—the number of Italian emigrants from 1895 to 1985 is estimated to be about 23 million and the number of immigrants from 2002 to 2009 is estimated to be about 3.5 million. Birth- and age-specific mortality rates are derived from census data (provided by the Italian Institute of Statistics), collected every 10 years; for intercensual years, we made a linear interpolation between the data from two consecutive censuses. Migration flows are estimated by combining estimates of the yearly number of migrants (provided by the Italian Institute of Statistics) and recent estimates of the age distribution of migrants [19]. Demographic data are shown in the electronic supplementary material.

To calibrate and validate the model, we make use of the following epidemiological data sets: serological data coming from three studies conducted in 1978 [16], 1997 [17] and 2004 [18]; yearly incidence from 1925 to 2009 (provided by the Italian Institute of Statistics); data on average age at infection in different years and age distribution of cases for the year 2009 [3].

In order to model the disease transmission process, the population of each age class *i* is divided into five epidemic classes: susceptible (*S _{i}*(

*t*), individuals who can acquire infection), susceptible vaccinated (

*F*(

_{i}*t*), individuals who are still susceptible after having received one or two vaccine doses), latent (

*L*(

_{i}*t*), individuals who have been infected and are not able to transmit the pathogen yet), infectious (

*I*(

_{i}*t*), individuals who have been infected and are able to transmit the pathogen) and removed (

*R*(

_{i}*t*), individuals who are immune to the disease, either because of vaccination or because they recovered from infection). The mean durations of both the latent (1/

*ω*) and the infectious (1/

*γ*) period are fixed at 7 days [12], whereas immunity is assumed to be lifelong.

As baseline choice, although it has largely been shown in the literature that mixing patterns are assortative by age [20,21], we make the simplest assumption of homogeneous mixing, i.e. all individuals in the population are exposed to the same force of infection. The main reason that led us to make this simplistic assumption is that, given the dramatic social and demographic changes observed during the past century (e.g. in household size and composition, school attendance rate and urbanization level), our knowledge on current age-mixing pattern cannot be considered representative of the whole studied period and no quantitative historical data on mixing pattern by age are available.

We do not assume any seasonal forcing even though within-year measles dynamics is well documented and follows a clear seasonal pattern, triggered by alternation of school term periods and summer holidays. Again, this assumption was made because, to our knowledge, current seasonal forcing cannot be considered representative of the whole period. Indeed, school attendance rates have dramatically increased during the past century and school term periods may have changed (e.g. because of changes in school organization and wars).

Although heterogeneous mixing and seasonality rate are not considered in the baseline model, the sensitivity of our results to different assumptions on contact patterns (in particular, we use age-specific contact matrices reported in [20] and [21]) and to seasonal variation of the transmission rate was considered.

The model accounts for vaccination of 1-year-old children and a second booster dose given to children at the age of 5 years. There are two parameters regulating the vaccination process: (i) vaccine efficacy, i.e. the probability that an individual becomes immune after the administration of a vaccine dose (which is assumed to be *f* = 95% as reported in the literature; see, for instance, [12]); and (ii) vaccination coverage, i.e. the fraction of individuals of a given age who are vaccinated for the first time, *v*_{1}(*t*), or receive the second dose, *v*_{5}(*t*). These data were obtained from the Italian Ministry of Health and the Italian Institute of Statistics and they are shown in the electronic supplementary material.

The model has only one free parameter, namely the disease-specific transmission rate *β*, and consists of the following set of ordinary differential equations:where *b*(*t*) is the yearly birth rate, *μ _{i}*(

*t*) are yearly age-specific mortality rates,

*N*(

*t*) represents the total population,

*n*= 100 years is the maximum age considered in the model and

*δ*is the Dirac delta function.

_{ij}As the age structure of the Italian population in 1901 can be considered in equilibrium (see the electronic supplementary material), initial conditions of the system in 1901 are defined as the equilibrium solution of the system with birth and mortality rates by age fixed to those observed in 1901 and a certain transmission rate *β*, and initialized with 10 infected individuals.

Posterior distribution of *β* was explored by Markov chain Monte Carlo (MCMC) sampling applied to the likelihood of serological data coming from two studies conducted in the pre-vaccination (1978) and vaccination (2004) era, respectively. A log scale was used for sampling as parameter *β* is positive definite [22].

The fundamental quantity regulating disease dynamics is the basic reproduction number (*R*_{0}), which represents the average number of secondary infections generated by a typical index case during the entire period of infectiousness [12]. *R*_{0} can be computed analytically for this model as the dominant eigenvalue of the next generation matrix associated to the transmission model [23], and can be approximated by *R*_{0} = *β*/*γ*.

Starting from results presented in [7,24], it is straightforward to show that changing the birth rate by a given factor *ε* produces exactly the same dynamical transitions as changing the transmission rate by the same factor, that is, changing birth rates from *b* to (1 + *ɛ*)*b* is equivalent to changing the transmission rate from *β* to *β*(1 + *ɛ*). This argument can be readily generalized to the case where a fraction *v* of newborns is vaccinated with vaccine efficacy *f*. Therefore, we define *β*(1 + *ɛ*)(1 − *fv*) to be the effective transmission rate of the disease. We analyse changes in measles epidemiology in terms of changes of the effective transmission rate.

Details on model formulation, parameterization and calibration can be found in the electronic supplementary material.

## 3. Results and discussion

### (a) Validation of the demographic model

A prerequisite for analysing the effects of demographic transition on measles incidence is the availability of a reliable and robust demographic model. The analysis of historical demographic data shows that the age distribution of the Italian population before the 1901 census was in a stationary equilibrium (see the electronic supplementary material). Afterwards, at least four main events have characterized the demographic dynamics in Italy: the drop in birth rates during the two world wars, the ‘baby boom’ in the 1950s, and another prolonged drop in birth rate since the seventies. These changes, coupled with the drop in mortality rates, especially in children and the elderly, have led to the progressive aging of the Italian population. Only in recent years has immigration contributed to sustaining population growth. As shown in figure 1*a*, the population model, informed with longitudinal data on birth rates, age-specific mortality rates and emigration/immigration rates, well explains the demographic dynamics of the Italian population since the beginning of the past century. The agreement with available historical data on the age structure of the Italian population is excellent (figure 1*b*) and the coefficient of determination *R*^{2} is larger than 0.99 (*p*-value < 0.0001) in all tested years. In the electronic supplementary material, we compare observed and predicted population size over time and analyse the effects of emigration and immigration fluxes.

### (b) Calibration of the disease transmission model

Although the model has only one degree of freedom, the fit to serological data is excellent (figure 2*a*,*b*): indeed, for the majority of data points, the observed average seroprevalence falls in the 95% CI of model predictions and for all age groups observed, and predicted 95% CI intersects. The transmission rate was estimated to be *β* = 677.9 per year (95% CI: 657.7–698.1 yr^{−1}), resulting in *R*_{0} = 13 (95% CI: 12.6–13.4) in 1901. This estimate is in quite good agreement with the estimated range 11–12 obtained in [12] by analysing historical measles data for England and Canada.

Remarkably, in the electronic supplementary material, we show that the classic equilibrium equations (stationary age structure and endemic equilibrium) produce a worse fit to the 1978 serological data (2004 data are not suitable for this analysis as they are affected by the effects of time-varying vaccination) with respect to the proposed non-stationary model. In fact, in general, serological data depend on a long history of infections, and, specifically for the 1978 serosurvey data, on infections occurred from 1965 to 1978. This period was characterized by one of the most remarkable drops in birth rate in the history of the Italian population (about 35% from 1965 to 1978) and, consequently, a proportional drop in the effective transmission rate. Therefore, older age classes (e.g. 1965–1972 birth cohorts) have experienced a higher force of infection, at least in their first years of life, with respect to younger age classes (e.g. 1973–1978 birth cohorts). This is the reason why classical equilibrium equations, not accounting for demographic changes, underestimate seroprevalence in older age classes and overestimate seroprevalence in younger age classes, thus resulting in an overall worse fit. All in all, this supports the use of the proposed model to analyse the role of demographic transition on measles incidence.

### (c) Validation of the disease transmission model

Once parameterized, we validate the model against a set of data not used for model calibration.

In the first half of the past century, measles incidence in Italy underwent a dramatic decreasing trend, even well before the introduction of vaccination (figure 3*a*). Although incidence data reported to the statutory notification system suffer from underreporting (as well as being variable over time and by region [16]) and have a crude temporal resolution (yearly data), a comparison between the incidence predicted by the model, and that reported to the statutory notification system, adjusted for underreporting, shows a high correlation (Spearman's rank correlation: *ρ* = 0.99, *p* < 0.0001; coefficient of determination: *R*^{2} = 0.94, *p*-value < 0.0001; figure 3*a*). Both the observed and the predicted incidence were averaged over a moving window of 15 years. Estimates of underreporting (7.4%) were obtained by analysing monthly incidence data from 1980 to 1996 (details in the electronic supplementary material). Poor validation before the 1950s can be explained by contact structure and seasonality being different from recent years and on which we have no information. By restricting the comparison to years after 1950, the coefficient of determination *R*^{2} increases to 0.99, *p*-value < 0.0001.

Estimates of the average age at infection are available only in three time windows, all lying in the vaccination era. Whenever such data are available, model predictions for the whole of Italy lie inside the minimum and maximum average age at infection observed in the 21 Italian Regions (figure 3*b*). It is worthwhile to stress that the large regional variability may depend on the observed different levels of vaccine coverage [25]. Moreover, our results show a bump in the average age at infection during the World War I, which is in line with the recorded historical data on measles mortality [26].

The model was also validated against serological data collected in 1997 [17]; the comparison between observed and predicted seroprevalence by age groups is shown in figure 3*c*. The model slightly underestimates seroprevalence in children aged 2–4 and 10–14 years. However, for all other age groups, the observed and the predicted 95% CI intersect.

Detailed age-specific incidence data are available for 2009: as shown in figure 3*d*, age distribution of cases peaked in adults (individuals aged 25–44 years) both in model predictions and in the surveillance data, and the predicted distribution of cases by age is in good agreement with the one observed (Spearman's rank correlation: *ρ* = 0.98, *p* = 0.0004; coefficient of determination: *R*^{2} = 0.97, *p*-value < 0.001).

### (d) Assessing the relative weights of demographic transition and vaccination in the decrease of measles incidence

The observed trend of the birth rate over time yields an estimate of an almost always decreasing pattern in the effective transmission rate during the pre-vaccination era. Specifically, the effective transmission rate drops from 677.9 per year (95% CI: 657.7–698.1 per year) in 1901 to 236.4 per year (95% CI: 229.3–243.4 per year) in 1982, almost three times lower than in 1901 (figure 4*a*). The introduction of vaccination in 1983 contributed to accelerate this decreasing trend (figure 4*a*). In particular, together with the effect of demographic changes, the immunization programme contributed to lowering the effective transmission rate so that the threshold for measles elimination was reached (figure 4*a*). Without the introduction of the immunization programme in 1983, the effective transmission rate would have remained nearly constant to those values that were predicted for the early 1980s (figure 4*a*).

Figure 4*b* shows the effects of demographic transition and vaccination on measles incidence. The average incidence was estimated to be 25.2 (95% CI: 25.1–25.3) per 1000 individuals in 1901 and dropped to 10.3 (95% CI: 10.2–10.4) per 1000 individuals in 1982. Without the introduction of measles vaccine, in the past 30 years, the incidence would have remained nearly constant at values observed in the early 1980s. We also found a strong positive correlation between birth rate and both the predicted and observed measles incidence (Spearman's rank correlation: *ρ* = 0.97, *p* < 0.0001 and *ρ* = 0.88, *p* < 0.0001, respectively). In addition, both predicted and observed measles incidence appear to be strongly positively correlated with the effective transmission rate (Spearman's rank correlation: *ρ* = 0.97, *p* < 0.0001 and *ρ* = 0.89, *p* < 0.0001, respectively).

We estimated that the number of avoided cases per vaccine dose (i.e. the difference between the number of cases predicted by the model without and with vaccination, divided by the number of vaccinated individuals, cumulated over time) readily increased over time and nowadays is close to 1 (figure 4*c*), strongly supporting the cost-effectiveness of measles vaccination. Our estimates comply with those reported in [27] where authors estimate the number of avoided infections per vaccine dose to be about 0.11 in Palermo (Sicily) from 1993 to 1997. Such a low value might be a consequence of the low vaccine coverage and force of infection (as a consequence of demographic factors) over that period.

As highlighted in figure 4*d*, a remarkable increase in the ages of cases over time is predicted in the pre-vaccination era: in particular, the average age at infection is estimated to have increased from 2.8 years (95% CI: 2.7–2.9 years) in 1901 to 5.9 years (95% CI: 5.8–6.1 years) in 1982. By looking at the vaccination era, we found a drastic increase in the average age at infection: the estimated average age at infection in 2009 was 17 years (95% CI: 16.7–17.3 years).

The impact of decreasing incidence and increasing age at infection in shaping age-specific measles seroprevalence can be appreciated by looking at figure 4*e* where, according to model simulations, we predict a progressive reduction of seropositive individuals in children and adolescents up to the start of the national measles immunization programme in 1983. This decline is well marked especially in individuals aged less than 10 years; in particular, as shown in figure 4*f*, during the pre-vaccination era variations of up to 50% in seropositive individuals are predicted in children aged 0–3 years and up to 40% in the other age classes. In later years, the effect of the immunization campaign prevails over that of demographic transition in shaping seroprevalence profiles by age. In particular, by looking at predictions for 2009, it emerges that low vaccine coverage from 1983 to the end of past century has generated a protection gap particularly relevant among individuals aged 10–30 years (less than 80% of individuals protected; see figure 4*e*) who are currently the most affected age group.

### (e) Sensitivity analysis

The effects of different assumptions on mixing patterns by age on model results were evaluated. The model was calibrated with two different contact matrices by age available for Italy, namely the Polymod matrix [20] and the synthetic contact matrix reported in [21]. Overall, results show that mixing patterns are likely to have changed so much during the course of the past century that our knowledge on the current age-mixing pattern cannot be considered representative for the whole period. Nonetheless, the use of contact matrices may contribute to improving predictions for the most recent period, characterized by higher average age at infection, and thus contacts between individuals in different age groups become more relevant. Details on this topic can be found in the electronic supplementary material.

It is well documented that measles dynamics follows a clear seasonal pattern. Our analysis, based on the analysis of monthly incidence from 1980 to 1996, supports the already proposed hypothesis (see for instance [7,8]) that this pattern is mainly driven by the school calendar. However, we also show that the within-year pattern does not influence the global trend of disease incidence, whose analysis is the main aim of this work. We also found that the seasonal reporting rate may vary between 1.5 and 25%, a result consistent with available estimates on its geographical variability [28]. An in-depth discussion on this topic can be found in the electronic supplementary material.

We found that results are insensitive to length and distribution of latent and infectious periods, unless we consider unreasonably long durations. Details on this topic can be found in the electronic supplementary material.

## 4. Conclusion

In this work, we aimed to answer the following two questions: (i) what are the causes of the observed decreasing trend in the incidence of measles cases before the introduction of measles vaccine in 1983? (ii) What is the net effect of vaccination on measles incidence? For the first question, our results clearly support the hypothesis that demographic transition is mainly responsible for the observed decreasing trend in incidence before the introduction of the measles vaccine. For the second question, our results show that incidence would have remained almost constant since the early 1980s as a result of nearly constant birth rates. Thus, our results support the hypothesis that the decreasing trend of incidence observed in the vaccination era is mainly attributable to the direct effect of vaccination.

More in detail, we found that demographic transition is responsible for a threefold decrease in the effective transmission rate during the pre-vaccination era, which in turn contributed to a proportional decrease in measles incidence (from about 25.2 per 1000 individuals in 1901 to 10.3 per 1000 individuals in 1982). Demographic transition also determined a drastic increase in the average age at infection (from 2.1 years in 1901 to 5.5 years in 1982) during the pre-vaccination era and remarkable differences (up to 50% in young age classes) in terms of measles seroprevalence by age. As for the vaccination era, we found that the number of avoided infections per vaccine dose readily rose to values close to 1 and vaccination greatly contributed to making the final step in decreasing the effective transmission rate to values below the critical threshold for disease persistence. Conversely, the average age at infection increased to values between 15 and 20 years in 2009 and, as an effect of initially low vaccine coverage before 2000, a protection gap particularly relevant in individuals aged 10–30 years (less than 80% of individuals protected) has emerged. These results support the hypothesis that measles has become no more endemic in Italy and explain why its epidemiology is characterized by series of sporadic spatially localized outbreaks occurring in age classes (mostly 10–30 years old) that are still unprotected.

From the methodological point of view, it is interesting to observe how (unexpectedly) well a non-stationary age-structured model, informed with longitudinal demographic and vaccination data and with only one free parameter (namely, the disease-specific transmission rate), is able to describe the observed measles epidemiology. However, in order to predict the future dynamics of measles and to evaluate the effectiveness of immunization programmes, the model needs to be improved in several aspects. First of all, demographic projections are required, accounting for changes in birth rates, immigration trends and mortality rates by age. Moreover, as highlighted by the conducted sensitivity analysis, the model should be enriched by accounting for structured contacts [29–31], for instance by age [20,21]. In fact, homogeneous mixing models, such as the one introduced in this paper, tend to overestimate the overall number of cases [32]. While this may not be an important issue when analysing historical data, where cases occur almost only in children and the majority of still susceptible children get infected during the same epidemic season, accounting for heterogeneous contacts by age (or social context) becomes more relevant during the disease elimination process, which is instead characterized by smaller numbers of cases, shift of the age at infection from childhood to older ages, spatially localized outbreaks and small clusters of secondary cases in high schools and universities. Another important feature that is needed to decipher the impact of demographic transition and immunization programmes, but that is missing in this model, is the spatial structure. Accounting for differences in population density [10] and patterns of human mobility [33–35] could help to explain large regional differences in terms of case incidence and average age at infection. However, retrieving longitudinal/historical mobility data to determine the strength of the connections between spatially distinct populations is a very hard task.

In general, our results recommend the introduction of a new class of models that deviate the age distribution from stationary equilibrium solutions by accounting for demographic changes. This will make it possible to deeply understand epidemiologic features of infectious diseases and, from the public health perspective, to better evaluate the short- to medium-term effects of mass vaccination programmes for endemic diseases, notably in developing countries that are undergoing dramatic demographic changes.

## Funding statement

This work was partially supported by the EU FP7 Epiwork project for research funding. The funders had no role in study design, data collection and analysis, decision to publish or preparation of the manuscript.

## Acknowledgements

The authors thank the anonymous reviewers for their comments, Laura Fumanelli and Alessia Melegaro for comments on an early version of the manuscript.

- Received October 14, 2013.
- Accepted November 26, 2013.

- © 2014 The Author(s) Published by the Royal Society. All rights reserved.