Coupled network analysis revealing global monthly scale co-variability patterns between sea-surface temperatures and precipitation in dependence on the ENSO state

The Earth’s climate is a complex system characterized by multi-scale nonlinear interrelationships between different subsystems like atmosphere and ocean. Among others, the mutual interdependence between sea surface temperatures (SST) and precipitation (PCP) has important implications for ecosystems and societies in vast parts of the globe but is still far from being completely understood. In this context, the globally most relevant coupled ocean–atmosphere phenomenon is the El Niño–Southern Oscillation (ENSO), which strongly affects large-scale SST variability as well as PCP patterns all around the globe. Although significant achievements have been made to foster our understanding of ENSO’s global teleconnections and climate impacts, there are many processes associated with ocean–atmosphere interactions in the tropics and extratropics, as well as remote effects of SST changes on PCP patterns that have not yet been unveiled or fully understood. In this work, we employ coupled climate network analysis for characterizing dominating global co-variability patterns between SST and PCP at monthly timescales. Our analysis uncovers characteristic seasonal patterns associated with both local and remote statistical linkages and demonstrates their dependence on the type of the current ENSO phase (El Niño, La Niña or neutral phase). Thereby, our results allow identifying local interactions as well as teleconnections between SST variations and global precipitation patterns.


Introduction
The heterogeneous latitudinal distribution of incoming solar radiation at the Earth's surface leads to a meridional imbalance with excess energy in the tropics and a corresponding deficit in the polar regions [1]. Meridional heat transport balancing this gradient is realized by two main mechanisms: via ocean currents transporting warm surface water towards the poles and cold subsurface water towards the tropics, and by large-scale atmospheric circulation patterns [2,3]. While the oceanic mechanism mostly dominates in the tropics, meridional heat transport in the mid-to-high latitudes contains a considerable atmospheric component [4].
In general, the ocean and atmosphere interact via different types of energy fluxes and thereby exchange heat via conduction and convection [5]. This oceanatmosphere coupling has been found to be essential for explaining various key climate phenomena, ranging from the seasonal cycle characteristics in the tropics over the El Niño-Southern Oscillation (ENSO) a e-mail: reik.donner@h2.de (corresponding author) to various patterns of decadal-scale climate variability [6,7]. Among others, dynamical patterns of seasurface temperatures (SST) generate energy and moisture fluxes from the ocean into the troposphere, air pressure variations and, hence, circulation patterns, and thereby emerging spatio-temporal precipitation patterns around the globe [4].
To better understand ocean-atmosphere coupling mechanisms, there have been many studies investigating which climate variability mode can trigger which spatio-temporal response patterns in different variables, with a particular focus on precipitation. For instance, several papers have investigated the influence of El Niño (i.e., the warm phase of ENSO) across Southeast Asia, Indonesia, Australia and the United States [8][9][10]. For Europe, the corresponding effect of the North Atlantic Oscillation (NAO) has raised particular interest [11,12], whereas for Australia, the precipitation responses to the changing states of the Interdecadal Pacific Oscillation (IPO) and Southern Oscillation Index (SOI) have been studied [13]. In addition to these mostly regionally focused studies, the effects of the ENSO, North Pacific Oscillation (NPO) and NAO on global precipi- tation patterns have been explored, revealing distinct regional response patterns to all these modes of climate variability [14]. However, especially precipitation extremes are most substantially affected by ENSO [15], which can be seen all over the world, including India, Africa, South America, the Pacific Rim, North America, and, more weakly, Europe [14].
The aforementioned studies have in common that certain index variables representing different largescale spatio-temporal climate variability patterns have been taken as references for studying linear responses in terms of correlation or regression maps. Notably, this viewpoint reduces the complexity of the problem under study considerably, but leads to a systematic loss of information on which parts (i.e., regions) of the addressed climatological field actually correlate strongly with, e.g., rainfall at a specific point in space. This general conceptual limitation can be overcome by considering the recently developed methodological framework of functional climate network analysis [16][17][18]. This approach takes the established framework of complex network theory as a basis for developing and employing a suite of characteristics describing the placement of strong mutual statistical associations among climate time series that may represent one or more climatological fields.
Complex network representations of climate variability can be used in parallel for data representation, analysis and visualization [19]. Thereby, they can help uncovering spatio-temporal patterns at multiple scales, ranging from local properties to global phenomena, which potentially improve our understanding of the overall dynamical organization of climate variability [20,21]. In the last 15 years, various studies have employed this modern framework to gain insights into essential patterns of climate dynamics. Statistical association between different climate time series, representing functional relationships between different parts of the climate system, has been typically characterized by either correlations or conceptually related nonlinear dependency measures [22][23][24] or based on the synchronous co-occurrence of extreme events [25][26][27][28]. In most cases, the individual time series correspond to nodes on some regular spatial grid covering either the entire globe or a specific region under study. Functional climate networks have been successfully used to hindcast extreme events, such as extreme precipitation in South America [29]), or to predict the occurrence of El Niño episodes [30]. In the context of the present work, we particularly highlight results on differential spatial organization patterns of climate variability associated with different ENSO phases [31,32] and the associated discrimination of different El Niño and La Niña flavors [33,34].
Most of the aforementioned climate network studies have exclusively focused on the dynamics within a single climatological field. Recently, some studies have extended this framework to the study of interdependences among different variables or atmospheric layers [35][36][37]. For this purpose, the idea of correlation-based single-variable networks has been thoroughly gener-alized into so-called coupled climate networks. This paper takes up a corresponding approach by focusing on the co-variability between SST and precipitation (PCP) during different seasons and ENSO phases, which is motivated by the known distinct teleconnection patterns triggered by El Niño and La Niña [38,39]. Thereby, we attempt to characterize the timedependent ENSO teleconnection patterns with global precipitation variability from a complex network perspective.
The present work shall serve as an initial case study highlighting what kind of information on SST-PCP coupling can be obtained with coupled climate networks by exploring relatively basic characteristics of the associated network connectivity structures. It is anticipated that follow-up work will continue along the outlined lines of research, for example, to differentiate between "local" (short-range) and "remote" (long-range) statistical SST-PCP linkages, which would be necessary for an appropriate process identification and associated attribution of the mutual connectivity patterns described in the course of the present work.
This paper is structured as follows: Sect. 2 describes the utilized data sets and methods. The obtained results are presented in Sect. 3. Here, we will first take a look at the local connectivity (i.e., the number of strong statistical associations) within one field with respect to the other variable (i.e., from SST to PCP and vice versa) as quantified by the so-called n.s.i. cross-degree. We will perform this analysis separately for strong absolute, positive and negative correlations, respectively. Subsequently, we employ the concept of network communities to the two coupled subnetworks, highlighting densely correlated regions across the two variables of interest. Finally, we provide a discussion of our findings in Sect. 4 and summarize them in Sect. 5.

Data
In this study, we utilize monthly globally distributed SST values from the Extended Reconstruction Sea Surface Temperature (ERSST v3) data set with a spatial resolution of 2 • × 2 • [40] and associated precipitation data with a spatial resolution of 2.5 • × 2.5 • provided by the Global Precipitation Climatology Project (GPCP) version 2 [41]. For the corresponding time series at each respective grid point, we first remove the mean annual cycle that otherwise would dominate the co-variability patterns among the studied series. In our analysis, we focus on the common time period of both data sets from 1979-2015. All grid points containing missing values (for example, due to temporary sea ice cover) have been removed from the data prior to our further analysis. As a result, the total number of considered SST and PCP grid points is N s = 9456 and N p = 10, 368, respectively.
As already emphasized in the introduction, the aim of this study is to analyze the spatial placement of strong correlations between ocean (SST) and atmosphere (PCP) variability during different seasons, i.e., boreal winter (DJF), spring (MAM), summer (JJA) and autumn (SON). In addition, we further differentiate the time series into El Niño, La Niña and neutral ENSO phases based on the Oceanic Niño Index (ONI), whereby the two former are identified if the values of the 3-month running mean SST anomalies in the tropical Pacific Niño 3.4 region exceed 0.5K (fall below −0.5K) for a minimum of 5 consecutive months [42]. As a result, we associate 156 months with El Niño, 120 with La Niña, and 168 with neutral ENSO phases. In the following, we will consider the corresponding data subsets and, due to their relatively small sizes, focus on linear (Pearson) correlation as a statistical association measure, since more complex measures typically cannot be estimated reliably from the given relatively small sample sizes.
We note that we use here monthly data instead of daily ones for different reasons: First, the two considered variables of interest exhibit very different spectral properties, with relatively smooth variations (weak high-frequency variability) of SST opposing strong high-frequency variability in PCP. Second, global precipitation data at daily timescales would be largely dominated by zero values hampering their appropriate statistical analysis in the context of the present work. Taking both aspects together, we have reasons to believe that comparing SST and PCP time series at daily resolution would provide results that are markedly affected by synoptic (weather) timescales rather than unveiling statistically robust interdependencies on longer timescales.
While providing a reasonable degree of robustness (in a statistical sense) to the obtained results, our choice of using monthly data implies that we cannot perform an evolving network analysis (i.e., studying spatiotemporal co-variability structures of SST and PCP for sliding windows in time) in a statistically meaningful way. Accordingly, comparing patterns between, for example, different El Niño episodes is not possible for simple statistical reasons (each season comprises only three, and each ENSO episode just a few consecutive months, which does not allow for reliable correlation estimates based on such short segments of monthly data). Hence, the choice of monthly timescales taken for the reasons already discussed above does not allow us studying differences in spatio-temporal co-variability patterns of SST and PCP among different time intervals with the same type of ENSO state, but only provides an integral picture averaging over all periods on record with either El Niño, La Niña or neutral ENSO conditions during a certain season. We are aware of the fact that this perspective does not take the important aspect of ENSO diversity [43] (i.e., different El Niño or La Niña flavors) into account, which however cannot be circumvented within the specific analysis framework employed in this work.

Coupled climate network analysis
We consider the spatial fields of the two climate variables under study (tailored to the specific situation of interest as described above) as being described by two sets of univariate time series X (s) In what follows, we will use the superscript index s (p) to refer to properties of the SST (PCP) field.

Network construction
In a functional climate network, each time series is represented by a node embedded on the Earth's surface at the spatial position of the respective grid or measurement point. In our case, the corresponding links indicate strong linear (Pearson) correlations between pairs of such series. Since we have to consider all time series (respectively, their subsets corresponding to the same season and type of ENSO phase) from both fields, we first estimate the full N × N lag-zero correlation matrix between all N = N (s) +N (p) time series individually for each combination of season and ENSO phase. This correlation matrix contains all pairwise correlations among SST grid points, among PCP grid points, and between pairs of SST and PCP time series, which can be conveniently represented by a block structure Here, the two sub-matrices P (ss) (of size N (s) × N (s) ) and P (pp) (N (p) ×N (p) ) represent the correlation matrices of the SST and PCP fields, respectively, which consist of elements In the above equations, σ (•) and C (••) denote the respective standard deviation vectors and covariance matrices of the two considered fields defined in the standard way. In full analogy, we identify the sub-matrices P (sp) and P (ps) as containing all correlation coefficients between time series in SST and PCP (respectively, PCP and SST). Since we use only lag-zero correlations in this study, the two latter matrices can easily be transformed into each other by simple transposition.
To transform the correlation matrix P into the binary adjacency (connectivity) matrix of the associated network representation (describing whether (1) or not (0) a link between two grid points exists), we employ a fixed yet individual threshold to each of the submatrices, retaining only those pairs of grid points for which the correlation exceeds the respective threshold. Here, for the submatrices describing only one field (SST or PCP, respectively), we select these thresholds T (ss) and T (pp) so that a fraction ρ (s) = ρ (p) = 0.01 (i.e., 1%) of all pairs of grid points (the so-called (internal) link density) are represented by links in each of the two subnetworks. This threshold implies that we only consider correlations above the empirical 99th percentile of all correlations between all time series within each field. For the linkages between both fields (i.e., the so-called bipartite or cross-links between SST and PCP nodes in the coupled climate network representation), we proceed in a similar way by choosing a threshold T (sp) = T (ps) so that a fraction ρ (sp) = ρ (ps) = 0.005 < ρ (s) , ρ (p) (i.e., the corresponding cross-link density) is retained.
Taken together, we obtain the coupled network's adjacency matrix as where Θ(•) denotes the Heaviside function. Note that we use here the so-called extended adjacency matrix, which also contains self-loops since each SST (PCP) grid point is trivially correlated with itself with the maximal possible correlation coefficient of 1. This setting is motivated by the further use of area-weighted network measures as discussed below [44]. Regarding the block structure of the extended adjacency matrix, it may be noted that the two diagonal blocks correspond to adjacency matrices of two single-variable (unipartite) networks, while the off-diagonal blocks describe bipartite (two-variable) networks exhibiting links only between different types of nodes (corresponding here to the two different climatological fields under study). As a justification of our choices of internal and crosslink densities, we note that commonly, the correlations within one climatological field tend to be stronger than between two different fields. Therefore, taking just a global common threshold to all submatrices might effectively eliminate most cross-links, which are, however, the properties of interest in this work. To account for this, we enforce the presence of both, internal and crosslinks, by imposing the aforesaid constraints via the different link densities, ensuring that the number of cross-links is still generally smaller than the number of internal links within each of the two considered fields. This setting is relevant to be able to properly interpret the resulting mathematical structure as the adjacency matrix of two coupled yet still identifiable subnetworks.
In summary, the analysis as described above depends on two parameters, the internal link and cross-link densities, which determine the thresholds to the elements of the different correlation matrices forming the blocks of P. In this study, we present only the results for one specific choice of both values for the sake of brevity, which have been guided by choices from previous works on functional climate networks [33,34,37]. One may expect the presented results to change quantitatively as any of the two parameters (or even both) are gradually varied, while moderate variations of both link densities do not alter the presented results qualitatively. As systematic exploration of this aspect is however beyond the scope of the present work.

Area-weighted cross-degree
By computing certain statistical properties on the adjacency matrix, we can obtain different network measures that help understanding different aspects of the underlying network structure at either local or global scale [45,46]. In case of climate data, we, however, need to account for the spatially heterogeneous spacing between neighboring grid points, which commonly decreases towards the poles, to avoid an over-representation of the high latitudes in the climate network characteristics. This can be done by assigning a weight to each node corresponding to the spatial area it represents, i.e. [21] w n = cos λ n , with the latitudinal position λ n of node n on the grid. In a more general context, this idea of area-weighting has been formalized in terms of node-weighted network characteristics (so-called node-splitting invariant or shortly n.s.i. measures) [44].
In this study, we solely focus on the simplest possible network measure, the n.s.i. version of the cross-degree [35] as previously studied in [36,37], which characterizes the areal share of the entire globe that one node in a given field is connected to in the other field: where V i and V j denote the respective sets of nodes describing the two subnetworks (i.e., SST and PCP field).

Community structure
Communities are groups of nodes within a network that are relatively tightly connected with each other while being only weakly connected to the rest of the network [47]. This general feature can be characterized by different more specific properties (for example, the corresponding network modularity [47]) and has recently found its first applications in the context of functional climate network analysis [48][49][50][51]. Together with the practical relevance of communities as serving as meaningful subsystems that are solely determined by the overall connectivity structure [52], the broadness of this concept has made community detection an active area of research within complex network science [53,54]. As a result, there exist a great variety of algorithms for community detection based on different optimization targets. While most commonly employed methods do not distinguish between different types of nodes and, hence, implicitly assume unipartite network structures, there also exist algorithms that are specifically tailored for the case of bipartite networks [55,56]. Coupled network structures like those explored in the present work, comprising both uni-and bipartite constituents, have received less specific attention, but can generally be studied using the same methodologies as standard unipartite networks.
In the present work, we apply the concept of communities to the giant adjacency matrix A + containing all links within the SST and PCP fields as well as all cross-links between both fields. Thereby, a network community can combine nodes representing different variables, which is different from most previous applications of the network community paradigm. We will refer to such multi-variable dense connectivity structures as two-variable communities in the remainder of this work. As all other concepts employed in this work, it can be readily generalized to the case of multi-variable communities in coupled networks comprising more than two distinct climatological fields.
To obtain as robust communities as possible, we utilize the Infomap algorithm [57]. This algorithm is particularly well suited for the purposes of this study, since it attempts to divide weakly modular regions into small communities. In the context of functional climate networks, such small communities are more likely to represent physically meaningful spatial co-variability patterns because they are more strongly connected internally. Moreover, recent applications of Infomap in a climate network context have demonstrated that it provides certain climatologically meaningful features like the spatial coherency and connectedness of the identified communities [50,51].

Results
In the following, we will present the results of our coupled climate network analysis between SSTs and precipitation (PCP) during different seasons and ENSO phases.

Area-weighted cross-degrees
First, we investigate the n.s.i. cross-degree patterns formed by links from nodes in the SST field towards the PCP field (SST → PCP) and vice versa (PCP → SST). For example, a large n.s.i. cross-degree of SST → PCP indicates a grid point in the SST field that is strongly correlated with a large spatial area in the precipitation field.

Absolute correlations
We start by constructing coupled climate networks based on the largest absolute correlation values, indicating the generally strongest linear statistical associations between variables. Figure 1 shows the corresponding area-weighted cross-degree patterns for SST → PCP.
During boreal fall (SON), El Niño years are characterized by strong cross-connectivity (i.e., strong effects on global precipitation) arising from two regions in the eastern equatorial Pacific and between the western equatorial and southern subtropical Pacific (Fig. 1A). Similar yet weaker influence patterns are also observable during La Niña years (Fig. 1E), while neutral ENSO years are characterized by the absence of the western equatorial Pacific center of influence, which is replaced by two other patterns in the Indian ocean and equatorial Atlantic (Fig. 1I).
With the progressive evolution of El Niño or La Niña conditions, the boreal winter season (DJF) is characterized by particularly strong influence patterns. El Niño years primarily feature the main ENSO region in the central to eastern equatorial Pacific exhibiting particularly strong co-variability with the global precipitation field (Fig. 1B), while a second smaller pattern emerges over the westernmost Pacific. Interestingly, the corresponding area-weighted cross-degree patterns are even more strongly pronounced during La Niña conditions (Fig. 1F), highlighting in particular the eastern equatorial Pacific and the subtropical western to northern central Pacific. During normal ENSO conditions, we find the strongest connectivity in the central and western equatorial Pacific (Fig. 1J).
Boreal spring (MAM) experiences the greatest influence of SST on precipitation during El Niño conditions in a small equatorial band spanning from the central Pacific to the eastern Atlantic ocean (Fig. 1C). A second region with high cross-connectivity emerges in the southernmost Indian ocean. In comparison with this, during La Niña conditions (Fig. 1G), much larger regions are characterized by high cross-connectivity with the precipitation field, including again the ENSO region in the eastern equatorial Pacific, but also the southeastern Pacific (resembling the typical spatial location of the Amundsen Low, a key atmospheric circulation pattern in the Southern hemisphere extratropics). The corresponding patterns become more dispersed during neutral ENSO conditions, resembling closely those already observed during DJF (Fig. 1K).
Finally, the boreal summer season (JJA) experiences strong effects of eastern equatorial Pacific SSTs on precipitation after El Niño events (Fig. 1D), while only relatively weak patches of elevated cross-degrees are visible during La Niña and neutral ENSO years (Fig. 1H,L). Figure 2 presents the corresponding area-weighted cross-degree between the precipitation variability at one grid point and all grid points in the global SST field. In comparison with Fig. 1, the strongly connected patterns comprise much less grid cells than for the opposite direction (yet with overall larger cross-degree values) and are mostly confined to two small bands in the equatorial Pacific and close to the Antarctic coastline. Regarding the latter region, we may speculate about a possible association with the sea ice boundary, which might have distinct yet more regionally confined effects on precipitation formation via sharp local gradients in evaporation, albedo (i.e., radiation balance) and heat fluxes (possibly inducing local circulation systems). However, another possibility that we cannot rule out here would be inaccuracies of the SST data product in the corresponding regions. For the neutral ENSO phase, PCP patterns in the central equatorial Pacific and over the maritime continent are most strongly correlated with SST during the boreal winter season (Fig. 2J). During La Niña years, this season features a strong pattern over the central equatorial Pacific (Fig. 2F), but also another smaller patch over the northern part of Brazil. In turn, during El Niño years (Fig. 2B), we observe only a single pattern over the central equatorial Pacific, which is however considerably weaker than during La Niña conditions.
During the boreal spring of neutral ENSO years (Fig. 2K), the dominating cross-degree patterns spread over a larger region in the equatorial Pacific. During La Niña phases (Fig. 2G), the patterns closely resemble those during DJF but appear somewhat less prominent. Similarly, during El Niño years (Fig. 2C), the observed Notably, during both, boreal fall and summer, the cross-degree patterns during neutral ENSO years are very weak and mostly cover some tropical regions (Fig. 2I, L). The same applies to La Niña phases, while there are still some stronger patterns over the eastern equatorial Pacific yet much weaker than during DJF and MAM (Fig. 2E, H). In turn, during El Niño conditions, the main cross-degree patterns stay confined to the central-to-eastern equatorial Pacific ( Fig. 2A, D).

Positive correlations
While having previously considered strong absolute correlations, we have mixed up information on regions where, e.g., strong positive SST anomalies trigger wet conditions with such exhibiting precipitation deficits. Since both have distinctively different implications for the affected ecosystems and population, we now turn to a separate analysis for positive and negative correlations, respectively. For the positive correlations, Figs. 3 and 4 reveal strong similarities with the networks constructed based on absolute correlations.
The high degree of similarity between networks based on strong absolute and positive correlations is somewhat expected, since a vast part of the identified strongly influential SST regions (except for some regions in the Southern Hemisphere extratropics) are located in the tropics and subtropics. In those regions, we may expect positive SST anomalies to result in excess evaporation, leading to enhanced convective activity, cloud formation and, subsequently, excess precipitation. Via cascading moisture recycling, those processes can even affect significantly larger regions of the globe, as well known, for example, for the case of the Amazon rainforest receiving most of its precipitation via such a mechanism [58][59][60]. Another relevant process in this regard are recurrent intraseasonal wave phenomena in the tropics, such as the Madden-Julian Oscillation (MJO) [61] or the Boreal Summer Intraseasonal Oscillation (BSISO) [62], which exert a prominent modulation especially on monsoonal rainfall in various parts of the low latitudes.

Negative correlations
Finally, we study the coupled climate networks describing only the strongest negative correlations. During boreal fall, the neutral ENSO phase is characterized by elevated cross-degree patterns in the northern tropical Atlantic, western equatorial towards southwestern Pacific and central Indian ocean (Fig. 5I). The associated cross-degree patterns in PCP are rather diffuse and concentrate over the central equatorial Pacific (Fig. 6I). During La Niña (Fig. 5E), high n.s.i. crossdegrees for SST → PCP are found in the western equatorial to southwestern Pacific and weakly over the Indian Ocean, while the corresponding largest crossdegrees for PCP → SST are still concentrated over the central equatorial Pacific (Fig. 6E). El Niño phases (Fig. 5A) are characterized by two patterns of elevated cross-degrees over the eastern equatorial Pacific and between the eastern equatorial and southwestern Pacific. The corresponding PCP patterns (Fig. 6A) are mainly concentrated over the western equatorial Pacific and the maritime continent. Interestingly, the previously observed circumpolar Antarctic precipitation pattern is markedly correlated with SST regardless of the ENSO state.
In general, we emphasize that the cross-degree patterns obtained based on positive and negative correla- tions look surprisingly similar. This points to the fact that strong positive SST anomalies in certain regions can trigger marked atmospheric circulation changes (as known for ENSO related SST anomalies along with intensity changes and dislocation of the Walker circulation [63,64]), which may cause excess rainfall in some regions (positive correlation) while dry conditions in others (negative correlation). Note that the presented cross-degree patterns do not distinguish between the destinations of links emanating from the SST field. For the latter purpose, our analysis calls for further indepth investigations by either utilizing measures that capture spatially explicit information, like link distance or directionality related network properties [50,51,65], or systematically identifying pairs of source and target regions of different types of links (representing positive and negative correlations, respectively) [51,66]. Such more detailed analyses are, however, beyond the scope of the present work.

Communities
As a final step of our analysis, we take a look at the twovariable community structure of the coupled SST/PCP network. Since we expect the most interesting differences between the three types of ENSO phases to arise during boreal winter (where El Niño and La Niña peak), we present here exclusively the results for this season. Moreover, to allow for some interpretable visualization, we disregard all communities that include less than 1% of the nodes of the coupled network. Figure 7 shows the identified two-variable communities for the different ENSO phases. We first recognize that these communities in most cases correspond to spatially contiguous regions in the SST field (with a few exceptions where individual communities contain more than one substantially large region), whereas they appear more fragmented in the PCP field. In general, in the SST field, many of the larger two-variable communities involve those nodes that have been previously identified as high cross-degree patterns, suggesting that the associated communities indeed also include regions in the PCP field (cf. the different colors in Fig. 7).
When comparing especially the respective contributions of the SST field to the individual communities, we observe that El Niño phases exhibit a larger number of relevant communities than La Niña or neutral ENSO phases. This observation could have two possible reasons. On the one hand, the linkage structure within the SST field (representing connections that either directly link different SST grid points or are indirectly mediated via PCP grid points) is more fragmented during El Niño phases. Indeed, it is well known that El Niño leads to large-scale synchronization of climate variability by triggering coherent climate responses in distinct regions across the global surface temperature and precipitation fields. On the other hand, and in line with the latter fact, the higher number of communities could also be interpreted as reflecting a generally more homogeneous network structure, for which already very minor differences in the link placement may cause the Isomap algorithm to split larger groups of nodes into different communities.

Discussion
To study ocean-atmosphere coupling from a complex network perspective, we have generated coupled network representations of global SST and precipitation fields and evaluated the associated n.s.i. cross-degree and two-variable community patterns.
In boreal fall, global precipitation patterns are especially linked with SSTs in the Indian Ocean, the Caribbean Sea, the Barents Sea and a region known for explosive cyclogenesis [67] near Newfoundland. In the subsequent winter and spring, the n.s.i. cross-degree for absolute correlations (Figs. 1 and 2) shows that even during the neutral phase of ENSO (which is often considered less relevant in terms of teleconnectivity), SSTs in the central equatorial Pacific are strongly linked to precipitation. Interestingly, the corresponding SST variability is accompanied by both, positive and negative correlations with precipitation patterns, whereby positive correlations can be primarily found in the eastern central Pacific and negative ones in the western central Pacific [64]. In boreal summer, these pattern are less pronounced and more homogeneously distributed, which could point to a stronger relevance of local-scale (convective) processes than remote teleconnectivity. All aforementioned regions are moreover also associated with negative correlations. Positive SST correlations with precipitation can be found in the central Pacific, while negative ones arise particularly over the maritime continent and Central America. Moreover, SST also strongly affects precipitation close to the Antarctic coastline. During La Niña phases, the patterns observed for neutral ENSO phases are further enhanced. This particularly applies to the cross-degree values in the equatorial Pacific and especially in boreal winter and spring. A few main patterns can be identified, including an elongated pattern in the central and eastern tropical Pacific, an arc-shaped pattern in the western Pacific and some patterns close to the Californian coast and in the Amundsen Sea in the South Pacific [68,69]. Like during the neutral phase, the western Pacific SST anomalies are linked with precipitation via negative correlations and their eastern Pacific counterparts via positive correlations. Strongly affected precipitation variability can be localized in a large area in the central Pacific (with positive correlations) and over the maritime continent and the Caribbean Sea (negative correlations). Precipitation is mainly associated with one large community involving the SSTs in the central-toeastern equatorial Pacific.
During El Niño periods, the Walker circulation is markedly weakened, which causes warmer SSTs in the eastern equatorial Pacific and a cooling in the western equatorial Pacific. These changes in the Walker circulation result in excess precipitation in the central-toeastern tropical Pacific and less rain over the maritime continent [64]. However, especially the SST anomalies in the eastern equatorial Pacific are strongly linked to precipitation. The precipitation anomalies are correlated with SST anomalies primarily in the central eastern Pacific, the eastern equatorial Pacific and over the maritime continent. Like during La Niña and the neutral ENSO phase, the western Pacific anomalies are related to negative correlations with precipitation, while the eastern Pacific anomalies exhibit positive correlations. Larger communities in the SST field are spatially reorganized and preferentially located in the western Pacific, while associated precipitation components commonly link to more than one SST pattern.
By studying the n.s.i. cross-degree not only for absolute correlations, but also separately for strong positive and negative correlations, we have demonstrated that it is possible to distinguish different patterns of variability from each other. Moreover, two-variable community analysis can be employed to classify the nodes of our coupled SST-PCP network representation into distinct groups of strongly co-varying grid point time series. Accordingly, the proposed two-variable community analysis provides additional information as compared with using only the n.s.i. cross-degree, which can be exploited to further characterize the main regions in both climate fields that mutually influence and interact with each other.
Interestingly, SST → PCP patterns that are markedly visible in the neutral ENSO phase are weakened during La Niña and even suppressed during El Niño. This particularly applies to the Indian ocean (with variability largely characterized by the Madden-Julian Oscillation and Indian Ocean Dipole [61,70]), the Mediterranean Sea to Caspian Sea, some region near Newfoundland, and the Caribbean Sea (reduced during El Niño, suppressed during La Niña). In turn, during El Niño phases the natural links between SST and precipitation are interrupted and dominated by the SST variability of the tropical Pacific (and southern Pacific for La Niña).
Our analysis has also recovered known teleconnections between SST and PCP triggered by El Niño, which, however, mostly correspond to wet anomalies over the Philippines, Uganda and near the US-Mexican border. Other regions previously identified to often exhibit dry conditions during El Niño phases (including Indonesia, northern Australia and the Amazonas region) do however not match well with the patterns obtained by our analysis. One reason for this could be the strong links from the circum-Antarctic ocean region that might hide other (weaker) statistical links. During La Niña, nearly all regions that are known to be markedly influenced by ENSO teleconnections exhibit elevated cross-degree, including wet regions in north-eastern South America and over the maritime continent as opposed to dry areas over Taiwan and Uganda.
From the methodological perspective, we have restricted ourselves to a specific setting involving linear correlations and a particular choice of thresholds. According to our experiences from other studies, we expect that the main patterns discussed above should not be altered qualitatively when employing different thresholds, while a detailed investigation of the associated impacts would be beyond the scope of the present study. Moreover, we emphasize that we have chosen to employ monthly data and focus on lag-zero (i.e., quasi-instantaneous) interrelationships between SST and PCP. Another reasonable extension of the presented analysis would therefore be allowing for timelagged correlations. While the choice of monthly scales in this work has put a clear focus on any kind of SST-precipitation relationship explained by persistent synoptic weather regimes, long planetary wave related connections or processes mediated by secondary actors may manifest at timescales considerably larger than one month. Identifying such processes along with associated spatial patterns would require a systematic variation of the mutual time lags between SST and precipitation field, which has been beyond the scope of our present study.
A natural extension of the present study to be considered in future work consists of a distinction between short-range (local) SST-precipitation connections and long-range teleconnectivity (due to phenomena like planetary waves and atmospheric moisture recycling). One straightforward way would be including information on the average spatial distance captured by the (cross-) links emanating from each node into the analysis. It can be expected that the thus obtained additional information will provide deeper insights into the coexistence of many strong positive and negative correlations in certain parts of the SST field, which clearly point to an effect of teleconnections either fostering or suppressing rainfall in different regions.
From the climatological perspective, it could furthermore be useful to employ a different classification of ENSO phases according to the state of eastern tropical Pacific SST anomalies during the corresponding peak season (boreal winter) only. In such case, it would be reasonable to not only study the period from the boreal fall before to the boreal summer after an El Niño or La Niña event. By contrast, it would be also relevant to consider the boreal summer season preceding such events. Unfortunately, such an analysis would be complicated by the existence of paired (multi-year) ENSO phases, as well as close successions of El Niño and La Niña phases in subsequent years. Moreover, another relevant yet highly intermittent actor in the SSTprecipitation interplay is strong volcanism, which has recently been shown to potentially foster El Niño conditions [71] and also strengthen the synchrony between ENSO variability and interannual changes of the Indian summer monsoon [72,73]. While sufficiently strong volcanic eruptions are, however, rare and may therefore not have affected the results of our presented analysis, further investigating them in terms of long equilibrium runs or ensemble simulations of state of the art general circulation models remains an interesting avenue for follow-up studies.

Summary
We have used the coupled climate network approach to systematically study global ocean-atmosphere interactions. Specifically, we have investigated how the covariability patterns between sea surface temperature (SST, ocean) and precipitation (atmosphere) change at monthly timescales during different seasons and phases of the El Niño-Southern Oscillation (ENSO). Unlike complementary multi-layer network approaches that can be used for inter-comparing spatio-temporal correlation structures of different climate fields that normally have to be provided at the same spatial grid, the coupled network approach takes up information on covariability both within and among two (or more) spatial fields of time series, thereby integrating unipartite (within-field) and bipartite (between-fields) correlation information and avoiding the restriction to equal grids [35,37].
In the present work, we have focused on two relatively basic structural network characteristics, the n.s.i. crossdegree field (measuring bipartite connectivity between SST and precipitation and thereby including both, local and remote connections) and the two-variable community structure (combining within and between field correlations) between SST and precipitation for absolute, as well as for positive and negative Pearson correlations. By studying the n.s.i. cross-degree, we have demonstrated that especially the tropical oceans affect global precipitation in a spatially coherent way, yet differently during different phases of ENSO. By comparing positive and negative values of correlations, we have found that during the neutral phase of ENSO, the western (west of 180 • E) Pacific, Indian ocean and Atlantic SSTs (in boreal fall) have negative effects on precipitation and vice versa, whereas the eastern regions like the centralto-eastern equatorial Pacific rather experience a positive SST effect on precipitation. Similar observations can be made during La Niña phases. As expected, the cross-degree patterns during El Niño phases reveal positive correlations between SST and precipitation especially in the eastern equatorial Pacific. Only during boreal fall, western equatorial Pacific SSTs also exhibit pronounced negative correlations with global precipitation patterns.
Beyond the cross-degree fields, two-variable communities point to strongly interacting regions between SST and precipitation, ignoring how each individual node interacts with the whole climate system. Communities highlight the main areas in which SST and precipitation are strongly connected among each other. Accordingly, during neutral ENSO phases the main interaction zone between SST and precipitation is located in the central equatorial Pacific. In turn, during La Niña phases, the main community in the SST fields covers much of the central-to-eastern equatorial Pacific, while the corresponding precipitation region mostly coincides with the central equatorial Pacific. These patterns are thus in general agreement with our findings for the n.s.i. crossdegree. During El Niño phases, two separate patches in the eastern and central equatorial Pacific are identified as the most relevant regions in terms of SSTprecipitation co-variability.
Our analysis has demonstrated how the interaction between two spatial fields of different climate variables can be studied with complex network tools during different ENSO types and seasons. Besides further integration of more sophisticated coupled network characteristics, including such taking up spatial information on the linkage patterns as well, the present work opens up several lines for follow-up research. On the one hand, one possible subject of future research could be extending the present approach to more than two variables, or considering two variables but focusing on two or more distinct key regions to capture their specific interactions. On the other hand, for studying climate dynamics one should also consider the effect of different timescales, i.e., how the interactions among different climate variables evolve across time and scale and how this is related with the network topology [74][75][76], including the community structure. For the latter purpose, one could use different methods (e.g., wavelet analysis) to decompose the time series into different timescales and study the interaction between two or more climate variables at separate timescales or between different scales to capture cross-scale variability [77].
Beyond its specific climatological focus, our study provides another example for the potential usefulness for diagnosing and attributing spatio-temporal patterns of climate variability across variables and regions. Future studies shall further take up these potentials of the employed methodology and investigate to which extent the gained knowledge can be further utilized for improving climate model diagnostics and statistical climate forecasts. 1306A), the Belmont Forum/JPI Climate project GOTHAM (Globally Observed Teleconnections and Their Representation in Hierarchies of Atmospheric Models, Grant No. 01LP 16MA) and the JPI Climate/JPI Oceans project ROADMAP (The Role of Ocean Dynamics and ocean-atmosphere Interactions in Driving ClimAte Variations and Future Projections of Impact-Relevant Extreme Events, Grant No. 01LP2002B). The authors gratefully acknowledge the European Regional Development Fund (ERDF), the German Federal Ministry of Education and Research and the Land Brandenburg for supporting this project by providing resources on the high performance computer system at the Potsdam Institute for Climate Impact Research. All computations have been performed using the Python package pyunicorn [78].
Funding Open Access funding enabled and organized by Projekt DEAL.

Author contribution statement
NE, CC and RVD designed the study. NE performed the numerical analyses with the help of CC. CK and RVD interpreted the results. NE, CK and RVD drafted the manuscript. All authors reviewed the manuscript and approved it for submission.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecomm ons.org/licenses/by/4.0/.