Revisiting interannual to decadal teleconnections influencing seasonal rainfall in the Greater Horn of Africa during the 20th century

This paper revisits teleconnections for three major rainy seasons at the Greater Horn of Africa for the period 1901–2013. Sea surface temperature‐based climate indices known to influence Short, Kiremt, and Long Rains are used in a comprehensive statistical analysis to detect non‐stationary behaviour in teleconnections and to split them into interannual and decadal time scales. Physical mechanisms are proposed for identified significant non‐stationary teleconnections. Interannual variability in the October–December Short Rains is predominantly influenced by the Indian Ocean Dipole (IOD) with a percent variance explained (PVE) of up to 80% in recent years. However, abrupt shifts in this teleconnection occurred around 1918, 1951, 1987, and recently. The Short Rains also correlate strongly with El Niño‐Southern Oscillation (ENSO). However, ENSO influence on Short Rains is mediated by in‐phase occurrence of IOD. Decadal variations in Short Rains are more directly explained by low‐frequency variability in the Pacific Ocean (PO). A hitherto undocumented non‐stationary relation was found between Atlantic Niño 3 and June–September Kiremt Rains. The non‐stationarity seems to be related to a decadal regime shift in the West African monsoon system in the late 1960s. The variability of Kiremt Rains is also strongly associated with ENSO, although the recent increased correlation is not non‐stationary. Consistent with recent studies, the post‐1998 March–May Long Rains decline is strongly associated with decadal variability in the PO. The PVEs in the stable correlations using Pacific Decadal Oscillation and Interdecadal Pacific Oscillation indices range from 25 to 64%, mostly due to low‐frequency variability (>8 years). The results have ramifications for seasonal‐to‐decadal forecasting and can be taken to design modelling studies to gain insights into the physical mechanisms.

warning systems network and world food programme, 2016). These extreme events and rainfall variability in GHA are dominated by changes on the large scale, with a clear link to tropical ocean-atmosphere variability. Given the uncertainty of climate projections in the region, the inevitability of drought and associated socio-economic distress, understanding the mechanisms that produce this variability and developing both dynamical and statistical approaches for extended-range forecasts and climate information is urgently needed.
Besides orographic effects and other factors, the seasonal rainfall climatology over GHA is highly influenced by the north-south displacement of the tropical rainbelt, often, but inadequately termed inter tropical convergence zone (Nicholson, 2018). Much of the northern Ethiopian Highlands and northwestern parts of GHA have boreal summer monsoon regimes with maximum rainfall during the period of June-September (JJAS) known locally as "Kiremt Rains". This season accounts for some 50 to 80% of the rainfall over Ethiopian agricultural regions (Korecha and Barnston, 2007). The equatorial part of GHA has two rainy seasons; the bimodal patterns are experienced over many parts of GHA during the "Long Rains" of March-May (MAM) and the "Short Rains" of October-December (OND).
Although the Long Rains contribute a larger fraction to the annual total, the Short Rains show a larger degree of interannual variability relative to the Long Rains (Hastenrath et al., 1993;Black et al., 2003) and have a larger impact on the society through changes in the regional hydrological cycle (Behera et al., 2005). Several studies have investigated the relationship between the Short Rains and the El Niño-Southern Oscillation (ENSO; Ropelewski and Halpert, 1987;Ogallo, 1988;Hastenrath et al., 1993;Mutai et al., 1998;Nicholson, 2017). These statistical studies find a strong connection between ENSO and East African rainfall. Contrarily, atmospheric general circulation models suggest that the Indian Ocean (IO) sea surface temperatures (SSTs) exert a greater influence on the East African Short Rains than the Pacific SSTs (Goddard and Graham, 1999;Latif et al., 1999). This has resulted in considerable controversy concerning the primary driver of the Short Rains interannual variability. Consequently, other studies tried to compare the fractional influence of IO Dipole (IOD) and ENSO on Short Rains (Saji et al., 1999;Clark et al., 2003) and suggested that IOD is a major contributor to the variability of the Short Rains. In view of these finding, the focus has been shifted to the IO for understanding the variability of the Short Rains (Black et al., 2003;Yamagata et al., 2004;Behera et al., 2005;Ummenhofer et al., 2009;Bahaga et al., 2015;Nicholson, 2015;Endris et al., 2016). The conclusion is that rainfall over East Africa (Indonesia) is increased (decreased) during a positive IOD event and reverse ocean-atmosphere conditions happen during negative IOD events. The physical mechanism for the IOD influence on Short Rains is explained by a Gill-type response (Gill, 1980) to the western warm pool of the IO and the modulation of the IO Walker circulation (Webster et al., 1999;Hastenrath et al., 2011).
Focusing on the Ethiopian Kiremt Rains, some studies identified large-scale forcings important for seasonal prediction (Camberlin, 1997;Gissila et al., 2004;Korecha and Barnston, 2007;Diro et al., 2008;Segele et al., 2009a;Diro et al., 2011aDiro et al., , 2011bNicholson, 2014;Degefu et al., 2016). These studies have shown a significant simultaneous association between the summer rainfall over the Ethiopian Highlands and ENSO indices. El Niño (La Niña) episodes were shown to be associated with less (more) summer rainfall across much of the country. Other studies found an association between Ethiopian JJAS rainfall and SST anomalies over the southern Atlantic Ocean (AO), Gulf of Guinea, and IO (Segele et al., 2009b;Diro et al., 2011a;Williams et al., 2012) also suggesting that the transport of moist air from the Congo Basin plays a role in the interannual variability of boreal summer rainfall in eastern Africa. Although southern AO, Gulf of Guinea, and Congo air masses have long been assumed to be a source of moisture for summer rainfall over the Ethiopian Highlands, the role of the south AO on the interannual variability is still disputed. Segele et al. (2009b) and Segele et al. (2015) suggest that the south AO has an influence on the interannual variability of Kiremt Rains. In addition, Lyon (2014) discusses the climate characteristics of the Ethiopian severe drought of 1982-1984 during the June-July-August season and point to a role of AO forcing. On the other hand, there have only been a few studies conducted on decadal summer rainfall variability and its relationship to global SSTs. Jury (2010) found a modest association between areas of summer rainfall over Ethiopia and SST variability associated with the Atlantic Multidecadal Oscillation (AMO), while variability rainfall over the bimodal southern and eastern zone of the Ethiopian Highlands were associated with the cool phase of the Pacific Decadal Oscillation (PDO).
Unlike the Short Rains, understanding interannual variability and predictability of the Long Rains has been a challenge for researchers and forecasters, since no strong oceanic or atmospheric forcing has been identified (Hastenrath et al., 1993;Camberlin and Philippon, 2002;Liebmann et al., 2014;Nicholson, 2017). Moreover, the variability of MAM rainfall reflects the net influence of factors acting on intraseasonal, interannual, decadal, and multidecadal timescales (e.g., Pohl and Camberlin, 2011;Omondi et al., 2013;Lyon et al., 2014;Nicholson, 2017;Vigaud et al., 2017;Ummenhofer et al., 2018). Motivated by increased frequency of drought in the last decade, some studies emphasized the recent downward trend of the Long Rains over parts of GHA (Funk et al., 2008;Williams and Funk, 2011). Lyon and DeWitt (2012) further examined these persistent droughts and suggested that the decline started abruptly around 1999. They also demonstrate that the decline is part of a large-scale precipitation anomaly pattern that extends across the IO and Pacific Ocean (PO). Further studies provide evidence that the decline in GHA Long Rains is associated with a shift of warmer SSTs over the western tropical Pacific and cooler SSTs over the central and eastern tropical Pacific (Lyon and DeWitt, 2012;Lyon, 2014;Lyon et al., 2014;Yang et al., 2014;Vigaud et al., 2017;Ummenhofer et al., 2018) and the PO multidecadal variability (Lyon, 2014;Yang et al., 2014). Similarly, Vigaud et al. (2017) suggest that a more frequent dry regime in May since 1998/1999, associated with an earlier onset of the Indian summer monsoon and Somali jet, could account for a recent abrupt shift observed in GHA Long Rains. In contrast, other studies link the drying trend to anthropogenic and aerosol forced Indo-PO warming, resulting in the westward extension of the Indo-Pacific warm pool, causing a westward shift of the Walker circulation and a subsidence anomaly and drying over East Africa (Williams and Funk, 2011;Funk et al., 2014;Liebmann et al., 2014;Funk and Hoell, 2015;Hoell et al., 2017). Alternatively, Hoell et al. (2017) suggest that human-induced changes in tropical SST exacerbated the effect of natural Pacific decadal variability, enhancing GHA drying during the Long Rain season. Additional modelling studies examined several hypotheses related to the decline of the Long Rains and concluded that anthropogenic aerosol emissions played a role (Rowell et al., 2015;Tierney et al., 2015).
Despite the above-referenced numerous studies on teleconnections between GHA seasonal rains and remote climate anomalies, apparent disagreements among them regarding attribution arise from consideration of different seasons, datasets, and varying analysis periods due to nonstationarity of climate drivers (Nicholson, 2015(Nicholson, , 2017. Especially, few studies exist on the non-stationarity and multidecadal changes of teleconnections. These studies only deal with changes of Short Rains teleconnections over time (Camberlin and Philippon, 2001;Clark et al., 2003;Nakamura et al., 2009;Manatsa et al., 2012;Manatsa and Behera, 2013;Nicholson, 2015Nicholson, , 2017. Furthermore, there is, to the best of our knowledge, no study that assigns the GHA seasonal teleconnections to interannual or/and decadal timescales. Several such studies exists for the West African monsoon (Diatta and Fink, 2014;Suárez-Moreno et al., 2018). To fill these gaps, the present study provides a comprehensive statistical analysis considering all three major GHA rainy seasons. The specific objectives of the present study are: • To statistically test the non-stationarity of previously identified teleconnections for 11-and 31-year windows during the whole investigation period 1901-2013. • Here for the first time, to consider SST indices in all ocean basins and examine the evolution of a statistical relationship over both varying time window lengths for fixed starting years, as well as over changing starting years for fixed time window length in a triangle representation (Diatta and Fink, 2014), indicating potential non-stationary behavior of GHA rainfall teleconnections in all possible combinations, • To propose a physical mechanism for identified significant non-stationary teleconnections. • To identify and assign the teleconnections to periodicities below or/and above 8 years. • To revisit and attribute the recent increased frequency of Long Rains droughts in GHA as part of the natural decadal variability driven by Pacific multidecadal variability.
In Section 2, datasets and analysis methods are described. A brief discussion on interannual and decadal variability during the three major rainy seasons over GHA is given in Section 3. Section 4 presents multidecadal changes in the teleconnections with some of the physical mechanisms for non-stationary correlation also explored. Finally, the results are summarized and discussed in Section 5.

| Rainfall indices in the GHA
In order to re-assess the interannual to decadal changes in the association between seasonal rainfall over GHA and indices of remote ocean-atmosphere variability, we have used observed precipitation and SST datasets. To investigate long-term rainfall fluctuations in the investigation period, we have used version 7 of the Global Precipitation Climatology Centre (GPCC) gridded monthly precipitation (Schneider et al., 2014). This is a gauge-based, global, land-surface dataset at 0.5 × 0.5 spatial and monthly temporal resolution for the period 1901-2013. To validate GPCC and show GHA rainfall climatology in recent times, gridded monthly rainfall estimates from the Climate Hazards group InfraRed Precipitation with a Station dataset (CHIRPS; Funk et al., 2015b) were also utilized.
In this study, the standardized precipitation index (SPI; McKee et al., 1993) was used to compute rainfall indices. SPI is very flexible, allowing calculations aggregated over different spatial scales and temporal domains (World Meteorological Organization, 2012;Lavaysse et al., 2015). This study focuses on the monthly timescale and therefore the SPI was calculated using area-averaged, monthly accumulated GPCC and CenTrends precipitation data. From the monthly SPI values, the seasonal mean SPI time series were constructed for three major rainy seasons.
The GHA region comprises 11 countries: Kenya, Ethiopia, Eritrea, Somalia, Djibouti, Uganda, Tanzania, Rwanda, South Sudan, Sudan, and Burundi and includes areas with bimodal and unimodal rainfall distributions. The GHA region is divided into three sub-regions based on seasonal rainfall climatologies (Figure 1). Our analysis particularly focuses on three well-defined rainfall regimes in the respective sub-regions, which cover most parts of GHA and are named as follows: the Equatorial East African Short Rains SPI (EEA-SR-SPI) is represented by the area-averaged OND rainfall over the land grid boxes of the region 5 o S-6 o N, 35 -50 E (blue box in Figure 1). EEA covers most of Kenya, the southern parts of Ethiopia and Somalia, and northern Tanzania. It is characterized by a bimodal rainfall distribution with the major rainfall season in MAM. EThiopian Highland Kiremt Rains SPI (ETH-KR-SPI) is defined by the areaaveraged JJAS rainfall over the land grid boxes of the region 6 -14 o N, 35 -40 E (orange box in Figure 1). ETH covers most of the northern Ethiopian Highlands and some parts of Sudan. It has a rainfall maximum during boreal summer and a pronounced dry season in boreal winter. The GHA Long Rains SPI (GHA-LR-SPI) designated by the area-averaged MAM rainfall over the land grid boxes of the region 10 o S-12 o N, 30 -50 E (red box in Figure 1). Thus, throughout the paper, the terms Long Rains, Short Rains, and Kiremt Rains will correspond to specific locations within the study domain (i.e., the red, blue, and orange boxes; Figure 1). We have verified that the results presented in this article are robust with respect to small variations in size and position of these regions. The SPI indices were computed for the entire investigation period.
In addition, the GPCC dataset was compared with the recently released Centennial Trends (CenTrends) precipitation estimates (Funk et al., 2015a). The running correlation between GPCC and CenTrend for a 31-year sliding window shows a stable correlation of nearly 1, both for EEA-SR-SPI and GHA-LR-SPI, while the correlations for ETH-KR-SPI are above 0.8 with an increasing trend after 1960 (not shown). Such minor difference might be attributed to the presumably large number of Ethiopian stations that contribute to CenTrends dataset. The sensitivity of these precipitation products is also tested for multidecadal changes in the seasonal rainfall teleconnections. In general, the results from the CenTrends are consistent with GPCC.

| SST-based indices of remote climate states
In this study, various previously identified SST-based climate teleconnection indices were utilized to explore their impacts on interannual to decadal seasonal rainfall fluctuations over GHA. In addition, these indices were also used to investigate and understand the multidecadal changes in the teleconnections. Gridded monthly SSTs were obtained from the extended reconstructed SST version 4 (ERSSTv4) of the National Oceanic and Atmospheric Administration (NOAA) National Climate Data Center (NCDC) (Huang et al., 2015). For comparison, we have also used the Hadley Centre Sea Ice and SST dataset (HadISST; Rayner et al., 2003). The SST-based indices fall into the geographical sectors shown in Figure 2. A brief summary of how indices of remote climate anomalies are defined, created, and relevant previous works suggesting a climate indices-rainfall relationship is given in Table 1. Unless otherwise noted, all indices were calculated for the entire investigation period 1901-2013, and the 1950-1979 climatological base period was used for computing anomalies. For more details on the indices and their definitions, the readers are referred to the relevant literature cited in Table 1.

| Reanalysis datasets
We use the European Centre for Medium-Range Weather Forecasts (ECMWF) 20th Century Reanalysis product (ERA-20C; Poli et al., 2016). Monthly averaged fields at 1 × 1 horizontal resolution were considered. The analysis spans 111 years, from January 1900 to December 2010. This product is used as an input for our composite analysis in Section 4.4.

| Statistical methods
To quantify the linear association between remote indices of ocean-atmosphere climate anomalies and GHA rainfall shown in later sections, several statistical methods were used. The Pearson's linear correlation and percentage of variance explained (PVE) are employed to evaluate the link between the remote climate indices and GHA seasonal rainfall. The triangle representations used in Fink et al. (2010) and Diatta and Fink (2014) were chosen to illustrate the correlation coefficient over all possible combinations of starting dates and time window lengths. Besides the simple linear correlations, a partial correlation technique (e.g., Yule, 1907) is used to examine relationships between three variables, for example, the partial correlation between dipole mode index (DMI) and Short Rains, while excluding the influence that arises owing to the correlation between N34 and Short Rains. Moreover, composite analysis techniques are performed to investigate potential physical mechanisms.
One aspect of this analysis is the use of both unfiltered and filtered data. The low-and high-pass Butterworth filters (Butterworth, 1930) with a half power period of 8 years were applied to both time series before the calculation of linear correlations. The low-pass filtered time series are thus representing fluctuations of 8 years and longer (i.e., decadal timescales), and the high-pass filtered data are reflecting interannual variability. The critical values of the significance levels for Pearson's linear correlations are calculated by an F-test (Taubenheim, 1969). We have then used both the standard degree of freedom (hereafter SDOF) based on the total number of data pairs that enter the correlations and the effective degree of freedom (hereafter EDOF) as described in Fink and Speth (1997), which takes autocorrelations into account. We have also applied the non-stationarity significance test proposed by Van Oldenborgh and Burgers (2005). This method uses Monte Carlo techniques to generate synthetic time series under constant regression coefficients and approximates the weather noise with a normal distribution. One thousand time series of climate indices are generated with the real observations replaced by synthetic observations with the same 113-year regression with the rainfall indices. Then the running correlations of the synthetic times series with observed rainfall indices are calculated using a window length of 11 and 31 years. The null hypothesis of constant correlation is tested, and if rejected, the teleconnection is considered to be non-stationary over the 113-year period.

| GHA RAINFALL INTERANNUAL AND DECADAL VARIABILITY
In this section, the 1901-2013 year-to-year and lowfrequency evolution of the three major GHA seasonal rainfall indices are examined. Figure 3a shows interannual variability and low-pass filtered SPI of the Short Rains based on GPCC data. Compared with the Long Rains (Figure 3c), there is recent evidence of strong interannual variability of Short Rains, with the years of 1961 and 1997 standing out as the wettest years ( Figure 3a). Moreover, there is a clear phase shift in the decadal variability. For instance, the dry Short Rains decadal signals during 1935-1943and 1981-1989occurred during wetter Long Rains, while during 1921-1927-2010 wetter Short Rains occurred with persistent Long Rains drought.
For ETH Kiremt Rains, Figure 3b indicates large interannual variations but also characteristics that are reminiscent of the Sahel multidecadal variability, with comparatively wet conditions during the 1940s, 1950s, and 1960s, then a downward trend with increased drought in the 1980s and 1990s and a recent recovery, albeit between 2009 and 2013 JJAS rainfall was quite low (Figure 3b). Although JJAS rainfall over ETH shows multidecadal behaviour similar to West Africa, the involved mechanisms responsible might vary substantially, since features such as the Somali low-level jet play a role (Segele et al., 2009a).
Unlike the Short Rains, interannual variability of the Long Rains shows a clear persistent decline in recent decades ( Figure 3c). Moreover, the low-frequency SPI variations suggest that the drying trend since about 1998 is part FIGURE 2 Locations of remote climate anomalies used in this study as predictors for GHA seasonal mean rainfall. The indices and corresponding geographical areas are detailed in Table 1 [Colour figure can be viewed at wileyonlinelibrary.com] of the decadal variability of climate over GHA. Furthermore, there are other periods where the time series of Long Rains shows upward and downward trends. This is consistent with recent studies that the post-1998 Long Rains decline in MAM is a manifestation of decadal variability (Lyon, 2014;Yang et al., 2014). In general, it is important to note that in marked contrast to the declining trend in the Long Rains after 1980, the Short Rains show a weak upward trend (Figure 3a,c; Rowell et al., 2015;Nicholson, 2017).

| Multidecadal changes associated with Short Rains teleconnections
In this subsection, the dominant teleconnections for the Short Rains and their multidecadal changes are analysed. Figure 4 shows the time evolution of the correlation coefficients between mean EEA-SR-SPI and four SST-based remote climate indices for a 31-year sliding window. It is clear from Figure 4 that there are non-stationary correlations to different degrees. One good example is the eastern IO index (EIO), which reflects the time variations of the eastern node of the DMI. It changed from insignificant correlations in the early 1940s to significantly negative correlations of about −0.7 in the late 1960s. Based on the non-stationarity significance test, the changes in the running correlation between the EIO index and the Short Rains are significant. A similar change in the link over time, but opposite in correlation, is that between EEA-SR-SPI and DMI. For the 31-year windows centred on years between mid-1920 and 1967, the linear correlation coefficient rose from statistically insignificant values of 0.3 to significant values of around 0.75, but for the DMI did not pass the significant test of non-stationarity. Besides, the change in the correlation with the western IO index (WIO) shows a fairly similar trend compared with DMI except that the evolution of correlations in WIO are lower than DMI after the 1950s (Figure 4).
On the other hand, the correlation coefficient between EEA-SR-SPI and N34 changes from a statistically insignificant correlation of nearly 0.2 in mid-1920s to a significant correlation of 0.65 in mid-1930s, then remains significant throughout most of the investigation period, despite an abrupt decrease below the significance level around 1976 ( Figure 4). This sudden change was presumably related to a Pacific or global climate shift (Trenberth and Hurrell, 1994;Zhang et al., 1997;Meehl et al., 2009). The recent increase in N34 and EEA-SR-SPI relation could be an interference of ENSO in the EEA-SR-SPI and IOD relation as discussed by Goddard and Graham (1999) and others, finding that the Pacific forces the changes in the IO that are linked to the interannual variability of Short Rains. Another possible      W. Short Rains (Ropelewski and Halpert, 1987) and Kiremt Rains (e.g., Korecha and Barnston, 2007) N3 Trenberth (1997) El Niño 3 as standardized SST anomalies averaged in region bounded by 5 explanation is that increased coupling between the tropical Pacific and IO during the ENSO cycle since the mid-1970s resulted in frequent co-occurring IOD and ENSO events (Annamalai et al., 2005). To disentangle this competing influence, we applied a partial correlation analysis between DMI and Short Rains anomalies after removing/partialling out the Niño-3.4 influence. Positive and significant correlations across East Africa indicate that excessive rains during OND are preferentially linked to concurrent IOD ( Figure S1a). However, a similar analysis with Niño-3.4 index and rainfall anomalies for the same season after removing/partialling out IOD indicate that the ENSO anomaly leads to a spatially non-coherent negative correlation over parts of East Africa ( Figure S1b). Corresponding modelling studies (e.g., Goddard and Graham, 1999;Bahaga et al., 2015) and observation studies (e.g., Black et al., 2003;Behera et al., 2005) revealed viable evidence of the dominant role of the tropical IO in modulating the Short Rains interannual variability (Hastenrath et al., 2011;Nicholson, 2015). Although there are still some disagreements among scientist, it has now become clearer that the association between the Short Rains in East Africa and ENSO is actually more the result of an indirect forcing by the PO on the IO (Black et al., 2003;Annamalai et al., 2005). Collectively, the centred correlations of all remote climate indices were insignificant for the period up to the late- 1930s. This weaker Short Rains correlation in first half of the 20th century might be linked to SST data quality issues. However, this is notably altered in the early 1950s and after, particularly with EIO ( Figure 4; Camberlin and Philippon, 2001;Manatsa and Behera, 2013;Nicholson, 2015). In addition, all remote indices indicate a change in the relation with EEA-SR-SPI around the mid-1970s, most likely related to a "climate shift" in this period, with global consequences (Meehl et al., 2009). The possibility that such abrupt changes in the Short Rains teleconnections are due to inclusion/exclusion extreme years of 1961 and 1997 (Figure 3a) was tested by replacing the mean (i.e., a zero anomaly) and redoing the correlations. Overall, such replacement does not affect the significant shift in the mid-1970s and multidecadal changes in the relationship between the Short Rains and DMI for short time windows (not shown). This is consistent with Nicholson (2015).

| Changes associated with Short Rains teleconnections in triangle representation
It is evident from Figure 4 that the correlations are nonstationary to different degrees. To further examine such nonstationary associations, the statistical relationship over both various time window lengths for fixed starting years and starting years for fixed time window length in triangle representations are explored. The first analysis is performed for DMI and EIO indices representing IOD and related to interannual variability in the tropical IO. Figure 5 demonstrates the linear correlation coefficients between EEA-SR-SPI and the two indices for the period 1901-2013. The triangle representation of the correlations can be interpreted as follows: first, moving from the starting year on the left ordinate horizontally to the right represents correlations at the starting year with increasing time window length. Second, moving from the top vertically downward yields correlations for starting years with a fixed moving time window denoted on the lower and upper abscissa. All correlations with 5% significant level based on SDOF (EDOF) are indicated by contours (stippling).
No significant connection was found between the DMI and the Short Rains for time series of less than 20 years for starting years between 1901 and 1951, with an exception of starting years near 1918 (Figure 5a). However, a relatively strong coupling between the two occurred after 1951 until 1978. Another abrupt shift in correlation occurred around 1983, signalling a significant change in the relationship. However, the relation reversed again near 1987 to the recently observed even stronger coupling. Our nonstationarity test applied to an 11-year window length further confirmed that the change in Short Rains and IOD link is significant at 1% significance level. Previous studies also mention a similar change in the relation between the Short Rains/IOD, particularly during 1918during , 1961during , 1983during , and 1997during (Clark et al., 2003Manatsa and Behera, 2013;Nicholson, 2015Nicholson, , 2017. These periods have already been confirmed as times when significant shifts occurred in the IO (Manatsa et al., 2012). However, a significant positive correlation for windows longer than 30 years is seen for starting years after 1918. East Africa receives above (below) normal rainfall during + (−) IOD years (Figure 5a). The PVE for such long sample periods is greater than 16%; particularly, in recent years with a short length of the time window the PVE increased to 80%. An analogue pattern is identified for high-pass filtered time series (Figure 5e). A similar analysis is shown for low-pass filtered data with a half power period of 8 years (Figures 5c). For this, window lengths equal or longer than 40 years are used. The correlation coefficients must be interpreted with caution due to the very few cycles of decadal variations that can occur in a 40-to 113-year-long time series. Positive correlations significant at 5% confidence level between low-pass filtered time series of DMI and the Short Rains are revealed in Figure 5c. This decadal link between DMI and the Short Rains is most likely a signal of decadal timescale PO influence on the eastern IO (Figures 5d, 6c). In fact, some previous studies have suggested that the PO decadal variability is Linear correlation coefficients between EEA-SR-SPI and standardized EIO index (right panels) and DMI (left panels) for the period 1901-2013 using a combination of all starting years (ordinate) and time window lengths (abscissa). (a, b) The unfiltered components, (c, d) low-frequency components, that is, longer than 8 years with only window lengths of ≥40 years displayed, and (e, f ) high-frequency components, that is, shorter than 8 years. Significance at 5% level based on the SDOF and EDOF is indicated by contours (stippling). All correlation coefficients are calculated with linearly detrended time series [Colour figure can be viewed at wileyonlinelibrary.com] related to the interdecadal variability of Short Rains (Omondi et al., 2013;Bahaga et al., 2015). A comparable but negative correlation between EIO index and EEA-SR-SPI is presented in Figure 5b. However, no significant correlations are found for time windows less than 60 years and for starting years until 1951. For a shorter time window, a shift in the relationship occurred around 1931, where the relationship changed to an insignificant positive correlation, then reversed to a significant negative correlation around 1951 and after 1991 (Figure 5b). The DMI (EIO) show positive (negative) significant correlations to low-pass filtered EEA-SR-SPI (Figure 5c, d), but no significant connections are found for WIO (Figure 6d), which signals that low-frequency association between DMI ( Figure 5c) and EEA-SR-SPI is most likely explained by interaction between western PO and eastern IO on a decadal timescales (Xie and Annamalai, 2002;Black et al., 2003). These studies suggest the dynamics behind ENSO and IO teleconnections, an interaction which goes beyond the scope of this study.
In order to investigate the variation of EEA-SR-SPI teleconnection to N34 and WIO, further analyses are performed. Figure 6a shows the correlation patterns between the N34 index and EEA-SR-SPI. N34 is positively correlated with East African Short Rains for long sample periods. The PVE is greater than 16%, which is low but statistically significant. This suggests that N34 plays a role in the interannual variability of this Short Rains, but that the N34 influence is mediated through a concurrent warming in the WIO as extracted by partial correlation analysis ( Figure S1b). Figure 6a further reveals multidecadal changes in the relationship between N34 and EEA-SR-SPI. For time windows of less than 20 years prominent changes occurred for starting years around the 1920s, in the late 1940s, and around the early 1970s. Another abrupt shift in correlation occurred in recent decades. These changes coincide with the negative phase of Interdecadal Pacific Oscillation (IPO, e.g., Figure 12) and might point to a multidecadal change in Pacific SSTs in 1998-1999, which is consistent with Figure 8c of Lyon et al. (2014); that is akin to the Pacific climate shift during 1976-1977(Trenberth and Hurrell, 1994Deser et al., 2004;Meehl et al., 2009). Thus, the statistically significant correlation between N34 and EEA-SR-SPI largely stems from the low-frequency fluctuations (Figure 6c). To further substantiate the role of IPO on Short Rains variability, the correlation between EEA-SR-SPI and IPO index is shown in Figure 7. For time windows of less than 20 years, a clear decadal change in the relation is evident around the early 1920s, the late 1940s, around 1970s, and in recent decades (Figure 7a), consistent with previously reported changes in the PDO (e.g., Mantua et al., 1997). In addition, while there are some differences, this correlation pattern is similar to N34, with N34 having the most pronounced signal (Figure 6a). This is further confirmed by a significant association between IPO and EEA-SR-SPI in the low-pass filtered data (Figure 7b). The evidence presented here indicates strong connections between the East African Short Rains with EIO, N34, and IPO on the decadal timescales, corresponding to a decadal interactions between the western Pacific and eastern IO (Xie and Annamalai, 2002;Black et al., 2003). Furthermore, the result found here reaffirms that of previous studies which showed the association between the East African Short Rains and ENSO on the interannual timescale is actually more the result of an indirect forcing by the PO on the IO. On the other hand, the Short Rains decadal variability is mainly driven by decadal variability in the PO (Figures 6c, 7b). Overall, we have shown that the western IO is a major driver of East Africa Short Rains interannual variability, while the low-frequency variability in eastern Indian and the PO are highly associated with Short Rains decadal variability.

| Multidecadal changes associated with Kiremt Rains Teleconnections
In this study, we have revisited the relation between SSTbased remote climate anomalies and Kiremt Rains over the Ethiopian Highlands. An in-depth analysis was executed to examine the non-stationary nature of Kiremt Rains teleconnections. Figure 8 shows the 20th century time evolution of the correlation coefficients between mean ETH-KR-SPI and three SST-based remote climate indices for a 31-year sliding window. It is evident from Figure 8 that there are nonstationary correlations to different degrees. One salient example is the Atlantic Niño 3 (ATL3) index. It changed from significant negative correlations of nearly −0.6 in the late 1940s to a nearly significant positive correlation of about 0.35 in the recent decade. Based on the nonstationarity significance test, the changes in the running correlation between the ATL3 index and ETH-KR-SPI are significant. This suggests a change from a simultaneous occurrence of a warm tropical east Atlantic and dry Kiremt Rains in the early part of the investigation period to an opposite behaviour in recent decades (Figure 8). The latter is consistent with Segele et al. (2015) based on modelling experiments for the period of 1971-2000, who showed that warmer AO sea surface temperature anomalies (SSTA) generally brought wetter conditions to northern Ethiopia. In addition, Srivastava et al. (2018) found similar multidecadal evolution in the teleconnection between ATL3 and West Africa Monsoon (WAM) rainfall. They showed that after 1980, El Niño (La Niña) tend to occur simultaneously with negative (positive) ATL3 and this in turn reduced (enhanced) WAM by increased south-westerly moisture flux.
A second insignificant but non-stationary correlation is that between the AMO index and ETH-KR-SPI. For 31-year windows, it changes from significant negative correlations of nearly −0.35 in the early 1940s to the significant positive correlation of about 0.5 in the late 1960s (Figure 8). For the given window length of 31 years, the index cannot represent the AMO which fluctuates on timescales of about 60 years (Enfield et al., 2001;Diatta and Fink, 2014). Thus, this nonstationary correlation might reflect the changes in the covariability between interannual SST anomalies in the North Atlantic and JJAS rainfall over the Ethiopian Highlands. As indicated by Figure 3, summer rainfall over the Ethiopian Highlands has a fingerprint of Sahel multidecadal variability. This is consistent with the changes in the running correlation between AMO and ETH-KR-SPI (Figure 8), suggesting that a fraction of the JJAS multidecadal rainfall variability is associated with SST changes in the Northern Hemisphere AO. Jury (2010) also found a modest association between areas of the summer rainfall over Ethiopia and SST variability associated with the AMO. Nonetheless, this requires further modelling experiments.
The sliding correlation coefficient between the Kiremt Rains and N3 (used here instead of N34 due to somewhat stronger correlations) also varies over the study period between statistically insignificant values of nearly −0.32 to the robust and significant correlation value of nearly −0.67 for 31-year windows centred around the late-1970 ( Figure 8). Nonetheless, it may be noted that the significant but lower correlation before 1930s may be related to data quality issues both in SST and Ethiopian rainfall data sets. In general, the observed relations reveal that the interannual variability of the Kiremt Rains is strongly associated with eastern Pacific SST anomalies. However, the time evolution of correlation coefficients did not pass the significance test of non-stationarity. It is evident that the occurrence of drought (excessive rainfall) over this region is associated with El Niño (La Niña) events ( Figure 8). This is consistent with the findings reported in previous studies (e.g., Korecha and Barnston, 2007;Segele et al., 2009aSegele et al., , 2009bDiro et al., 2011a). In addition, it is important to underline increasing influence of ENSO on Kiremt Rains over Ethiopian Highlands in the recent period (Figure 8). This is consistent with Suárez-Moreno et al. (2018) who suggested that the nonstationary teleconnection between WAM and ENSO is likely due to the varying amplitude of the latter. Some recent decades have seen strong ENSO variations, which would suggest an enhancement of teleconnections in our case. However, our results show that such an enhancement cannot be statistically proven in the period under study.

| Changes associated with Kiremt Rains teleconnections in triangle representation
The sliding correlation in Figure 8 suggests instabilities of the links between Kiremt Rains and remote climate indices. In this subsection, we performed further in-depth analyses to examine the non-stationary in triangle representations.  Figure 9a demonstrates the correlation patterns between the N3 index and ETH-KR-SPI with the unfiltered data. No significant connection could be derived between the N3 and the Kiremt Rains for time windows less than 20 years for starting years between 1901 and 1934, with an exception for starting years near 1908. However, a relatively strong coupling between the two occurred after the starting year 1935, and between 1980 and 1990. These shifts were associated with warm IPO phases, and ENSO variability was in phase with the IPO/PDO ( Figure S2). For N3, a significant positive correlation for time windows longer than 30 years is seen for all starting years. Ethiopian Highlands receives above (below) normal rainfall during La Niña (El Niño) years (Figure 9a). This correlation is significant, stable with PVEs larger than 50%, and mainly based on high-frequency variability (<8 years). A similar pattern with a pronounced background decadal change is identified for high-pass filtered time series, exhibits that N3 explains the interannual variability of the Kiremt Rains (Figure 9e). Moreover, there is also significant correlation for low-pass filtered time series (Figure 9c), suggesting a potential decadal link between ENSO and the Kiremt Rains. Similarly, a linear correlation between the ATL3 SSTs and the Kiremt Rains for all combinations of starting years and time window lengths is illustrated in Figure 9b. On the interannual timescales, ATL3 is negatively correlated with the Kiremt Rains for starting years between 1931 and 1960 (Figure 9b), although the negative correlation coefficients are smaller and less significant. This correlation changes from negative to positive in the late 1960s, while being insignificant after this reversal. What causes the significant non-stationarity correlation between ATL3 and ETH-KR-SPI?
In order to answer this question and further investigate and understand the mechanism by which ATL3 influences the interannual variability of the Kiremt Rains, we performed a composite analysis. Similar to the approach followed in correlations between ATL3 and ETH-KR-SPI. The years taken for composites can be seen in Table S1. Figure 10 , 1934-1964) (Dry, 1934-1964) (Wet, 1969 (Dry, 1969-1999) precipitation over Africa, HadISST, and 200-hPa wind from ECWMF ERA-20C reanalysis representing patterns for: (a) wet years in 1934-1964, (b) dry years in 1934-1964, (c) wet years in 1969, and (d) dry years in 1969. Wet Kiremt Rains in the 1934-1964 period occurred simultaneously to WAM dipole years with a cold ATL3 (Figure 10a), and vice versa (Figure 10b). On the other hand, Kiremt wet years co-occurred with warmer ATL3 during 1969-1999, but without coherent dipole pattern (Figure 10c). A corresponding regime shift in 1968 within the WAM system has apparently changed the influence of the Atlantic Niño on Kiremt Rains at same time (e.g., Diatta and Fink, 2014;Nicholson et al., 2018). Notably, most wet years did occur with weak La Niña in conjunction with warmer ATL3 but most dry years were strongly associated with El Niño (Figure 10c, d). This is consistent with Figure 8, showing a clear change in the relation to ENSO. In addition, it can clearly be seen that there are easterly wind anomalies at 200 hPa in excess Kiremt Rains years and a westerly anomalies in deficit Kiremt Rains composites ( Figure 10). The corresponding 850-hPa wind composites are shown in Figure S3.
As to be expected, they show that easterly (westerly) lowlevel anomalies over northern tropical Africa are related to dry (wet) Kiremt seasons. Corresponding to the well-known shift in the relation between ATL3 and West African rainfall in the late 1960s (e.g., Nicholson et al., 2018), namely, a dipole (monopole) response to warm ATL3 before (after) this date ( Figure S3b and c), the low-level winds did change direction during warm ATL3 years between these periods. This zonal wind change, the reason of which remains elusive, went along with the changed response of the Kiremt Rains to ATL3 forcing. To further examine the evolution of the linear statical relationship between IPO and the Kiremt Rains, correlation in triangle representation are shown in Figure 11. The IPO index is negatively correlated with ETH-KR-SPI for long sample periods starting between the 1900s and 1950s and more recently also for shorter sample period lengths ( Figure 11). Significant association between IPO and ETH-KR-SPI in the low-pass filtered correlation pattern (Figure 11b) is similar to N3 (Figure 9c), with N3 having the most pronounced signal. This might be an indication of interdecadal modulation of the IPO upon the impact of ENSO on the Kiremt Rains. For instance, the frequent dry Kiremt Rains during 1982-1984 appears to have occurred when El Niño events coincide during warm phases of IPO (Figure 3; Lyon, 2014). Thus, the robust and significant correlation between IPO and the Kiremt rainfall index in a longer time window may reflect ENSO-induced variability of the IPO or superposition of these two variability modes (Figures 3, 9a, e, 11-12;Newman et al., 2016). In summary, despite the increase in ENSO influence on Kiremt Rains in recent decades, our simple statistical analysis of centurylong observations do not exclude that this increase is due to yet short sampling period and that the underlying teleconnection is not changed.

| Multidecadal changes associated with GHA Long Rains teleconnections
As mentioned in the Introduction, understanding interannual variability and predictability of the Long Rains has been a challenge for researchers and forecasters. In this study, we have considered all SST-based remote climate anomalies suggested by previous studies and explored (Table 1). Consistent with Hastenrath et al. (1993) and Camberlin and Philippon (2002), we did not find any significant correlation between indices representing ENSO interannual variability and Long Rains over GHA (not shown). Recently, Funk et al. (2014) suggested western Pacific SST gradient and SSTs in the central IO (Table 1) as a predictor of the first and the second principal components of the interannual variability of MAM Long Rains. As a consequence, we have analysed the evolution of the concurrent correlation coefficients between mean GHA-LR-SPI and Western Pacific Gradient (WPG) and Central Indian Ocean (CIO) (cf. Table 1) for a 31-year sliding window. However, no significant correlations are found (not shown). This is consistent with Liebmann et al. (2014) who have highlighted that the sources of interannual variability of the Long Rains have been difficult to pin down, in part because of MAM rainfall appears to be only weakly constrained by SST anomalies on interannual time scales.

| Changes associated with Long Rains teleconnections in triangle representation
Consistent with previous studies, we did not find any SSTbased indices significantly correlated to Long Rains and hence it is difficult to identify the sources of interannual variability. In this section, we show further analyses to examine Long Rains decadal variability and teleconnections. As stated in the Introduction, there has been more frequent drought in Long Rains regions over GHA in recent years. So far there are three hypotheses developed in order to explain the recent decline. One of those suggestions attributed this with decadal variability in the PO; however, no specific indices of decadal ocean variability were suggested. Here, we use the definition of Deser et al. (2004) to create PDO and IPO time series (Table 1) and present further evidence on the significant decadal relationship between indices of Pacific multidecadal variability and the MAM Long Rains. Figure 12 shows detrended PDO and IPO time series and detrended GHA MAM rainfall anomalies from GPCC. The time series cover our investigation period and all indices have been smoothed using a 9-year running average. The detrended Long Rains time series exhibits significant correlation with both PDO and IPO indices ( Figure 12). This is in agreement with Lyon (2014) and Yang et al. (2014). Figure 12 reaffirms that during warm PDO phases, for instance, 1924PDO phases, for instance, -1945PDO phases, for instance, and 1977PDO phases, for instance, -1998, GHA received above normal precipitation, while during cold phases from 1914-1925, 1947-1976, and 1999-present precipitation over the GHA was below normal.
To further elucidate the decadal relationship between PO SST anomalies and the MAM Long Rains anomalies, Figure 13 shows the regression pattern of SST anomaly onto the GHA Long Rains anomaly relative to 1950-1979 climatological base in MAM over the period of 1901-2013. The long-term trends for the whole investigation period in both SST and precipitation were removed before the regression calculation. To show teleconnections at a decadal time scale a 9-year running average has also been applied to both the observed SST and the precipitation, and dotted shading stands for the regression coefficients significant at the 5% level. The regression coefficients spatial pattern shows a typical horseshoe shape characterized by cold SST anomalies in the western North Pacific encircled by warm SST anomalies in the eastern part of the basin extending to the central and eastern tropical Pacific (Figure 13), corresponding the positive phase of PDO spatial pattern (Deser et al., 2012;Farneti, 2017). As can be clearly identified from Figure 13, the wet phase of the Long Rains is associated with a pattern of the positive phase of PDO, while the dry phase of the Long Rains is related to a pattern of the negative phase of PDO  Lyon, 2014;Yang et al., 2014). Overall, the regression analysis further supports Figure 12 and previous suggestion that the post-1998 MAM Long Rains decline is a part of decadal climate variability over GHA, remotely driven and associated with the negative phase of the PDO.
To illustrate the evolution of Long Rains decadal teleconnection, we have used correlations in triangle representation for PDO ( Figure 14a) and IPO (Figure 14b) and the GHA-LR-SPI. For all starting years, the correlations between the GHA Long Rains and Pacific SSTAs represented by PDO and IPO are positive and significant based on both SDOF and EDOF tests. Consistent with the relation shown in Figures 12 and 13, there are relatively high correlation coefficients greater than 0.6. PVE ranges from 25 to 64% and the association stems from low-frequency variability (>8 years) for both PDO and IPO (Figures 14a, b). Thus, our analysis support that the recent Long Rains decline is part of natural decadal fluctuation and strongly associated with multidecadal variability in the PO. In addition, the significant correlation between the rainfall and indices of decadal variability in the PO hint at useful predictors of rainfall variability on decadal timescales.

| DISCUSSION AND CONCLUSION
The main aim of this study was to revisit previously known teleconnections influencing three major rainy seasons over GHA, the Long Rains, the Short Rains, and the Kiremt Rains over the 1901-2013 period. In addition, the evolution of statistical relationships between remote climate anomalies and rainfall indices are examined both over varying time window lengths for fixed starting years and varying starting years for fixed time window length in triangle representation. Focusing on the stationarity, linear correlations were studied for unfiltered, low-pass (larger than 8 years), and high-pass (less than 8 years) filtered detrended time series. In addition, a composite analysis was used to investigate the physical mechanisms explaining the non-stationary behaviour of the correlation between ATL3 and Kiremt Rains. To this end, monthly precipitation datasets were used to create standardized precipitation indices over homogeneous rainfall regimes in the respective sub-regions. SST datasets were used to construct indices for remote climate anomalies (Table 1). The major findings related to the three rainy seasons are: • Consistent with previous studies, the interannual variability of the Short Rains is more closely linked to the IOD than to ENSO. The East Africa region receives above (below) normal rainfall during positive (negative) IOD events. We found abrupt changes in the relation between the Short Rains and IOD, particularly around 1918IOD, particularly around , 1951IOD, particularly around , 1983IOD, particularly around , and 1987. Our non-stationarity test applied for 11-year window length further confirmed that this non-stationarity in the link is significant at the 1% level. This is consistent, but expands the findings of Camberlin and Philippon (2001), Clark et al. (2003), and Nicholson (2015). Moreover, ENSO shows strong correlation with the Short Rains year-to-year variability, yet the ENSO influence is mediated by the co-occurrence of in-phase oscillation of the IOD. Our result further reveals that the low-frequency variations of the Short Rains are most likely explained by decadal variability in the PO, as it is demonstrated by significant correlations to the low-pass filtered time series of ENSO and IPO indices. • We found a significant non-stationary correlation between ATL3 and the Kiremt Rains over the Ethiopian Highlands. It changed from significantly negative to positive after in the 1960s. This is in agreement with a WAM regime change in 1968 (e.g., Diatta and Fink, 2014;Nicholson et al., 2018). The composite analysis also revealed that in the 1934-1964/1969-1999 period, wet Kiremt Rains over the Ethiopian Highlands tend to simultaneously occur with cold/warm ATL3 region ( Figure 10; Diro et al., 2011a;Segele et al., 2015). Our analysis also reveals that interannual variability of the Ethiopian Kiremt Rains is strongly correlated with the ENSO indices. This association is significant and stable with PVEs larger than 50% in recent periods and mainly based on high-frequency variability (<8 years). From simple time series analysis, a statistically significant larger influence of ENSO on Kiremt Rains in recent decades cannot be inferred. • On the other hand, the IPO index is negatively correlated with Kiremt Rains for long sample periods. The negative (positive phase) of IPO corresponds to frequent occurrence of wet (dry Kiremt Rains), and might reflect the interannual ENSO events occur in phase with the IPO. For instance, the frequent dry Kiremt Rains during 1982-84 appears to have occurred when El Niño events coincide during warm phases of IPO (Figure 3; Lyon, 2014). This might further suggest a potential modulation of ENSO forced interannual teleconnections by the background Pacific decadal variability (Wang et al., 2014;Newman et al., 2016). • The recent Long Rains decline over GHA is part of natural decadal fluctuations and significantly associated with indices representing multidecadal variability in the PO. Both the PDO and IPO show a positive and significant correlation to the Long Rains. PVE ranges from 25 to 64% and the association stems from low-frequency variability (>8 years) for both PDO and IPO. The result presented is in agreement with Lyon (2014) and Yang et al. (2014). This significant correlation between the Long Rains and IPO/PDO might offer a useful inference in predicting the Long Rains at the decadal timescale.
Reasons for non-stationary behaviour in interannual teleconnections can be manifold. Apart from data inhomogeneities, changes in the background climate state at interdecadal time scales are another potential cause. For example, Suárez-Moreno et al. (2018) suggest this is one reason for non-stationary teleconnections between the Eastern Mediterranean SSTs and West African monsoon rainfall. Alternatively, the amplitude of ENSO is suggested as a reason for variations in the link to West African monsoon rainfall. Changes in the East African Short Rains teleconnections have been discussed by Camberlin and Philippon (2001), Clark et al. (2003), and Nicholson (2015). Camberlin and Philippon (2001) suggested the change in large-scale circulation anomalies and inclusion/exclusion of extreme years as justifications for the changes in the apparent strength of the teleconnections. Here we propose that the multidecadal background state of the PO (i.e., IPO, PDO) is one possible explanation for the non-stationary of Short Rains interannual teleconnections. The correlation between the Short Rains and N34/DMI appears to be modulated and enhanced during the cold phase of IPO, while the impact is insignificant during the warm IPO phase. This is consistent with the damping effect discussed in Wang et al. (2014). Some of the abrupt changes found in the present study may also coincide with 1976-1977 climate shift (Trenberth and Hurrell, 1994).
The present comprehensive study stresses the usefulness of SST-based statistical indices for statistical prediction models of GHA rainfall. But at the same time, it also highlights that non-stationary behaviour as well as correlations between indices on interannual and decadal times scales needs to be taken into account. A similar finding has been noted by Diatta and Fink (2014) for the West African monsoon. Further studies should consider individual months (Vigaud et al., 2017) and the spatial coherence of the signals within the three regions considered (e.g., Camberlin and Philippon, 2002). Furthermore, our results can be taken to design modelling studies to gain insights into the physical mechanisms responsible for the internanual and decadal teleconnections and their non-stationarities.