## Abstract

Due to the physical coupling between atmosphere and ocean, information about the ocean helps to better predict the future of the atmosphere, and in turn, information about the atmosphere helps to better predict the ocean. Here, we investigate the spatial and temporal nature of this predictability: where, for how long, and at what frequencies does the ocean significantly improve prediction of the atmosphere, and vice versa? We apply Granger causality, a statistical test to measure whether a variable improves prediction of another, to local time series of sea surface temperature (SST) and low-level atmospheric variables. We calculate the detailed spatial structure of the atmosphere-to-ocean and ocean-to-atmosphere predictability. We find that the atmosphere improves prediction of the ocean most in the extratropics, especially in regions of large SST gradients. This atmosphere-to-ocean predictability is weaker but longer-lived in the tropics, where it can last for several months in some regions. On the other hand, the ocean improves prediction of the atmosphere most significantly in the tropics, where this predictability lasts for months to over a year. However, we find a robust signature of the ocean on the atmosphere almost everywhere in the extratropics, an influence that has been difficult to demonstrate with model studies. We find that both the atmosphere-to-ocean and ocean-to-atmosphere predictability are maximal at low frequencies, and both are larger in the summer hemisphere. The patterns we observe generally agree with dynamical understanding and the results of the Kalnay dynamical rule, which diagnoses the direction of forcing between the atmosphere and ocean by considering the local phase relationship between simultaneous sea surface temperature and vorticity anomaly signals. We discuss applications to coupled data assimilation.

## 1. Introduction

The ocean and atmosphere are coupled by numerous processes involving exchange of mass, heat, and momentum (Frankignoul 1985; Kushnir et al. 2002; Sobel 2007). This implies that knowing the state of the atmosphere helps to better predict the future state of the ocean, and vice versa. This fact has motivated the improvement of ocean observing systems, as well as coupled atmosphere–ocean models and data assimilation (Sluka et al. 2016; Penny et al. 2019). In this study, we employ a model-independent approach to investigate the improvement in prediction of the atmosphere and ocean by including information from the other component. To later interpret the predictability results, we briefly review the literature on the dynamic coupling and predictability of the atmosphere and ocean.

Climate and its variability are strongly influenced by the ocean, and in particular by SSTs, which also play a key role as a source of potential predictability for climate fluctuations. The large-scale structure of SST anomalies depends not only on large-scale atmospheric circulation and its ensuing heat fluxes but also on heat transport by currents and vertical mixing (Ekman currents and pumping). Ekman pumping is especially energetic at subsynoptic scales (Frankignoul 1985; Deser et al. 2010). The coupling between SST anomalies and the overlying atmospheric circulation varies geographically. It is known that in the extratropics, it is primarily the atmosphere that drives SST rather than vice versa (Frankignoul 1985). This atmosphere-forced variability is an important source of low-frequency variability in the climate system (Hasselmann 1976; Frankignoul and Hasselmann 1977). In modeling studies, Luksch and von Storch (1992), Luksch (1996) found that much of the low-frequency SST variability in the North Pacific and the North Atlantic could be explained by wind anomalies, mainly through anomalous heat fluxes and Ekman transport. Several studies have examined how predictable SST is from atmospheric forcing, including Scott (2003) with an idealized stochastic model. Model results have also shown that tropical SSTs are highly predictable when atmospheric fluxes are prescribed (Shukla and Kinter 2006).

On the other hand, SSTs substantially drive low-level atmospheric flow in the tropics. Two major classes of models exist for low-level tropical flow (Sobel 2007). Matsuno–Webster–Gill models assert that the heating due to deep convection drives the surface winds (Gill 1980). Lindzen–Nigam models, on the other hand, assert that the lower troposphere in the tropics is well mixed due to convection below the trade inversion, such that low-level atmospheric temperature gradients strongly resemble the temperature gradients at the ocean surface (Lindzen and Nigam 1987). Contemporary evidence suggests that zonal surface winds are well explained by convective heating, but that in regions of strong meridional SST gradients the Lindzen–Nigam model is successful (Sobel 2007).

SST gradients and their associated turbulent fluxes also affect the pattern of moisture convergence, and thus the pattern of diabatic heating anomalies, which can have significant effects on large-scale circulation. This effect is more pronounced in the tropics, because divergent flow grows larger compared to rotational flow closer to the equator (Shukla and Kinter 2006). The fact that baroclinic instability is less significant in the tropics implies that longer-range predictability can be obtained from boundary conditions such as SST (Charney and Shukla 1981). Shukla (1998) showed that the tropical atmospheric flow and rainfall is so strongly determined by SST that it can be forecast as long as SST can be forecast, unlike the extratropical atmosphere. In the extratropics, the influence of SST anomalies on the atmosphere is more difficult to identify, but a large number of studies have shown that there are small, yet discernible, effects (Kushnir et al. 2002). A reason for the weaker extratropical response is that in the tropics vertical advection dominates in the atmospheric response to heating, while horizontal advection dominates in the extratropics (Hoskins and Karoly 1981); thus it is easier for heating to influence the free atmosphere in the tropics (Thomson and Vallis 2018a).

Additional studies of the spatial structure of atmosphere–ocean coupling have been done using a dynamical rule, first proposed by Kalnay et al. (1986). This rule determines whether the atmosphere or ocean is the dominant local driver based on the phase relationship between SST and 850-hPa vorticity anomalies (see Fig. 1). Using this rule, Peña et al. (2003) showed that the ocean generally drives the atmosphere in the tropics, and the atmosphere drives the ocean in the extratropics. By testing the Kalnay rule on different reanalyses and CMIP5 models, Ruiz-Barradas et al. (2017) found that it was robust and that it can be used for identifying model deficiencies. Other studies have used the relationship between SST and rainfall to determine the local driver (Wu and Kirtman 2007; Kumar et al. 2013a). Kumar et al. (2013a) found that the SST drives the atmosphere most strongly in the tropical eastern Pacific.

There is substantial regional variation in predictability, especially between the tropics and extratropics, but also on smaller scales. Analyses that group different regions together can thus miss the detailed spatial structure of predictability (DelSole and Tippett 2007). In this paper, we study the predictability of SST from the low-level atmosphere, and vice versa. We investigate the spatial variation of predictability, maximum lead times, seasonality, and frequency decomposition, and provide dynamical interpretations. We believe this to be the first work to provide a global picture of the predictability of the atmosphere from the ocean and the predictability of the ocean from the atmosphere. We also use it to determine where the atmosphere and ocean are driving and compare it to the Kalnay dynamical rule. The dynamical rule is based only on Ekman pumping and the direct, linear (Kushnir et al. 2002) impact of ocean temperature anomalies on the atmosphere. Here, we study the predictability of SST and the atmosphere using a more general statistical method, which allows the identification of coupling through other mechanisms as well. Furthermore, this work has possible applications to coupled data assimilation, which we will discuss.

## 2. Methods and data

### a. Granger causality

Granger causality is based on modeling two sets of time series as stochastic autoregressive processes, and quantifying the improvement of prediction skill beyond what could be obtained from only the historical record of the one signal by including the information from the other. Granger causality has been applied to some climate phenomena, serving as a more robust alternative to more commonly used lagged correlation analyses (McGraw and Barnes 2018).

According to the Wiener–Granger definition of causality, given variables *X* and *Y*, “*X* causes *Y*” if, in an appropriate statistical sense, *X* assists in predicting the future of *Y* beyond the degree to which *Y* already predicts its own future (Wiener 1956; Barnett and Seth 2014). Granger (1969) formalized this definition in terms of linear autoregression.

A *p*-order vector autoregressive model, VAR(*p*), for an *n*-dimensional stationary stochastic process **U**, sampled at discrete time indices *t*, is defined as

where the *n* × *n* matrices $Ai$ are the optimal regression coefficients, and $\epsilon t$ is the independent and identically distributed (iid) residual noise vector for time *t*. Assuming that **U** is split into two jointly distributed multivariate processes **X**_{t} and **Y**_{t}, the VAR(*p*) model (1) can be written as

The **Y** component from the full regression model (2) is

Here, the state of **Y** depends on the past data of both itself and **X**, provided $Ayx$ ≠ 0 and $Ayy$ ≠ 0.

One can investigate the influence of **X** on **Y** by using a reduced model where **X** is removed from the set of information and **Y** is solely predicted by its own past,

The terms $Ayy,i\u2032$ and $\epsilon y\u2032$ are the reduced (optimal) regression coefficients and reduced residual, respectively. Recalling the Wiener–Granger definition, if the full regression model (3) yields a significantly better prediction than the reduced regression model (4), then the null hypothesis (i.e., no Granger causality) is rejected and **X** is identified as a Granger cause of **Y**.

One way to quantify the changes in prediction error of **Y** between the full (3) and reduced (4) VAR models is the log-likelihood ratio

where $\Sigma yy=cov\u2061(\epsilon y,t)$ and $\Sigma yy\u2032=cov\u2061(\epsilon y,t\u2032)$ (Barnett and Seth 2014). In the case that **Y** is one-dimensional, $\Sigma yy$ and $\Sigma yy\u2032$ are the mean-square (prediction) errors of their respective VAR models, a measure of forecast skill (Barnett and Seth 2017). In the case that **Y** is multidimensional, the determinants of the covariance matrices are used, called the “generalized variance” (Barrett et al. 2010), and $F$ can be expressed as a sum over combinations of the univariate “cause” and “effect” variables (Barrett et al. 2010). Note that $F$ can be tested for statistical significance using the *F* or *χ*^{2} tests (Barnett and Seth 2014).

Similarly, one can investigate the other direction, $FY\u2192X$, for a Granger causal effect of **Y** on **X**; unlike correlations, $FY\u2192X\u2260FX\u2192Y$ in general. Note that $F$ is always nonnegative, since the inclusion of additional explanatory variables in the full linear model can only decrease the prediction error. A larger value of $F$ indicates a larger relative improvement of prediction skill of **Y** by including **X**. A comparison of the magnitude of different $F$ values is meaningful because of their asymptotic equivalence to transfer entropy for a large class of processes (Barnett and Seth 2014; Barnett et al. 2009).

In this paper, we will use Granger causality in order to determine whether including information from the atmosphere or ocean significantly improves prediction of the other component. We will refer to the Granger causality from the atmosphere to the ocean, $FAtmos\u2192SST$, as “atmosphere-to-ocean predictability” (and vice versa). Although we cannot rigorously infer physical causality from statistical predictability, the predictability in the system does have physical origins, which we endeavor to explain in section 3.

### b. Lead times

To determine the dependence of Granger causality on lead time, we consider shifted signals. By backshifting the predictor (“cause”) signal by a time *τ* > 0, the fitting procedure for the full model does not “know” about the *τ* most recent data points of the predictor available at any given time. This corresponds to fitting a modified version of (3), where the indices of the predictor variable are shifted:

If introducing this shift of *τ* days diminishes the improvement in prediction due to a certain variable to the point where it loses statistical significance, we say that this variable does not improve prediction at a lead time of *τ* + 1 days (the extra 1 day appears because we are considering one-step-ahead prediction). If there is still significance, we can infer that the predictor variable provides improved prediction of the predictand (“effect”) *τ* + 1 days in advance.

We will present results for the maximum lead time at which atmosphere-to-ocean and ocean-to-atmosphere predictability is significant in each grid cell. This is determined by evaluating whether the predictability is significant at a sequence of lead times, and then taking the maximum lead time for which all the previous lead times are also significant. For example, if the predictability is significant at a lead time of 15 days and at 30 days but not at 45 days, we take 30 days as the maximum lead time.

### c. Frequency decomposition

Granger causality can also be decomposed by frequency *ω*, yielding a function *f*_{X→Y}(*ω*). The latter can be interpreted as the proportion of spectral power of **Y** at frequency *ω* that can be attributed to interaction with **X** (Barrett and Barnett 2013). We can then consider the Granger causality limited to the band of frequencies between *ω*_{0} and *ω*_{1}:

When *ω*_{0} = 0 and *ω*_{1} is the Nyquist frequency (half of the sampling rate of the signal), this is equal to the regular Granger causality $FX\u2192Y$. Thus, the time-domain Granger causality can be considered an average over all frequencies of the frequency-domain Granger causality (Barnett and Seth 2014).

### d. Implementation

In this study, we consider the Granger causality between SST and a five-dimensional signal composed of SST anomalies and five atmospheric variables that characterize the low-level atmosphere (Atmos): surface pressure anomalies and vorticity, divergence, air temperature, and specific humidity anomalies at 850 hPa. The five atmospheric time series are treated as a single signal representing the atmospheric state; that is, we look at whether the SST improves prediction of the multidimensional atmospheric state, and whether the multidimensional atmospheric state improves prediction of SST. We choose 850 hPa since it is in the free atmosphere (above the boundary layer), more likely to be useful for prediction (Thomson and Vallis 2018a). We calculate the Granger causality measure in both directions, (Atmos→SST) and (SST→Atmos), at every grid cell of the global oceans. We also calculate the statistical significance at all grid cells. If the statistical significance level at a grid cell is less than the chosen threshold, we cannot reject the null hypothesis of no improvement in prediction. In the subsequent figures, regions where the null hypothesis cannot be rejected are shown in white. When multiple significance tests are done at different lead times, in order to account for the problem of multiple testing we use the Benjamini–Hochberg procedure (Benjamini and Hochberg 1995). We use the Multivariate Granger Causality (MVGC) Toolbox for Matlab to carry out the Granger causality analysis (Barnett and Seth 2014). We provide our results in NetCDF format, as well as the open-source code for our analysis and associated plots.^{2}

We assign each grid cell a lag order *p* for computing the Granger causality in that cell. We do this by fitting the full VAR model with lag orders from *p* = 1 day to *p* = 45 days, and select the order that minimizes the Akaike information criterion (AIC) (Akaike 1974). The most common *p* selected by AIC is 4 days, and tends to be higher in the tropical oceans (typically between 4 and 8 days). We keep these same lag orders for each cell when considering the dependence of predictability on lead times in section 2b. AIC is commonly used in selecting the order of AR processes (von Storch and Zwiers 2002), achieving a compromise between underfitting the model by using too low an order and overfitting by using too high an order. The AIC can be understood as a comparative measure of the model’s prediction error, compensating for the fact that the training error will be a biased estimate of the prediction error (Efron and Hastie 2016). However, the predictability estimates are not very sensitive to the lag order selection; we have tested selection of lag orders with the Bayesian information criterion (which generally picks more parsimonious models than AIC) as well as a fixed lag order *p* = 5. Both of these give very similar results of Granger causality.

For the analysis by season, the same season in different years is assumed to be an independent realization of the same VAR process. Each VAR model is thus fit once based on the set of time series of that season in the years available.

### e. Data

We perform the analysis at daily average temporal resolutions over the oceans. The SST and atmospheric fields are obtained from the ERA-Interim reanalysis (Dee et al. 2011) for 1979–2017. Due to the spatial scale dependence of predictability (see section 2f), we use the reduced N128 Gaussian grid. The Gaussian grid has cells with almost equal areas, about 80 km × 80 km, as opposed to a regular latitude–longitude grid whose cell areas vary with a factor of cos(latitude) due to the spherical geometry. We use the full length of the dataset, meaning that the calculated predictability represents the overall behavior of the coupled atmosphere–ocean system over the entire time period. Note that although the reanalysis of the atmospheric fields is carried out in an uncoupled mode in ERA-Interim, the coupled dynamics of the real-world system are reflected in the time series through the atmospheric data assimilation and the use of observed SSTs (Kumar et al. 2013b; Ruiz-Barradas et al. 2017). Penny et al. (2019) have shown explicitly that data assimilation can cause an uncoupled atmosphere to synchronize to coupled dynamics.

Note that with smoothing methods (such as 4D-Var, used in ERA-Interim) information flows from the future to the past within the assimilation window (Carrassi et al. 2018), so it is important to ensure that the assimilation window does not overlap over multiple points in the time series used for the Granger analysis. Here, we use 24-h temporal resolution while ERA-Interim’s assimilation window is 12 h (starting at 0000 and 1200 UTC), and we average from 0000 UTC to the next 0000 UTC; thus, there is no overlap.

The anomaly time series are obtained by subtracting the first two harmonics (annual and semiannual cycles) from the daily time series as computed in each grid cell, as in Peña et al. (2004) and Ruiz-Barradas et al. (2017). We use the anomaly time series since otherwise there is a trivial component of predictability due to the regularity of the seasonal and subseasonal cycles.

### f. Methodological notes

We briefly address the use of Granger causality to measure predictability of the atmosphere and ocean. Linear autoregressive processes, especially AR(1) and AR(2), are widely used to model climatic time series (von Storch and Zwiers 2002). In particular, it has been suggested that many features of the atmosphere–ocean system at time scales from months to several years can be understood by modeling atmospheric time series as white noise and SST as an AR(1) process (Frankignoul and Hasselmann 1977). Moreover, a fairly general class of discrete-time stochastic processes, including nonlinear ones, can be approximated as VAR processes with sufficiently high lag order, even if this is not the most parsimonious model (Poskitt 2007; Barnett and Seth 2014, 2017). Additionally, Granger causality is theoretically invariant under a wide class of invertible digital filters (Barnett and Seth 2011), an important property for predictability measures (DelSole and Tippett 2007). The fact that this is a data-driven method, independent of model representations of physical processes except insofar as they are reflected in the reanalysis, is an advantage over GCM-based predictability studies (National Academies of Sciences, Engineering, and Medicine 2016).

Predictability depends on spatial scale, as larger structures tend to be more persistent than smaller ones (DelSole and Tippett 2007). For example, longer waves in the atmosphere are more predictable than shorter ones (Shukla 1985; Dalcher and Kalnay 1987), and large-scale SST anomalies are generally more persistent than those at smaller scales (Frankignoul 1985). Because we calculate predictability locally in each grid cell, at higher grid resolutions the forecasts should lose skill sooner due to the time taken for wave propagation and the advection of structures. Recent studies have also shown the significance of atmosphere–ocean interactions due to mesoscale ocean eddies, in which forcing can sometimes have the opposite direction than that observed on larger scales (Bishop et al. 2017; Hewitt et al. 2017). We have repeated our calculations with the same ERA-Interim reanalysis but on the reduced N48 Gaussian grid, which has cells of about 210 km × 210 km, or about 7 times coarser than N128, which we use for the main analysis. The results for predictability and maximum lead times are almost identical, except for isolated regions that have slightly longer maximum lead times, as expected from the above discussion. Thus, even though there is spatial scale dependence of predictability, it is not very strong at this range of resolutions. Last, our method is inherently local, as the Granger test is carried out individually in every grid cell. It thus cannot account for interactions between remote regions, such as those between the tropics and extratropics, as well as other teleconnections.

Globally averaged over the full time series of daily anomalies, and with no lead time, 94% of the variance in the SST anomalies is explained by the full VAR model for SST: 1 − [Σ_{SST}/Var(SST)] ≃ 0.94, where Σ_{SST} is the mean-square error of the full VAR model for SST. A total of 93% of the generalized variance in the atmospheric variables is explained by the full VAR model for the atmosphere: 1 − [|**Σ**_{Atmos}|/Var(Atmos)] ≃ 0.93, where Var(Atmos) is the generalized variance of the atmospheric variables.

For interpretations in terms of physical causality, correlation analyses are more liable to capture spurious correlations than methods that take into account the autocorrelation of time series, such as the Granger method (Dean and Dunsmuir 2016; Cryer and Chan 2008; BozorgMagham et al. 2015); examples of the advantages of Granger causality over lagged regressions in the climate context are given in McGraw and Barnes (2018).

## 3. Results

### a. Atmosphere-to-ocean and ocean-to-atmosphere predictability

The atmosphere-to-ocean predictability $FAtmos\u2192SST$ is shown in Fig. 2a. We observe that the atmosphere improves prediction of the ocean almost everywhere. Furthermore, there is a pronounced pattern of increase of atmosphere-to-ocean predictability along regions of large SST gradients. As per Frankignoul (1985), the wind-driven current contribution to the rate of change of an SST anomaly *T*′ is

where ** τ**′ is the anomalous wind stress, $T\xaf$ and

*T*′ are the mean and anomalous SST,

*ρ*is the density,

*f*is the Coriolis parameter, and $h\xaf$ is the mean mixed layer depth. Thus, regions of higher climatological SST gradients (see Fig. S1 in the online supplemental material) should display more of an atmospheric influence, reflected in the predictability. This is indeed seen in the regions of high SST gradients associated with the warm Gulf Stream, Kuroshio, Agulhas/Agulhas Return, and Brazil Currents. Along the currents themselves, where the SST is substantially driven by advection, $FAtmos\u2192SST$ is relatively small. This may be the reason for the small atmosphere-to-ocean predictability in the Antarctic Circumpolar Current region as well. Other features, such as the fact that the highest atmosphere-to-ocean predictability is seen south of the Kuroshio, even though the region of highest SST gradients is to its north, may be due to atmospheric forcing of SST anomalies through anomalous heat fluxes. Frankignoul and Reynolds (1983) found that the heat flux term is larger than the Ekman transport term in the North Pacific, and the sum of these two terms is larger south of the Kuroshio. Bishop et al. (2017) found that SST variability is primarily ocean-driven in regions of high climatological SST gradients;

^{3}however, they only included the sensible heat flux contribution to the ocean mixed layer temperature, while we consider also other terms that can contribute to SST anomalies, such as Ekman transport.

At very high latitudes (in the Arctic Circle, and some parts of the Southern Ocean) both the atmosphere-to-ocean and ocean-to-atmosphere predictability lose significance. This may be due to the presence of sea ice, which complicates atmosphere–ocean interactions (Zhang et al. 2018). It may also be caused by the poor predictability in the polar regions due to the scarcity of in situ observations and inadequate understanding and modeling of polar processes (Spengler et al. 2016), which in turn is reflected in the reanalysis fields.

Figure 2b shows the ocean-to-atmosphere predictability $FSST\u2192Atmos$. This type of predictability, arising due to the boundary conditions of the atmosphere, is often termed “boundary-forced predictability” (Shukla 1985). The ocean improves prediction of the atmosphere most in the tropical Pacific, but also over nearly all of the extratropical ocean. The statistically significant results over almost all of the extratropics is notable, as such an effect is notoriously difficult to isolate using GCM studies; in these studies, even when unrealistically large SST anomalies are imposed, the signal-to-noise ratio is often too low to make out a response (Thomson and Vallis 2018a). Our results are consistent with the results summarized in Kushnir et al. (2002), which conclude that SST anomalies do have a small effect on atmospheric circulation in the extratropics.

Figure 3 shows the zonal median Granger causality. Comparing the magnitudes of $FAtmos\u2192SST$ to $FSST\u2192Atmos$, we see that the atmosphere improves prediction of the ocean more than vice versa, except in a small band of latitudes around the equator. The large degree of hemispheric symmetry suggests that the continents are not critical in determining the large-scale patterns of predictability, except for the location of SST gradients, which appears to influence the atmosphere-to-ocean predictability, as described above.

### b. Dependence of predictability on lead times

Next, we investigate the dependence of atmosphere-to-ocean and ocean-to-atmosphere predictability on lead times (see section 2b). A time limit for predictability is inevitable in a chaotic system due to the sensitivity on initial conditions. The extratropical atmosphere is particularly chaotic, while some regions of the tropical atmosphere are so strongly dependent on SST that they do not behave chaotically given the SST as boundary condition (Shukla 1998).

Figure 4a shows the maximum lead time for which there is significant atmosphere-to-ocean predictability. The predictability is short-lived in most regions of the extratropical oceans, generally lasting fewer than 16 days. In the tropical oceans (between about 20°N and 20°S) there are regions with longer-lived atmosphere-to-ocean predictability of a few months. Thus, although the magnitude of the atmosphere-to-ocean predictability in the tropics is relatively small compared to the extratropics (see Fig. 2a), it is longer-lived. This is consistent with the finding that tropical SSTs are highly predictable when atmospheric fluxes are prescribed (Shukla and Kinter 2006).

Figure 4b shows the maximum lead time for ocean-to-atmosphere predictability. There is long-lived predictability in the tropical oceans, lasting several months almost everywhere between about 30°N and 30°S. In the tropical Pacific and the western tropical Atlantic there are regions of predictability longer than a year. The longer predictability in these regions is likely due to the longer decorrelation times of SST anomalies in these regions (see Fig. S2). This is consistent with Peña et al. (2003), who found a long-lived correlation between leading SST and vorticity in the tropical Pacific. It is also consistent with Shukla (1998), who found that the tropical atmosphere is highly predictable from SST, especially over the Pacific. Furthermore, due to the slower time scale of the ocean compared to the atmosphere, SST signals have higher autocorrelation (Frankignoul and Hasselmann 1977) and we should expect SST-driven regions to have longer predictability from SST than atmosphere-driven regions from atmospheric variables. There are, in addition, large parts of the subtropics that exhibit ocean-to-atmosphere predictability longer than two weeks, which could have implications for subseasonal-to-seasonal prediction (National Academies of Sciences, Engineering, and Medicine 2016).

Figure 5 shows the global mean atmosphere-to-ocean and ocean-to-atmosphere predictability at lead times up to one month, demonstrating that both generally decay monotonically with increasing lead time. Initially, the atmosphere provides more predictability to the ocean than vice versa, but this predictability also decays faster. By five days the ocean-to-atmosphere predictability becomes dominant, largely owing to the long-term predictability in the tropics.

### c. Spectral analysis

Figure 6 shows the frequency decomposition of the Granger causality (see section 2c). The atmosphere-to-ocean predictability is almost white at periods longer than a few months, and smaller for shorter periods. The ocean-to-atmosphere predictability is stronger at longer periods. Thus, the most predictable frequencies of the atmosphere from the ocean are low frequencies, and the most predictable frequencies of the ocean from the atmosphere are also low frequencies. To quantify the “average” frequency for the predictability, we define the mean frequency weighted by the spectral Granger causality:

The inverse of the weighted-mean frequency of the atmosphere-to-ocean predictability averaged over all grid cells is 8.6 days, while the inverse of the weighted-mean frequency of ocean-to-atmosphere predictability is 13.6 days. Figure S3 shows the spatial variation of the weighted-mean frequency of predictability. The weighted-mean frequency of atmosphere-to-ocean predictability is higher in the tropics (particularly in the tropical Pacific and Atlantic) and in the Southern Ocean, and lowest in the midlatitudes. On the other hand, the weighted-mean frequency of ocean-to-atmosphere predictability is lowest in the tropics and generally increases with latitude from the equator, except in the tropical Atlantic.

The low-frequency maximum of the atmosphere-to-ocean predictability is consistent with Hasselmann (1976) and Frankignoul and Hasselmann (1977), who found that a slow climate variable (such as SST) stochastically driven by a fast weather variable (such as the atmospheric variables) is most affected by the low-frequency part of the weather spectrum. It is also consistent with Ruiz-Barradas et al. (2017), who found that there is a greater proportion of coupled anomalies of vorticity and SST (indicating atmosphere–ocean interactions as per the dynamical rule; see Fig. 1) at lower temporal resolution. The near whiteness of the atmosphere-to-ocean forcing low frequencies has also been noted in previous studies (Frankignoul 1985). The strong ocean-to-atmosphere predictability at low frequencies can be understood from the ocean’s slower time scale and confirms the ocean’s role as a source of predictability of low-frequency variability of the atmosphere.

### d. Seasonality

We see a strong seasonality in the predictability. Both atmosphere-to-ocean and ocean-to-atmosphere predictability are stronger in the summer hemisphere, although the seasonality is more striking in atmosphere-to-ocean predictability. Figure 7 shows the summer and winter atmosphere-to-ocean and ocean-to-atmosphere predictability.

SST is more predictable from its own past in the winter than in the summer, with the reduced model having a root mean-square error of 0.13 K in the winter hemisphere and 0.17 K in the summer hemisphere with no lead time. This is likely because the decay of SST anomalies due to negative heat flux feedbacks takes longer in the winter due to a deeper mixed layer (Frankignoul 1985; Park et al. 2005). Thus, there is more “room to improve” by including the atmosphere in predicting the SST in the summer, which explains the higher atmosphere-to-ocean predictability. Similarly, the atmosphere is also more predictable in the winter (Shukla 1985), providing a possible explanation of the seasonality of the ocean-to-atmosphere predictability. We see that the ocean-to-atmosphere predictability in the tropical Pacific does not vary significantly by season, consistent with Thomson and Vallis (2018b). There are features of high ocean-to-atmosphere predictability around 20°S that are present only in DJF.

Note that in the seasonal analysis some regions lose statistical significance of ocean-to-atmosphere predictability; this is due to having fewer data (25% as many data points) to establish significance. This reinforces the weakness of the ocean-to-atmosphere forcing in these regions.

### e. Local driver

In each grid cell we compute the logarithm of the ratio of the Granger causalities:

A positive value indicates that the atmosphere is predominantly driving the ocean, and vice versa.

We compute (10) at daily resolution, and also for frequencies lower than 1 month^{−1} [i.e., with *ω*_{0} = 0, *ω*_{1} = 1 month^{−1} in (7)]. The statistical significance is not considered since the frequency-limited Granger causality cannot be tested with the standard significance tests. Figure 8a shows the results for daily resolution, showing that the atmosphere is the local driver over almost all of the extratropics, except in currents where the SST is mostly advection driven, as described in section 3a. The ocean is the main driver in the eastern tropical Pacific, and in smaller regions in the Indian and Atlantic Oceans. Figure 8b shows the results for frequencies lower than a month, showing that the ocean-driven regions greatly expand in the tropical oceans (and also in the Southern Ocean, but this is uncertain due to the weak significance here, as discussed in section 3a).

We now compare the results of the Granger analysis to those of Ruiz-Barradas et al. (2017), which provides an indication of whether the ocean or atmosphere is driving using the Kalnay dynamical rule (Fig. 1). Since the dynamical rule uses only the vorticity, to compare the two methods we use only vorticity at 850 hPa with the Granger method, instead of the set of five atmospheric variables used in Fig. 8 and in the rest of the paper. For the dynamical rule, we compute the log of the ratio of the anomaly counts where atmospheric forcing was dominant to that where oceanic forcing was dominant [see Ruiz-Barradas et al. (2017) for details].

Figure 9 displays comparisons between the two results at 5-day resolution. We do not consider statistical significance for the dynamical rule, since it is not clear how to determine it. There is general agreement that atmosphere-to-ocean predictability/forcing is dominant in most of the extratropics, while ocean-to-atmosphere predictability/forcing is dominant in the tropics between about 30°N and 30°S, and especially in the Pacific. The direction of the predominant forcing agrees in 77% of the grid cells.

## 4. Summary and discussion

We employ Granger causality analysis to determine the local atmosphere-to-ocean and ocean-to-atmosphere predictability between time series of SST and low-level atmospheric variables over the global oceans. We find that the atmosphere improves prediction of the ocean and the ocean improves prediction of the atmosphere in both the tropics and extratropics. Our finding of a statistically significant signature of the ocean on the atmosphere nearly everywhere in the extratropics is notable because it is difficult to demonstrate with GCMs. The atmosphere-to-ocean predictability is stronger in regions of high SST gradients, and lasts a few days in the extratropics, but up to several months in regions of the tropical Pacific and Indian Oceans. The ocean-to-atmosphere predictability is strongest in the tropical Pacific, and is long-lived across the tropical oceans, in many regions lasting longer than a year. Both the atmosphere-to-ocean and ocean-to-atmosphere predictability are larger at low frequencies, the latter more so. Both the atmosphere-to-ocean and ocean-to-atmosphere predictability are stronger in the summer hemisphere. We find that at daily resolution the ocean predominantly drives the atmosphere in the tropical Pacific, and the atmosphere is the primary driver in the extratropics. At frequencies lower than a month, the ocean is the driver across most of the tropical oceans. The patterns are broadly similar to the patterns of forcing determined by the Kalnay dynamical rule (Ruiz-Barradas et al. 2017).

The results highlight the regions where predictability can be obtained from the atmosphere and the ocean, which could be useful for identifying the regions where forecasts, especially at subseasonal-to-seasonal time scales, could be improved (National Academies of Sciences, Engineering, and Medicine 2016). This should also have applications for coupled atmosphere–ocean data assimilation, in particular for guiding variable localization. Variable localization is the problem of determining when to use the cross-covariance between variables in the atmosphere and ocean in the analysis. In some cases, including this cross-covariance may degrade the analysis due to spurious correlations. In other cases it should improve the analysis by employing more information about the relationship between the variables in the assimilation process, as in strongly coupled data assimilation (Sluka et al. 2016; Penny and Hamill 2017; Penny et al. 2019). It is thus important to localize variables based on their physical relevance or on their background error correlation (Yoshida and Kalnay 2018). Furthermore, the benefit (or detriment) of using the cross-covariance can be highly regionally dependent. Several studies have shown that assimilating the atmosphere into the ocean benefits the tropical oceans the most (Lu et al. 2015; Sluka 2018; Storto et al. 2018) and that assimilating the ocean into the atmosphere benefits the extratropics most (Sluka 2018). This may be due to the fact that in the weakly coupled case information is flowing from one component to the other through the dynamics (primarily from the atmosphere to the ocean in the extratropics, and vice versa), so that the additional information gained through the strongly coupled assimilation benefits primarily the other direction of information flow (Sluka 2018). However, other studies have shown the opposite or mixed results (Liu et al. 2013; Sluka et al. 2016), highlighting the need for more research.

Besides the applications to data assimilation, several other extensions of the current work could be undertaken. To more closely examine the processes that produce the atmosphere-to-ocean predictability, different terms of the SST anomaly equation as in Frankignoul (1985) could be regressed onto the atmosphere-to-ocean predictability. To isolate the effect of El Niño/La Niña, years in the different phases of ENSO could be analyzed separately, and similarly for other climate oscillations. Interactions between other parts of the Earth system could also be explored using this method, such as stratosphere–troposphere and land–atmosphere processes—for example, the relationship between surface temperature on land and rainfall (Trenberth and Shea 2005). A similar predictability analysis could be undertaken for the design of observing systems, as a model-independent and inexpensive alternative to observation system simulation experiments (OSSEs). Convergent cross mapping, which is designed for systems where separability of “cause” and “effect” variables is difficult (Sugihara et al. 2012; BozorgMagham et al. 2015), could also be used.

## Acknowledgments

We thank the anonymous reviewers of previous versions of this manuscript, whose comments have substantially improved the paper. We also thank Cesar B. Rocha, Jagadish Shukla, Travis Sluka, and Takuma Yoshida for helpful discussions. The late distinguished professor Fuqing Zhang provided unique insights. E.B. is supported by the University of Maryland Flagship Fellowship and Monsoon Mission II funding (Grant IITMMMIIUNIVMARYLANDUSA2018INT1) provided by the Ministry of Earth Science, Government of India. E.K. and S.M. acknowledge Lev Gandin funding (Grant 2956713) provided by G. Brin. AR-B gratefully acknowledges support from the U.S. National Science Foundation through Grant AGS1439940. The authors acknowledge the University of Maryland supercomputing resources (http://hpcc.umd.edu) made available for conducting the research reported in this paper.

## REFERENCES

*Monsoon Dynamics*, J. Lighthill and R. P. Pearce, Eds., Cambridge University Press, 99–110, https://doi.org/10.1017/CBO9780511897580.009.

*Time Series Analysis: With Applications in R.*2nd ed., Springer Texts in Statistics, Springer-Verlag, 491 pp., https://doi.org/10.1007/978-0-387-75959-3.

*Next Generation Earth System Prediction: Strategies for Subseasonal to Seasonal Forecasts*. The National Academies Press, 350 pp., https://doi.org/10.17226/21873.

*Issues in Atmospheric and Oceanic Modeling*, B. Saltzman, Ed.,

*Advances in Geophysics*, Vol. 28, Academic Press, 87–122, https://doi.org/10.1016/S0065-2687(08)60186-7.

*Predictability of Weather and Climate*, T. Palmer and R. Hagedorn, Eds., Cambridge University Press, 306–341, https://doi.org/10.1017/CBO9780511617652.013.

*The Global Circulation of the Atmosphere*, T. Schneider, and A. H. Sobel, Eds., Princeton University Press, 219–251.

*Modern Mathematics for the Engineer: First Series*, E. F. Beckenbach, Ed., McGraw-Hill, 165–190.

## Footnotes

Supplemental information related to this paper is available at the Journals Online website: https://doi.org/10.1175/JCLI-D-18-0817.s1.

© 2019 American Meteorological Society. For information regarding reuse of this content and general copyright information, consult the AMS Copyright Policy (www.ametsoc.org/PUBSReuseLicenses).

^{1}

Although it is the curl of the surface wind stress that induces Ekman suction and pumping, the anomalies of (surface) vorticity and curl of the surface wind stress are usually in the same direction.

^{3}

In that paper, the regions of highest gradients were identified with the location of the currents; however, the currents correspond to local maxima of SST, so along the currents themselves the gradients vanish.