A predictive study was conducted based on an analysis of all clinically and laboratory-confirmed measles cases reported by health districts in Côte d’Ivoire between 2010 and 2024.

Côte d’Ivoire is located in western sub-Saharan Africa and covers an area of 322 462 km². It is bordered by Burkina Faso and Mali to the north, Liberia and Guinea to the west, Ghana to the east and the Gulf of Guinea to the south. The country has four seasons, namely a long rainy season from March to June, a short dry season from July to August, a short rainy season from September to October and a long dry season from November to February. The total population living in Côte d’Ivoire is estimated to be 29 389 150 inhabitants, with a density of 91 inhabitants/km² and an average annual growth rate of 2.9%.¹⁰ The health system comprises 33 health regions and 113 health districts.

Measles surveillance data in Côte d’Ivoire from 2010 to 2024 were used to adjust the ARIMA model. To record these data, health centres report suspected measles cases to health districts, which in turn notify the central body (Directorate for the Coordination of the Expanded Program on Immunisation [DC-PEV]). Laboratory confirmation of the diagnosis is made by the Pasteur Institute of Côte d’Ivoire (IPCI) based on the detection of immunoglobulin M (IgM). All these data are stored in the Epi.info database of the DC-PEV and are validated through harmonisation sessions with different structures of the health system (the DC-PEV, the National Institute of Public Hygiene [INHP]), the IPCI and the Directorate of Health Information (DIS):

Raw data on measles surveillance in Côte d’Ivoire are collected daily by health districts. These data were aggregated on a half-yearly basis by grouping them by year and half-year and summing up the number of measles cases. This choice was made in order to smooth out daily fluctuations and, above all, the seasonal dynamics of measles in Côte d’Ivoire. The time series used for model adjustment is therefore semi-annual. The time series used for model adjustment is therefore semi-annual.

An exploratory analysis of the data was performed by visualising the trend and seasonality with graphs (time series curves, Autocorrelation Function [ACF] and Partial Autocorrelation Function [PACF]) and checking the stationarity of the series with the Dickey–Fuller test¹¹ after transforming the data into a time series using the ‘ts()’ function in R software.

This time series of measles cases was decomposed using the R function ‘decompose’, making it possible to identify the trend and seasonality. The number of differentiations was automatically calculated using the ndiffs and nsdiffs functions for non-seasonal and seasonal differentiations, respectively. It was subsequently applied two differentiations to remove the trend and retested the stationarity on the differentiated series using the Dickey–Fuller test (significance threshold 5%).

It was adjusted for seasonality (since seasonal cycles were present) by automatically identifying the best^12,13 ARIMA model. Model stationarity was checked using model residuals, and the best-fitting model was automatically identified using the auto.arima function.

As the model identified by auto.arima does not take double differentiation into account, we tested a set of models that do take double differentiation into account. The model with the lowest Root Mean Squared Error (RMSE) and whose Akaike Information Criterion (AIC) and residuals were white noise was selected.

The model was validated using a cross-validation approach adapted to time series, in particular, rolling-origin cross-validation, which was used to assess the model’s predictive performance. This method evaluates the stability of the model in different temporal contexts, such as the coronavirus disease 2019 (COVID-19) pandemic period.

The model’s performance was evaluated using RMSE, Mean Absolute Error (MAE) and Theil U test compared to a naive seasonal reference (snaive).

Calculating accuracy metrics between the training data and the validation data allowed us to evaluate the accuracy of the model based on the smaller values with the validation data test than the training data test (RMSE, MAE and Theil U test). The Theil inequality coefficient is most often used for model validation by evaluating the accuracy of the model’s prediction. It is a coefficient whose values range from 0, meaning a perfect prediction, to 1, representing maximum inequality.¹⁴

It was conducted to assess the robustness of the model in the face of potential structural breaks in the time series. This involved re-estimating the model over different periods, with a particular focus on the period of the COVID-19 pandemic.

A seasonal ARIMA model was estimated over the entire series by including a dummy variable coded 0 before rupture and 1 after rupture. This variable was introduced as an external regression (xreg) to represent the level change.

Based on the selected and validated ARIMA model, we made half-yearly predictions on the number of measles cases from 2025–2030 in Côte d’Ivoire, taking into account the uncertainty of the forecasts (confidence intervals). Taking into account the population data, from the general population census of Côte d’Ivoire, and especially on the projections of this population (2025–2030) transmitted by ANSAT to the Ministry of Health, both for the general population and for children aged 0–15 years.¹⁰ For these years (2025–2030), we calculated the annual incidence of measles. These predictions were made for the general population and for children aged 0–15 years.

All analyses were performed using R.4.2.1 software, including the prediction of the number of measles cases using the R ‘forecast’ function.¹⁵

This study used routinely collected aggregate surveillance data that did not contain any personally identifiable information. As such, the analysis posed minimal risk to individuals. Data were processed in accordance with national regulations and institutional guidelines on data protection and confidentiality. Given the use of anonymised secondary data, no formal informed consent was required. Ethical approval to conduct this study was obtained from the Research Ethics Committee of the National Institute of Public Health. The ethical clearance number is 001-2025.

The decomposition of the time series generated three main components, namely an ascending trend from 2017, seasonality with periodic fluctuations on a biannual basis and high residuals in 2020, 2021 and 2022, suggesting unusual events not explained by the other components. These findings are the same as the complete data (Figure 1a) and those of children aged 0–15 years (Online Appendix 1 Figure OA1-2a).

Lags 1 and 2 show significant positive autocorrelations, with the decrease in autocorrelations after a few lags indicating that the series may be stationary (Figure 1b and Figure 4b).

Automatic identification of the number of differentiations noted one non-seasonal differentiation (d = 1) and one seasonal differentiation (d = 1).

It should be noted that the series required two differentiations to be stationary (test of Dickey-Fuller, p = 0.03).

The ACF and PACF of the residuals show that most of the autocorrelations are not significant (Figure 1c).

Some residual correlations at specific lags might suggest incomplete adjustment; however, the Ljung–Box test indicates that the residuals do not exhibit significant correlation (p = 0.13).

It initially used the auto.arima function to select the model. However, the model proposed by auto.arima() (ARIMA(0,0,1)(1,1,0)[2]) did not take into account the two differentiations. We then compared five models based on double differentiation and low AR and MA orders (between 0 and 2) to avoid overfitting. (Online Appendix 1 Table OA1-1).

The ARIMA(1,1,1)(0,1,1)[2] model was chosen based on its smaller RMSE.

After adjusting the model over the entire series, there is a significant ACF peak (Online Appendix 1 Figure OA1-1). We then decided to increase q by +1 in order to capture this residual dependence.

Finally, it was applied this ARIMA(1,1,2)(0,1,1)[2] model to the entire data series.

The selected SARIMA model (ARIMA(1,1,2)(0,1,1)[2]) is good, considering the following observations:

Indeed, the model residuals do not show significant autocorrelation as indicated by the Ljung-Box test (p = 0.05698 (> 0.05)). The residuals are therefore similar to white noise, which indicates that the model captures the main structures of the data well.

The chosen model has a simple structure: ARIMA(1,1,2)(0,1,1)[2] where (1,1,2) indicates a first-order moving average in non-seasonal data with a non-seasonal difference of order 1 and (0,1,1)[2] indicates a first-order seasonal autoregressive component with a seasonal difference of order 1 and a semi-annual seasonality ([2]).

Both predicted and observed data diverge with seasonal variations. Moreover, observed data remain within the confidence interval of predicted data (Figure 2).

For the training data, the RMSE and MAPE are respectively 215.8% and 107.6%. Similarly, the results for the test data are respectively 169.2% and 22.2%, and the proposed model performs better than the basic one.

Figure 2 compares observed values with forecasts validated by cross-checking at one step. Most forecasts fall within a tolerance band of ±15% during stable periods (before 2018), while deviations increase during periods of high volatility, particularly from 2020 onwards, indicating reduced accuracy for sudden peaks.

Furthermore, the higher RMSE (approximately 290 cases) compared to the MAE (185 cases) reveals the presence of significant errors visible during the sharp peaks after 2020, highlighting the model’s difficulty in anticipating outbreaks.

It was therefore checked for structural breaks using the strucchange() package before performing an ARIMAX analysis.

This check noted that the second half of 2020 was the main structural break with the smallest Bayesian Information Criterion (BIC) value (413.1) (Online Appendix 1 Table OA1-2).

The inclusion of a break indicator variable in the ARIMA model (ARIMA(1,1,2)(0,1,1)[2] + break dummy) made it possible to capture a significant structural change that occurred in the second half of 2020. This coefficient associated with the intervention (xreg = 607.7) indicates an increase in the average number of measles cases per half-year after that date. Furthermore, the inclusion of this break improved the model’s performance (RMSE fell from 290 to 170, a reduction of more than 40%). Furthermore, the Ljung-Box test had a p-value of 0.2055; the residuals are therefore white noise. This confirms the relevance of the model with intervention in this time series of measles data.

However, in cross-validation, the naive seasonal benchmark outperformed the ARIMA model with intervention (Theil’s U = 1.28) (Online Appendix 1 Table OA1-3), indicating superior short-term predictive accuracy. Nevertheless, the ARIMA model was selected for inference and scenario analysis because it explicitly accounts for the structural break identified in 2020 and provides interpretable estimates of post-break level changes.

After 2025, the model predicts periodic fluctuation around a stable mean for the number of cases. The confidence interval increases over time, reflecting increasing uncertainty in predictions as the time horizon extends (Figure 3).

The final fitted model was used to generate forecasts for the period 2025–2030.

It is noted that the predicted incidence of measles will remain high in the general population, even if a downward trend will be observed. Indeed, the incidence will vary from 51.5 to 45.6 cases per million from 2025 to 2030 in the Ivory Coast (Table 1).

In children aged 0–15, the incidence will be higher in 2025 by around 146 cases per million and will increase by 18 points to 160 cases per million in 2030 (Table 2).

This study has certain limitations that must be taken into account when interpreting projections for the elimination of measles in Côte d’Ivoire by 2030.

Firstly, the analysis is based on systematically collected national surveillance data, which may be affected by under-reporting and differences in case detection between health regions and over time.

Secondly, the chronological approach assumes that historical trends, seasonality and the identified structural break will persist in the future. However, changes in the performance of vaccination programmes, the implementation of additional vaccination activities, population mobility or the capacity to respond to epidemics could significantly alter the future dynamics of measles.

Thirdly, the models do not incorporate certain key variables such as national vaccination coverage, disruptions to the health system, or behavioural factors influencing vaccine uptake. It should also be noted that the 95% confidence intervals for incidence remain wide, reflecting considerable uncertainty about future transmission patterns. Therefore, the results should be interpreted as scenario-based projections to support strategic planning and programme decision-making, rather than definitive predictions about measles elimination. Continued strengthening of surveillance systems and routine immunisation remains essential to validate and update these projections.

To reduce the impact of measles and achieve its elimination in the general population, and specifically in children aged 0 to 15 years, three levels of intervention can be considered, namely increasing vaccination coverage, strengthening surveillance and conducting certain studies; in particular, studies evaluating current vaccination strategies and the age of measles vaccination in Côte d’Ivoire.

Increasing measles vaccination coverage should be effective through the implementation of routine vaccination that guarantees the administration of both doses of measles vaccine and targeted campaigns, reaching rural areas and vulnerable populations. This requires controlling the denominator through well-conducted target enumeration. Strengthening the immunity of the most at-risk age group is an imperative that requires the organisation of a multi-age measles vaccination campaign that includes children aged from 9 months to 15 years.

Also, for these different activities, including awareness raising, critical periods (before and during the first semester) should be targeted to take into account seasonal variations. Back-to-school periods should also be considered in the implementation of activities.

Disease surveillance – particularly for measles – will need to be strengthened to improve reporting in all health areas and on a continuous basis through strengthening health infrastructure, staff capacity and educating parents.

Finally, studies should be carried out to elucidate the socio-demographic and immunological profile of people aged over 15 years who have measles in order to take targeted action against them.

Understanding and acting on the determinants of vaccine hesitancy to improve Vaccine Coverage (VC), and understanding the perception of the importance of measles among stakeholders in order to improve the quality of surveillance.