The influence of global climate drivers on monsoon onset variability in Nigeria using S2S models

Rainfall onset has a lot of implications on the sustainability of the socio-economic activities in Nigeria. This study assesses the skills of CMA, ECMWF, and UKMO sub-seasonal-to-seasonal (S2S) models in predicting monsoon onset and its variability over Nigeria. It also investigates the global drivers modulating the variability and their teleconnections with rainfall onset anomaly. All the models, their ensemble members, and the observations were subjected to quantitative statistical analyses from 1998 to 2016. Results show that the three models are able to simulate the Northwards migration of the onset dates adequately with inherent biases and unique characteristics. They are also able to capture the evolution and variability of the global drivers modulating the monsoon onset. While CMA and the ECMWF models improve progressively towards the Sahel, the UKMO model performance is best over the Gulf of Guinea. In addition, despite the fairly poor performance of the models in predicting the variability of onset dates over the Gulf of Guinea and the Sahel, there is a considerable improvement in the correlation skill of the models over the Savannah. Furthermore, results show that only the ECMWF model was able to produce the strength of both the African Easterly Jets (AEJ) and the Tropical Easterly Jet (TEJ) in spatio-temporal mode. These are two of the crucial global drivers modulating the dynamics of West African monsoon. It was also found out that most global drivers, especially the Inter-tropical Discontinuity (ITD) and the Sea Surface Temperature (SST) over the Central Pacific, exhibit direct teleconnection with the onset anomaly. This direct relationship is shown to be strongest over both the Gulf of Guinea and the Sahel. Although the CMA model might have the least skill, it, however, showed that all the S2S models, despite the inherent biases, are able to predict rainfall onset over Nigeria, within the sub-seasonal timescale. Finally, the results show that improvements in multi-model ensembles are valuable added information able to significantly improve model performance.


Introduction
The onset of the monsoon in West Africa is considered a crucial and integral part of monsoon dynamics, because it provides valuable information for decision-makers. In Nigeria, for instance, almost 68% of the population (about 122 million people) are peasant farmers and completely depend on rains for sustainable farming (Lawal et al. 2016). Furthermore, national sectors that are related to water management such as agriculture, health institutions (i.e., managers of waterborne diseases), water resource managers (i.e., dams), and hydro-electricity generation, etc., do greatly benefit from the prediction of the dates of monsoon onset (Bombardi et al. 2017). Apart from the onset, monsoon season also consists of two other elements-the cessation dates and rain days-but it has been generally accepted that monsoon onset is the most critical for agriculture and water resource managers (Kowal and Knabe 1972;Omotosho et al. 2000). As the science of the sub-seasonal-to-seasonal (S2S) project is improving, so also are the corresponding improvements in the weather forecast skills (Olaniyan et al. 2018) and the promotion of its uptake by the National Weather Centres (NWC). The knowledge of S2S has also been declared as a new frontier in weather forecasting (Vitart et al. 2012;Vitart 2014), and one of the main objectives of the S2S project is to determine the forecasting capabilities of the characteristics of atmospheric dynamics relevant to monsoon onset using S2S models (Vitart et al. 2017). Although, S2S forecasts are not seasonal forecasts; however, studies (e.g. Vellinga et al. 2013;Alessandri et al. 2015;Bombardi et al. 2017) have shown that it is possible to predict significant rainfall onset dates within the sub-seasonal time scales. Unfortunately, the aspect of utilizing S2S models to determine the onset dates, as well as utilizing the models to understand the roles played by atmospheric dynamics in the determination of onset dates and variability, is lacking over Nigeria. Therefore, this present paper will seek to address this aspect over Nigeria.
There has been no consensus among climate scientists on the concept of determining the most appropriate rainfall onset dates. Based on various assumptions, there are more than 16 scientific concepts in determining the rainfall onset date (Fitzpatrick et al. 2015(Fitzpatrick et al. , 2016. The concepts have been classified into two basic definitions-the regional and the local definitions. Most regional definitions are based on the abrupt movement of the Inter-tropical discontinuity (ITD) from the Gulf of Guinea to the Sahel (Sultan and Janicot 2003;Fontaine and Louvet 2006;Hagos and Cook 2007;Fontaine et al. 2008;Gazeaux et al. 2011). Regrettably, most studies on the regional definition have always focused on the Sahel when applying these definitions to onset dates. The local definition is basically used in the operational activities of the NWCs including the Nigerian Meteorological Agency (NiMet: www.nimet .gov.ng). It is determined using daily precipitation data at a meteorological station with inherent ignorance of rainfall from any neighboring stations (Fitzpatrick et al. 2015). This may, however, be as a result of few and sparsely dense data sets that are normally the bane of meteorological data set over the tropical region. Because of this, the previous works (e.g. Stern et al. 1981;Silvakuma 1988;Omotosho et al. 2000;Odekunle et al. 2005) may not have taken the neighboring rainfall teleconnection into cognizance when defining onset dates.
Furthermore, in spite of the huge importance of the rainfall onset to the West African economy, long-range seasonal climate forecasts for the timing of the onset dates are not always available (Ingram et al. 2002). This may be because most global prediction centers responsible for the seasonal climate forecasts do not make it freely available to NWCs in Africa. It often costs a lot of money for the forecast to be released for real-time operational use. We, though, acknowledge that, recently in 2018, the Copernicus project (www.ecmwf .int/data-catal ogues /c3s-seaso nal) has made the seasonal climate forecast freely available, but it still very laborious to access the data by most African NWCs. Most African NWCs may not be able to have the Internet resources required for the retrieval of these data even when using MARS (Meteorological Archival and Retrieval System (Raoult 2001).
In addition, studies that relate the skill of S2S models to the prediction of rainfall onsets, over West Africa, are very scarce. This is in spite of one of such studies over South America, East Asia, and Northern Australia (Bombardi et al. 2017). Bombardi et al. (2017) evaluated the skill of S2S precipitation forecasts from three global coupled models, namely: the China Meteorological Administration (CMA) model; the Japan Meteorological Agency (JMA) model; and the National Centre for Environmental Prediction (NCEP) model. Using hybrid methodology (combination of models with observations), (Bombardi et al. 2017) concluded that, despite the existence of some inherent complexities in the skill of S2S prediction, they are, however, skillful enough to identify onset dates at sub-seasonal lead times. Though, studies have shown the predictive skills of various S2S models over different regions of the world for different purposes: e.g., Tompkins and Feudale (2010) worked over some areas in Africa; (Lynch et al. 2014) presented their work over Europe; Tompkins and Giuseppe (2015) showed some results over Uganda and Kenya; (White et al. 2015) considered S2S skill over Australia. Recently, (Olaniyan et al. 2018) evaluated the predictive skill of the European Center for Medium-range Weather Forecasting (ECMWF) S2S model over Nigeria during the peak of the boreal summer monsoon from 1998 to 2012. They showed that the model is skillful in reproducing atmospheric dynamics that are modulating rainfall during the peak of the monsoon seasons. However, they concluded that rainfall prediction by the ECMWF-S2S model, in areas south of latitude 10°N over Nigeria, is better represented by the model's ability to reproduce atmospheric dynamics rather than the actual quantitative rainfall forecasts. Yet, none of these above-mentioned studies focused on the abilities of the S2S models to predict reliable rainfall onset dates and variability within the sub-seasonal timescales, especially over Nigeria. In addition, these studies did not investigate the understanding of the teleconnection between the atmospheric dynamics, described in detail by (Olaniyan et al. 2018), and onset variability over Nigeria. This present paper will address these knowledge gaps, because the demands for reliable seasonal and/or sub-seasonal climate forecasts within a timescale that is right for effective response are tremendously growing all over the world.
Besides a few shreds of evidence of the S2S models' abilities in predicting rainfall onset dates over different regions using different methods, there are still some concerns that are needed to be addressed based on the performance of S2S models in predicting monsoon onset over Nigeria. First, are the S2S models skillful in predicting the monsoon onset dates and its inter-annual variability?; second, what are the contributions of global climate driving mechanisms in modulating the monsoon onset variability?; third, are the S2S models skillful in producing these global climate driving mechanisms; and, finally, are the improvements in multi-model ensembles merely due to size or does the added information by the ensembles improve model performance. Therefore, it is the aim of this study to address these concerns in other to elevate the trust levels of the S2S models at predicting an event as sensitive as monsoon onset.
The general concept of this paper is introduced and discussed in this section. The configuration of the S2S models from three model centers including the methods of analysis and data sets used is described in the next section. In the subsequent section, results from various analyses are presented. The conclusion and the implications of our findings are discussed in the final section.

Brief descriptions of the S2S models used in this study
This study uses three S2S models, based on different model configurations, from three different global prediction centers. These global prediction centers are the China Meteorological Administration (CMA), the ECMWF, and the UK MET-OFFICE. The earliest mentioned center utilizes the Beijing Climate Center Climate Prediction System (BCC-CPS-S2Sv1) version 1. The configuration is based on lagged average forecasting (LAF) method using a fully coupled BCC Climate System Model BCC-CSM1.2. The S2S Forecasts are running every day since 1 January 1994 and end with a 60-day integration. Each forecast consists of 4 LAF ensemble members, which are initialized at 00 UTC of the first forecast day (Wu et al. 2010;Li et al. 2017). ECMWF utilizes the integrated forecasting System (IFS) version 41r1. The ECMWF-S2S ensemble hind-casts Variable Resolution Ensemble Prediction System (VarEPS) is based on IFS version 41r1. It runs on an octahedral grid with 51 member ensemble (Buizza et al. 2008;Vitart et al. 2012). Operationally, the system is composed typically of coupled land, ocean, and atmosphere components. The system provides daily ensemble forecasts of a wide variety of atmospheric variables (e.g., precipitation, 2 m temperature, sea surface temperature (SST), horizontal components of wind flow at 700 hPa level, mean sea-level pressure (MSLP), dew point temperatures, etc.) with daily and sub-daily temporal resolution of the order of 6 h. More details on the ECMWF-S2S are available at http://s2spr edict ion.net/. S2S forecasting system from the UK MET-OFFICE prediction center made use of the GloSea4 modeling technique. The technique is an ensemble prediction system that uses the HadGEM3coupled GCM to model interactions across all physical components of the climate system: ocean, atmosphere, land    (Olaniyan et al. 2018) surface, and sea ice (Arribas et al. 2011). Summaries of the above descriptions are presented in Table 1.

Data set and methods
Three data sets from 1998 to 2016 are analyzed and compared in this study ( Table 2). The first data sets are the observed daily rainfall data. They are retrieved from the archive of NiMet. The data sets are from the rain gauges of 46 networks of stations spread across Nigeria (Fig. 1a). These data sets are used to determine the onset date of the monsoon. Here, the onset date of the monsoon follows those of (Stern et al. 1981) and (Omotosho et al. 2000). Here, the number of days of dry spell is slightly modified from 7 to 10 days for compliance with operational practice at NiMet. Therefore, in this study, we defined conditions for determining the onset date of monsoon as follows: (1) rainfall accumulation in the first 5 days must be at least 25 mm; (2) the first day and at least 2 other days must be wet, i.e., greater than 0.3 mm; and, (3) There should not be 10 days in the next 30 days without rainfall, i.e., no 10 consecutive days with rainfall < 0.3 mm. The second and third sets of data are, respectively, the reanalysis and S2S reforecast data sets (Table 2). While the reanalysis is from the ECMWF's Era-Interim, the three global S2S reforecast data sets are from the CMA-S2S ensemble hind-casts produced by the BCC-CPSv1, the ECMWF-S2S ensemble hind-casts produced by the VarEPS, and the UKMO-S2S ensemble hind-casts produced by the GloSea4. They are all retrieved from the ECMWF-S2S database, as supported by the World Meteorological Organization (WMO) through the World Weather Research Program (WWRP) and World Climate Research Program (WCRP). The variables from these data sets are the rainfall, SST (Newman et al. 2003), the Inter-Tropical Discontinuity (ITD: represented by 15 °C dew point temperature line [Olaniyan et al. 2015[Olaniyan et al. , 2018Lawal et al. 2016]), zonal (u) and meridional (v) winds at 700 and 200 hPa levels to, respectively, track the African Easterly Jet (AEJ: Grist and Nicholson 2001;Afiesimama 2007) and the Tropical Easterly Jet (TEJ: Nicholson et al. 2007;Nicholson and Grist 2003;Thorncroft and Blackburn 1999;Druyan and Hall 1996), and thermal low (represented by the coincidence of the hottest 2 m temperature and the lowest MSLP: Lavaysse et al. 2006Lavaysse et al. , 2009Lavaysse et al. , 2010. Here, apart from the temporal and spatial structure of both AEJ and TEJ, the evaluation will also focus on the months when both jets coincide. Studies (Nicholson et al. 2007;Omotosho and Abiodun 2007;Besson and Lemaître, 2014;Rai and Dimri 2017) have shown that enhances upper level divergence, as a result of the presence of TEJ is crucial in driving the monsoon, which, in turn, determines wet and dry years. The gradient of the thermal low is also analyzed. The gradients are determined between the Sahara (average between latitude 18-20°N and longitude 0-20°E) and the Coast (latitude 4-6°N and longitude 0-20°E). This is done by subtracting the heat-low components of the Coast from that of the Sahara. All S2S reforecast data sets are initialized on the last Thursdays of every month regardless of whether the hind-cast is "on the fly" (i.e., hind-cast produced as the same date as forecast; e.g., CMA and ECMWF) or "fixed" (i.e., hind-cast produced at a fixed date; e.g., UKMO). Retrieved reforecast variables are rainfall, SST, the ITD, zonal (u), and meridional (v) winds at 700 and 200 hPa levels, 2 m temperature, and MSLP. Here, 3, 10, and 2 ensemble members from the hind-cast of CMA-S2S (except for 2015 and 2016), ECMWF-S2S, and UKMO-S2S models (except for 2016) are used, respectively. Rainfall forecast data sets produced by the three S2S models are used to evaluate the skills of the models in determining onset dates and variability; while the other data sets, such as SST, 2 m temperature, MSLP, and u wind components at 700 and 200 hPa levels, are used to evaluate the performance of the models in producing the global climate drivers that are modulating the monsoon onset variability. The Era-interim reanalysis data sets are used to determine both the intra-seasonal and interannual teleconnections between the global climate driving mechanism and the variability of the monsoon onset. In this study, the teleconnection between the global climate drivers and the onset is based on the intra-seasonal anomaly of the global climate drivers and the inter-annual anomaly of the onset. Here, only months of highest synchronization with the onset anomaly are considered and be discussed over each region.
Statistical measures, in the form of Taylor diagrams (Taylor 2001), are utilized to determine the intra-seasonal and the inter-annual variability of standardized anomalies of all parameters. Quantitatively, correlation coefficients (r), and the normalized standard deviation (NSD) between all the S2S models (their ensemble members and the ensemble mean) with reference to the observations (station rain gauge and era-interim) are also determined. We also applied synchronization (Misra 1991) technique to quantitatively determine how well the S2S models capture the timing and anomaly signals. This only applies to the inter-annual variability of the onset dates. As explained by (Lawal et al. 2016) and (Olaniyan et al. 2018), synchronization is the percentage of periods in which the observed signs of anomalies coincide with those of simulations and or reforecasts. Furthermore, some measures of statistical significance, such as p value (Mason 2008), are performed for the synchronization and correlation skills that were determined in this study. We estimated the level of significance, where p value = 0.05, for a two-tailed experiment to decide whatever the linear association that may exist between the correlated parameters are plausible. Nigeria is divided into three climatological zones based on common climatology (Fig. 1b). Sea surfaces are also classified to depict Atlantic and Nino indices (Fig. 1b, c;(Olaniyan et al. 2018)). Finally, we regridded all data sets used in this study to 1° × 1° for uniformity and ease of analyses.

Inter-tropical discontinuity (ITD)
The three models, CMA-S2S, ECMWF-S2S, and the UKMO-S2S, are able to reproduce the intra-seasonal variability of the evolution of the Inter-tropical discontinuity (ITD) position as in Fig. 2a. This is more evident in Fig. 2b-d. Here, all the ensemble members and the ensemble means are highly correlated with the observation (Era-Interim). However, the ECMWF-S2S have the highest correlation compared to UKMO-S2S and CMA-S2S models. The strong relationship in the variability is also shown in the NSD, although the ECMWF is shown to be closest to the observation with NSD ≈ 1; in general, the correlations between the observations, the ensemble mean, and the ensemble members in all the three models are statistically significant (evaluated p values ≪ 0.05).

Temperature and pressure gradients
The S2S models generally reproduce the intra-seasonal variability of both the temperature gradient and the pressure gradient, as shown in Figs. 3a and 4a over the study areas. For the temperature gradient, the evidence of the good variability with the observed is shown in Fig. 3b-d as the NSD of the three models range from about 0.8 (CMA) to almost 1.0 (ECMWF and UKMO). While both the CMA and UKMO are seen to show little bias in both the dry (November-March) and the wet (April-October) periods, Fig. 4 a Intra-seasonal variability of the Pressure Gradient (mb) over the sahara and the coastal area of Nigeria by the CMA ensemble mean (green lines), the ECMWF ensemble mean (red line) and observed Era-Interim (black line); b, c Taylor diagrams showing the normalized standard deviations and the correlation coefficients of CMA (green) and ECMWF (red) S2S ensemble simulations with observation respectfully (triangle-ensemble members; circle for negative correlations), ensemble mean-star, and Era-Interim (observation)-black semi-circle the ECMWF seems to be colder than the observed in both periods. This implies that both the CMA and UKMO have stronger temperature gradient than the ECMWF, especially during the peak of the wet periods. In addition, all ensemble members and their means showed a strong correlation with the observed. The three models also have more than 0.95 correlation with the observation, with the ECMWF having the highest value of 0.98. Furthermore, the significant test carried out on the correlation showed that the correlation from all the model ensemble members and ensemble mean is significant as the p values are ≪ 0.005.
In case of the pressure gradient, the evidence of the good variability with the observed is also shown in Fig. 4b-d where CMA and ECMWF models have NSD of about 1.2. In addition, all ensemble members and ensemble means of both models showed a strong correlation with the observed. The two models also have more than 0.95 correlation with the observation. Significant test carried out on the correlation also indicated that the correlation from both models is significant at p values ≪ 0.005.

Sea surface temperature (SST)
To a very large extent, all the three S2S models in this study are able to reproduce the intra-seasonal variability of the standardized SST anomalies, as shown in Fig. 5a. All ENSO phases over the central and eastern Pacific are also reproduced by all the models (Fig. 5b-d). Although some rare variability spikes can be noticed, the models (the ensemble members and the ensemble means) do not, however, deviate much from the observation. This is because the NSDs of models are approximately 1.0, as shown in Figs. 5b-d and 6b-d. There are also strong relationships between the observation and the three S2S models, though with varying skills over different parts of oceans. For instance, over the Central Pacific, the average spatial correlation of the ensemble members and the ensemble mean of the ECMWF and UKMO models are well above 0.95, while the CMA is just above 0.90 (Fig. 5b-d).
The spatial correlation of the ensemble members and the ensemble mean slightly reduced over the Eastern Pacific with a correlation value of about 0.8 for CMA model, 0.90 for the ECMWF model, and about 0.95 for the UKMO model ( Fig. 6b-d). The correlations over both oceans are While S2S models from CMA and the UKMO failed to adequately capture the NSD, the ECMWF-S2S model is able to reproduce the observed intra-seasonal variability of the standardized SST anomalies over the GOG. The poor representation of the variability of the CMA and UKMO models is depicted in Fig. 7b-d. Here, the NSDs of the two models are not up to 0.30 as compared to the ECMWF of about NSD of 1.0. The relationship between the observation, the ensemble members, and the ensemble mean of the CMA and UKMO models are also very weak. As shown in Fig. 7c-d, both models average spatial correlation range between 0.25 and 0.35. However, the ECMWF model displays a strong spatial correlation of 0.8 approximately, between the ensemble members, the ensemble, and observation. Here, also the tested p value is less than 0.05, showing that all correlations are statistically significant.

African easterly jets (AEJ) and tropical easterly jet (TEJ)
The reproducibility of the intra-seasonal variability of the observed upper air jets by the 3 S2S-models is not consistent. Each of the model ensemble members and the ensemble means, as shown in Fig. 8b-d, displays different unique characteristics. The inconsistencies in the models are shown in the strength, the spatial, and the temporal positions of the two jets. For instance in Fig. 8c, the ECMWF model is able to reproduce the strength of the AEJ as compared to the observation excluding 2010. Throughout the years considered, the probability that at least 10% of the ECMWF ensemble members will reproduce both the spatial and the temporal position of the AEJ is strong. However, also throughout the years, at least two of the ensemble members will produce AEJ in the month of May and June. The spatial and temporal position of the TEJ is also well reproduced by the ECMWF model. Although the ECMWF model may under-estimate the strength of the TEJ, especially in those years when it is not observed, it, however, reproduced all the years when the AEJ and TEJ coincide.
In general, as in Fig. 8b, the CMA model (including all the ensembles) poorly reproduced both the AEJ and TEJ in strength and in position. For instance, it does not adequately reproduce the position and the strength of the TEJ. However, for any years (1998, 2003, 2008, 2012,  The UKMO model (ensemble) also showed poor reproducibility of the Jets. It shows a very weak skill in producing the TEJ throughout the research periods excluding 1999 and 2003. Although it is able to reproduce the AEJ, however, it also inadequately reproduced the temporal and spatial positions of the Jet, almost through the period of years considered.

The monsoon onset: anomaly and global drivers teleconnection
The monsoon onset and anomaly As shown in Figs. 9a, 10a, and 11a), all the three models used reproduced the northward migration of the onset dates from the Coast to the Sahel. In terms of the inter-annual variability of the onset dates, each of the model displays distinct biases over each region (figure not shown). For instance, while the UKMO displays both negative and positive biases over all the regions, the biases of the CMA and EMWF models progressively improve from the Coast to the Sahel. In addition, the inter-annual variability of the observed onset dates as reproduced by the three models varies across the region. For example, as shown in Figs. 9b, 10b, and 11b, three ensemble members of the CMA model over the Guinea area, one member over the Savannah and all the ensemble members over the Sahel, showed relatively high variability when compared to the observation. The variability of the ECMWF ensemble members is generally normal as most members cluster around NSD value of 0.8 and 1.3 in the entire regions. The inter-annually variability as reproduced by the UKMO model is also normal as the average value of the ensemble mean ranges between 0.75 and 1.2.
The skill of all the models used, based on direct relationship with the observed onset dates, over the Gulf of Guinea and the Sahel is poor and statistically not significant. Over these areas, the correlation of all the ensemble members and ensemble means, as shown in Figs. 9c, d and 11c, d, is not more than 0.4. While the ECMWF have few members showing the negative correlation in all the regions, as shown in Fig. 11b, d, the negative correlation showed by both the UKMO and the CMA is only in the Sahel. Meanwhile, as inferred from Fig. 10b-d, all the three models have considerable improvement in skill over the Savannah with the UKMO (ensembles) models having the strongest skill. As seen over the same area, the relationship between the UKMO model and the observation has an average correlation of almost 0.65. This correlation is statistically significant as the derived p value is ≪ 0.05. Again, the relationship between the CMA ensemble members and the mean is weak as the correlation of most ensemble members is below 0.4. The derived p value is ≫ 0.05, indicating that it has no statistical significance. The least Skilful over the Savannah is the ECMWF with all the ten ensemble members. Here, apart from having two of the members negatively correlated, almost all the ensemble members are less than 0.4 in correlation.
As can be seen in Figs. 12, 13, 14a, c, e, the inter-annual variability of the monsoon onset anomaly as produced by the ensemble members including the ensemble mean of the three models considered, spread over the observed, with varying exceptions. For instance, over the Gulf of Guinea in 2000 and 2016, Fig. 12a showed that the ensemble spread from the CMA model is well outside the observation. This is also true for the ECMWF ensemble spread in 2016, as shown in Fig. 12c. Over the same area of the Gulf of Guinea in 1999, 2007, 2008, 2009 shows that the UKMO model ensemble spread is outside the observation. The exceptions over Sudan are in 2002, 2004, 2005by CMA model, while the exceptions are in 2003, 2004, 2008, and 2012 by the UKMO model. The correlation of all the models in simulating the inter-annual variability of onset anomaly, as shown in Figs. 12, 13, 14b, d, f, is similar in skills with how they simulate the inter-annual variability of onset dates. This is because all the three models exhibit the worst correlation skill over both the Gulf of Guinea and the Sahel with improved skill over the Savannah. For instance, in all the regions, the correlation skill of the CMA ensemble members does not exceed 0.45, although, with negative correlation on all the regions, the correlation skill of the ECMWF ensemble members ranges from 0.1 to 0.5. However, despite the poor skill over the Gulf of Guinea and the Sahel, the ensemble mean and one of the ensemble members show relatively strong correlation skill of more than 0.65 over the Savannah.
While most ensemble members from the three models may have low correlation skills, in terms of a direct relationship with the observation, the synchronization skill based on the sign consistencies, as shown in Fig. 15a-c, suggests different disparities over the regions. For instance, the synchronization skills of both the CMA and the ECMWF generally increase from the Sahel to the Gulf of Guinea. The synchronization of the ensemble mean from CMA ranges from 45 to 70% over (Sahel and Gulf of Guinea). However, instead of the sharp rise in the CMA, the ECMWF ensemble mean gradually increases with synchronization range from 65 to 70%. Meanwhile, the synchronization skill of the UKMO is generally strongest over the Savannah (65%) and weakest at the Gulf of Guinea (45%). In addition, there are three different patterns in the synchronization by the ECMWF ensemble members, as shown in Fig. 15b. The first pattern is shown by ensemble members (2, 6, and 10). Here, all the members showed the synchronization over the Sahel (52, 52, and 58%, respectively) to be the weakest and over the Savannah (80, 75, and 78%) to be the strongest. The second pattern as shown by ensemble members (3, 7, 8, and 9) revealed a general weak synchronization over the Savannah as compared to the Gulf of Guinea and Sahel. The third pattern between ensembles 1 and 4 is similar but reversed. For instance, as shown in Fig. 15a, while ensemble-member-1 showed synchronization skill increases toward the Sahel, ensemblemember-4 showed the synchronization skill decreasing to the Sahel. Ensemble members (1 and 2) of the CMA are also similar in pattern with ensemble members (1 and 4) of the ECMWF. However, here, ensemble-member-1 showed that the synchronization skill increases toward the Gulf of Guinea, while ensemble-member-2 showed the synchronization skill decreasing to the Gulf of Guinea. In addition, ensemble-member-3 of the CMA is similar in pattern as ECMWF ensemble-member-8 where synchronization skill over the Savannah as compared to the Gulf of Guinea and Sahel is showed to be the weakest. There is no definite pattern displayed by the UKMO ensemble members; however, Fig. 15c shows that, out of all the three models used, the ensemble-member-2 of the UKMO model has the strongest synchronization skill of more than 85%, but over the Savannah. It, therefore, implies that, over the Gulf of Guinea, for instance, both CMA and ECMWF could, on average, make seven correct forecasts out of ten, based on the onset anomaly. While the ECMWF could make the same correct forecasts over the Sahel, both the CMA and UKMO could only make four correct forecasts out of ten. Meanwhile, over the Savannah, both ECMWF and the UKMO could still approximately make seven correct forecasts out of ten. The monsoon onset anomaly and global drivers teleconnection The teleconnection between the inter-annual variability of the anomaly of the global drivers, with the interannual variability of the onset anomaly (i.e., early or late onset) over the Gulf of Guinea, Savannah, and the Sahel, as shown in Figs. 16,17,18 and 19 is also presented in Table 3. The results from the Table showed different unique characteristics of the anomaly of each of the global drivers in different months, with the variability of the onset anomaly over each region.
For instance, over the Gulf of Guinea, although, with different correlation and synchronization skills, the table shows three direct and a "1 year lag" link (teleconnection) between the standardized anomaly of global drivers and the onset standardized anomaly. Over this area, the January SST anomaly over the Central Pacific, the Gulf of Guinea, and the Pressure Gradient anomaly are shown to be directly linked with the anomaly of the onset dates. Although, there are disparities in the correlation and synchronization skill of each global driver; however, the anomaly of the ITD position in August is also directly linked with the onset anomaly but with a year lag. This implies that the late onset (positive anomaly of the onset) in one particular year is linked with above normal in the position of the ITD in August of the previous year. Similarly, the SST anomaly over the Central Pacific in February and the Gulf of Guinea in May is directly linked with the onset anomaly over the Sahel. In addition, linked but with a year lag is the anomaly in the position of ITD in March. The teleconnection of each driver with the onset over this zone has different skills in synchronization and correlation. Here, while other drivers showed a positive correlation, the link between the Pressure Gradient anomaly and the rainfall onset anomaly is negatively correlated. It is the only driver over this region, with a "reversed 1 year lag" teleconnection. This implies that late onset (positive anomaly of the onset) in one particular year is linked with below normal in the magnitude of the Pressure Gradient the previous year's March. Throughout the years under study, the result reveals no link between the anomaly of the Temperature Gradient and the onset anomaly except over the Savannah. Here, the link between anomaly of the Temperature Gradient and the onset anomaly is direct but also lag by a year, with about 78% synchronization and correlation of 0.5. In addition, linked directly is the SST anomaly over the Central Pacific in April. Furthermore, Fig. 8a shows that both the AEJ and TEJ coincided 11 times, during the period under study. The coincidence of these jets in a particular year occurs with the onset anomaly over different regions with varying degrees but strongest over the Gulf of Guinea. For instance, the late onset over the Gulf of Guinea in a particular year coincides by about 70% with the coincidence of the two jets the previous year.
The existing year lag in the teleconnection between some of the global drivers especially the ITD, and the anomaly of the onset, may not be unconnected with the two-way interactions between soil moisture and the atmosphere. Several studies have documented the modulation of surface water and heat fluxes with the soil moisture variability in form of feedback mechanism Suarez 1996, 2001;Douville et al. 2001;Berg et al. 2017). For instance, according to Koster and Suarez (2001), a period of a particular rainfall high has the tendency of producing a positive anomaly in soil moisture, and the dissipation of that anomaly through evaporation and other processes may take weeks to months. In effect, soil has retention capability; this simply implies that the soil can ''remember'' periods of dry and wet events that caused an anomaly long after these conditions have been forgotten by the atmosphere. Overall, all correlations below 0.5 are statistically not significant having the evaluated p value ≫ 0.05.

Conclusion
This paper assesses the skills of three sub-seasonal to seasonal (S2S) models in predicting monsoon onset and its variability. It also examines the global drivers modulating the variability and their teleconnections with rainfall onset anomaly. The results showed that each of the models used exhibits unique and different characteristics over each classified region in Nigeria. For instance, all three models are able to simulate the Northwards migration of the onset dates adequately with inherent biases. While the biases of both the CMA and the ECMWF models improve progressively towards the Sahel, the bias of the UKMO model over the Gulf of Guinea is considerably small (± 10 days). In case of the onset anomaly, results showed that despite the poor performance of the models over the Gulf of Guinea and the Sahel, there is a considerable improvement in the correlation skill of the models over the Savannah. Overall, as the CMA Fig. 11 Inter-annual variability of rainfall onset dates (Julian dates) over a the Sahel Area of Nigeria from rain gauge observation(black line) the CMA ensemble members and ensemble mean (solid and broken green lines), the ECMWF ensemble members and ensemble mean (solid and broken red lines) and the UK-METOFFICE ensemble members and ensemble mean (solid and broken blue lines); b-d Taylor diagrams showing the normalized standard deviations and the correlation coefficients of CMA (green), ECMWF (red) and UK-METOFFICE (blue) S2S ensemble simulations with observation respectfully (triangle-ensemble members; circle for negative correlations), ensemble mean-star, and observation-black semi-circle model is shown to have the least skill in all the regions, the UKMO model is shown to have the best correlation skill especially over Savannah. The probability of making a correct forecast based on onset anomaly, as suggested by the synchronization skill and as displayed by different model ensemble-member, differs over each region. However, on average, the ECMWF model could make seven correct forecasts out of ten overall regions; the CMA, four and seven Results further showed that all the three models reproduced the evolution and variability of the global drivers modulating the monsoon onset. However, only the ECMWF Fig. 13 Inter-annual variability of rainfall onset dates standardized anomaly over the Savannah Area of Nigeria from rain gauge observation (blue bar) with a the CMA ensemble members and ensemble mean c the ECMWF ensemble members and ensemble mean and e the UK-METOFFICE ensemble members and ensemble mean. b, d, f Taylor diagrams showing the normalized standard deviations and the correlation coefficients of CMA (green), ECMWF (red), and UK-METOFFICE (blue) S2S ensemble simulations with observation respectfully (triangle-ensemble members; circle for negative correlations), ensemble mean-star, and observation-black semi-circle model was able to reproduce the spatio-temporal characteristics of both the AEJ and the TEJ. Again, the understanding of the teleconnection between most global drivers and the anomaly of the onset may require further investigation. For instance, results show a year lag teleconnection between the anomaly of ITD position in August and March and onset anomaly over the Gulf of Guinea and the Sahel respectfully. As explained in "The monsoon onset anomaly and global drivers teleconnection", the lag may be as a result of the retention capacity of soil moisture, which studies have shown to have a feedback effect on rainfall variability. Therefore, further research may still be needed. This is especially to understand how the soil moisture feedback mechanisms influence and modulate onset anomaly in Nigeria. However, despite the lagging displayed by some drivers, teleconnection of the SST over the Central Pacific is direct. It also has good correlation and synchronization skills, especially with the onset anomaly over both the Gulf of Guinea and the Sahel. The results showed that the CMA model has the least skills. However, all the S2S models, despite the inherent biases are able to predict rainfall onset over Nigeria, within the sub-seasonal timescale. Finally, the results show that improvements in multi-model ensembles could be a value-added information, able to significantly improve model performance. Fig. 15 Synchronization (%) of the simulated inter-annual rainfall onset dates anomalies between a the CMA Ensemble members and the Ensemble mean (b), the ECMWF Ensemble members and the Ensemble mean, and (c) the UK-METOF-FICE Ensemble members and the Ensemble mean. All models are with reference to observation over the climatological zones of Nigeria (GoG, Savannah, and Sahel) Fig. 16 Inter-annual variability of standardized SST anomaly over the central pacific in the month January (red), February (green), and April (yellow) with the rainfall onset standardized anomaly over the Gulf of guinea (a), Savannah (b), and Sahel (c) Areas of Nigeria  Fig. 17 Inter-annual variability of standardized SST anomaly over the Gulf of Guinea in the month January (red) and May (green) with the rainfall onset standardized anomaly over the Gulf of guinea (a), Savannah (b), and Sahel (c) Areas of Nigeria Fig. 18 Inter-annual variability of the ITD position in the month March (red), and August (green) with the rainfall onset standardized anomaly over the Gulf of guinea (a), Savannah (b), and Sahel (c) Areas of Nigeria  Fig. 19 Inter-annual variability of the pressure gradient standardized in the month January (red box) and March (green box); and temperature gradient in January (yellow box) with the rainfall onset standardized anomaly over the Gulf of guinea (a), Savannah (b), and Sahel (c) Areas of Nigeria Table 3 Summary of the teleconnection between the inter-annual variability of the anomaly of the global drivers with the inter-annual variability on the onset anomaly over the Gulf of Guinea, Savannah, and the Sahel NA not applicable