Abstract

A dynamically consistent state estimate is used for the period 1992–2011 to describe the changes in oceanic temperatures and heat content, with an emphasis on determining the noise background in the abyssal (below 2000 m) depths. Interpretation requires close attention to the long memory of the deep ocean, implying that meteorological forcing of decades to thousands of years ago should still be producing trendlike changes in abyssal heat content. Much of the deep-ocean volume remained unobserved. At the present time, warming is seen in the deep western Atlantic and Southern Oceans, roughly consistent with those regions of the ocean expected to display the earliest responses to surface disturbances. Parts of the deeper ocean, below 3600 m, show cooling. Most of the variation in the abyssal Pacific Ocean is comparatively featureless, consistent with the slow, diffusive approach to a steady state expected there. In the global average, changes in heat content below 2000 m are roughly 10% of those inferred for the upper ocean over the 20-yr period. A useful global observing strategy for detecting future change has to be designed to account for the different time and spatial scales manifested in the observed changes. If the precision estimates of heat content change are independent of systematic errors, determining oceanic heat uptake values equivalent to 0.1 W m−2 is possibly attainable over future bidecadal periods.

1. Introduction

The major observational obstacle to understanding the role of the ocean in climate is the extreme brevity of the instrumental record in a system having some memory exceeding several thousands of years. Datasets depicting the global interior ocean state begin with high accuracy altimetry only in 1992. The Argo array became quasi global in the mid-2000s. Assuming that these technologies continue to be supported (by no means clear), the community will ultimately have comparatively long records at least of the phenomena visible in upper-ocean hydrographic profiles and sea surface elevation.

Even in this recent period, major spatial and temporal inhomogeneities exist in these and related data. The main purpose of this paper is to examine the nature of the thermal variability in the deep ocean (below about 2000 m). At the present time, the Argo array (Roemmich et al. 2009), supplemented by elephant seal data (Roquet et al. 2013), is confined to the upper 2000 m and with the bulk of the extant values above 1000 m. Altimetric data respond to motions over the entire water column, although the partitioning of the motions they represent remains the subject of considerable debate. Most of the available abyssal measurements are sparse deep CTD profiles (Fig. 1) from hydrographic programs, sometimes designed to depict special regions (e.g., the Kuroshio or the Nordic Seas).

Fig. 1.

Hydrographic data reaching to 2000 m (a) between 1992 and 2000 and (b) between 2001 and 2011. (c),(d) The corresponding distributions reaching at least to 3600 m in the same two intervals.

Fig. 1.

Hydrographic data reaching to 2000 m (a) between 1992 and 2000 and (b) between 2001 and 2011. (c),(d) The corresponding distributions reaching at least to 3600 m in the same two intervals.

An important wider issue, which instigated this study, is the nature of a practical set of future observations capable of providing a basis for the understanding of ongoing ocean changes. At the end, some comments will be made about this problem, drawing on the results of the present analysis.

Much of the recent literature focuses on the ability to detect past and ongoing trends in ocean temperatures and heat content. The reality and magnitude of such changes is not the central goal here; rather it is to characterize the extent to which more general variability can be detected using the much more dense observational system of the last 10–20 yr. On the other hand, some order of magnitude numerical values are helpful for context.

Consider, for example, that greenhouse gas warming of the ocean is widely believed to be of order 1 W m−2 (e.g., Hansen et al. 2005) or less.1 The volume of the ocean is about 1.3 × 1018 m3. Using a mean density of 1038 kg m−3, the total mass is about 1.34 × 1021 kg, and with a heat capacity of roughly 3.8 × 103 J kg−1 °C−1, the global heat capacity is approximately 5.4 × 1024 J °C−1. A heating rate of 1 W m−2, if maintained for 20 yr, produces an energy content change of about 2.2 × 1023 J for a change in global ocean mean temperature of about 0.04°C. If the heating were confined to the upper 700 m, then based on a mean ocean depth of about 3700 m the temperature change is increased to about 0.2°C, and if all were confined to the region below that depth, the temperature change would be about 0.05°C (see Table 1). Recent observationally based estimates (Church et al. 2011) produce estimates closer to 0.5 W m−2, exacerbating the detection problem. (That the atmospheric radiation budget includes such poorly determined elements as changes in aerosols and cloud distributions is a major impetus to determining actual ocean heat storage changes.) Alternatively, a 1 mm yr−1 thermally induced change in global-mean sea level, if sustained for 20 yr, is consistent with a full-ocean volume mean temperature change of about 0.03°C, although important spatial variations exist in the sea level response to a fixed temperature change.

Table 1.

Approximate oceanic temperature changes implied by a 1 W m−2 heating rate over different times and depths, as well as the temperature change equivalent of a 1 mm yr−1 global-mean sea level (GMSL) change.

Approximate oceanic temperature changes implied by a 1 W m−2 heating rate over different times and depths, as well as the temperature change equivalent of a 1 mm yr−1 global-mean sea level (GMSL) change.
Approximate oceanic temperature changes implied by a 1 W m−2 heating rate over different times and depths, as well as the temperature change equivalent of a 1 mm yr−1 global-mean sea level (GMSL) change.

Another important question, pursued elsewhere, is whether available observations alone are capable of determining mean ocean temperatures, and whether the related heat content changes with time to accuracies and precisions useful at these levels. Estimating the global average change is especially challenging and here is only a by-product.

Historically, deep hydrographic measurements (below a few hundred or perhaps 1000 m) have been both difficult and expensive to acquire (see Abraham et al. 2013). The consequence has been sampling by a few, rare (in a multidecadal or centennial context), fragmentary top-to-bottom hydrographic stations and sections. Systematic global surveys did not begin until the era of the World Ocean Circulation Experiment, circa 1990. Figure 1 displays all of the oceanic temperature data (all CTD values) below 2000 m and below 3600 m since 1992 and used here [taken from the World Ocean Database 2009 of the National Oceanic and Atmospheric Administration (NOAA)]. Elephant seal temperature data do exist below 2000 m but are rare and are not included. By some standards (e.g., paleoceanography; see Huybers and Wunsch 2010), an impressive amount of data does exist: an evaluation of their adequacy can only be made in the context of the signal-to-noise structure and magnitudes at depth. Determining time changes with these datasets involves segregating them by interval with a consequently great reduction in the numbers available in any year or multiple of a year. To convey some of the observational difficulties, Fig. 2 displays the space–time standard deviation as a function of depth (not area weighted) as well as the standard deviation of the annual cycle. Accurate removal of the annual cycle and the temporal mean from individual data points is a major problem in the upper ocean, but not discussed here.

Fig. 2.

The standard deviation of temperature in the grid cells, in space, and time over 20 yr in the state estimate domain as a function of depth (solid line). In the absence of the eddy field, this curve is a very optimistic basis for determining average temperatures. Not area weighted. The variance includes the spatial time-mean contribution and which strongly dominates. The ability to remove it accurately is an issue in computing time changes from direct point observations. The dashed line is the global-mean standard deviation of the annual component.

Fig. 2.

The standard deviation of temperature in the grid cells, in space, and time over 20 yr in the state estimate domain as a function of depth (solid line). In the absence of the eddy field, this curve is a very optimistic basis for determining average temperatures. Not area weighted. The variance includes the spatial time-mean contribution and which strongly dominates. The ability to remove it accurately is an issue in computing time changes from direct point observations. The dashed line is the global-mean standard deviation of the annual component.

Some basic elements of the sampling problem are compiled in Table 2. About 52% of the ocean lies below 2000 m and about 18% below 3600 m. By defining a volume as having been “probed” if at least one CTD station existed within a roughly 60 × 60 km2 box in the interval 1992–2011, a minimal measure of sampling can be obtained. Thus, about ⅓ (11% of the total volume) of water below 2000 m was sampled during that time. Of the 16% lying below 3600 m, about 17% was measured.

Table 2.

Middle column refers to the volume of the ocean (total, below 2000 m, and below 3600 m) and the volume probed at least once by CTD (i.e., at least one CTD station existed within a roughly 60 × 60 km2 box) between 1992 and 2011 for various conditions (km3). Right column is the percentage of total volume. SH refers to the Southern Hemisphere.

Middle column refers to the volume of the ocean (total, below 2000 m, and below 3600 m) and the volume probed at least once by CTD (i.e., at least one CTD station existed within a roughly 60 × 60 km2 box) between 1992 and 2011 for various conditions (km3). Right column is the percentage of total volume. SH refers to the Southern Hemisphere.
Middle column refers to the volume of the ocean (total, below 2000 m, and below 3600 m) and the volume probed at least once by CTD (i.e., at least one CTD station existed within a roughly 60 × 60 km2 box) between 1992 and 2011 for various conditions (km3). Right column is the percentage of total volume. SH refers to the Southern Hemisphere.

As a consequence of this undersampling, even with the improvements in the last 20 yr, many papers have been published that simply assume no significant changes have taken place in the deep ocean over the historical period. Shifts in the deep-ocean properties may indeed be so slight that their neglect in discussions of heat uptake and sea level change is justified. The history of exploration suggests, however, that blank places on the map have either been assumed to be without any interesting features and dropped from further discussion, or at the other extreme, filled with “dragons” invoked to explain strange reports.2 It is also physically possible that in a search for abyssal trends that the higher-frequency, higher-wavenumber noise is negligible compared to the signals. In that view, the existing reports of deep trends based upon hydrographic lines (Roemmich and Wunsch 1984; Bryden et al. 1996; Joyce et al. 1999; Purkey and Johnson 2010; and others) are adequate. Recently, Balmaseda et al. (2013) offered estimates of abyssal changes with claimed accuracies of order of 0.01 W m−2 (0.0004°C temperature change equivalent over 20 yr) below 700 m. If that accuracy has in fact been obtained, the sparse coverage, perhaps extended to the scope of the World Ocean Circulation Experiment (WOCE) hydrographic survey, repeated every few decades, would be sufficient.

Given the combination of the high societal stakes in the accurate estimation of global heating rates and sea level rise, and the fundamental science questions of the understanding of oceanic variability, direct confirmation or refutation of this sufficiency hypothesis is essential. If it proves false, discussion can take place concerning the design of an adequate system.

2. A framework: The state estimate

Apart from the large-scale hydrographic survey done as part of WOCE (see Talley 2007), most direct ocean measurements have been made at the sea surface (altimetry, sea surface temperature, drifters) or obtained from XBTs (some reaching to O(750) m) and more recently from Argo floats, profiling primarily to 1000 m and more recently to 2000 m [e.g., von Schuckmann and Le Traon 2011; an extended listing of the available datasets is in Table 1 of Wunsch and Heimbach (2013a)]. Hypothetically, a highly accurate estimate of, for example, heat and salt content changes in the upper ocean, coupled with altimetric, meteorological, and so on measurements would allow inference of the deep-ocean changes as residuals in the data from the subtraction of upper-ocean contributions. The strategy used here is to exploit both this idea, and the deep data that do exist, through the vehicle of a constrained general circulation model. How well the upper ocean is determined, and thus the accuracy of the abyssal residuals so calculated, is still not so clear.

The Estimating the Circulation and Climate of the Oceans (ECCO)3 “state estimate” has been described in a number of places (e.g., Wunsch and Heimbach 2007, 2013a,b). In summary, it is a weighted least squares fit of a general circulation model [an evolved version of the Massachusetts Institute of Technology general circulation model (MITgcm); see Marshall et al. (1997) and Adcroft et al. (2004) for early forms] to the quasi-global datasets (which include the atmospheric forcing) using Lagrange multipliers. The estimate has 1° zonal resolution and a meridional resolution ranging from about 0.25° near the equator and poles to 1° at midlatitudes. An initial (then adjusted) meteorological forcing is derived from the Interim European Centre for Medium-Range Weather Forecasts (ECMWF) Re-Analysis (ERA-Interim; Dee et al. 2011). A numerical algorithm for fitting using the Lagrange multipliers is sometimes known as the “adjoint method” or in meteorology as four-dimensional variational data assimilation (4DVAR). The specific estimate used is labeled version 4, release 1, and in contrast to earlier estimates includes a full sea ice model (Losch et al. 2010; Fenty and Heimbach 2013) and extends to the North Pole [see Forget et al. (2013, unpublished manuscript) for full details].4

Very recently, Abraham et al. (2013) have published a useful discussion of the methods used both historically and today for direct ocean temperature measurements including, especially, the ongoing debates about systematic errors in the different datasets. The present state estimate uses all of the post-1991 data types they discuss, but combines them also with the continuous high-density altimetric height and other measurements as well as with the best initial estimate we could obtain of the air–sea heat transfers. Thus, the direct thermal measurements are combined with numerically much more numerous estimates of atmospheric heat transfers, implied sea level shifts, and other data.

Note that over the great volume of the oceans, the ECCO state is in slowly time-evolving, geostrophic, hydrostatic balance that, unlike most “data assimilation” products, satisfies the model equations without any artificial sources or sinks or forces. The state estimate is from the free running, but adjusted, model and hence satisfies all of the governing model equations, including those for basic conservation of mass, heat, momentum, vorticity, and so on, up to numerical accuracy.

Data assimilation schemes running over decades are usually labeled “reanalyses.” Unfortunately, these cannot be used for heat or other budgeting purposes because of their violation of the fundamental conservation laws; see Wunsch and Heimbach (2013a) for discussion of this important point. The problem necessitates close examination of claimed abyssal warming accuracies of 0.01 W m−2 based on such methods (e.g., Balmaseda et al. 2013).

As with other extant estimates, the present state estimate does not yet account for the geothermal flux at the sea floor whose mean values (Pollack et al. 1993) are of the order of 0.1 W m−2 and that are minute relative to the surface heating. But they are not negligible compared either to the vertical heat transfer into the abyss from above (measured, e.g., by κT/∂z, where κ is a vertical diffusion coefficient; cf. Emile-Geay and Madec 2009)5 or the change in atmospheric radiative forcing. Absence of this abyssal heating is one of the reasons we do not emphasize what prove to be weak trends in the state estimate.

The methodology used by Kouketsu et al. (2011, hereinafter K2011) is analogous to that employed here, although some of their inferences are different. Those differences and their possible causes are discussed later.

A total change in heat content, top to bottom, is found (discussed below) of approximately 4 × 1022 J in 19 yr for a net heating of 0.2 ± 0.1 W m−2, smaller than some published values (e.g., Hansen et al. 2005, 0.6 ± 0.1 W m−2; Lyman et al. 2010, 0.63 ± 0.28 W m−2; or von Schuckmann and Le Traon 2011, 0.55 ± 0.1 W m−2; note that differing averaging periods were used), but indistinguishable from the summary Fig. 14 of Abraham et al. (2013). Perhaps coincidentally, it is similar to the 135-yr 700-m depth ocean rate of 0.2 ± 0.1 W m−2 of Roemmich et al. (2012). On multiyear time scales accessible with a 20-yr record, the present estimate is sensitive in the upper ocean to the prior estimates of atmospheric heat transfers. In contrast, the abyssal ocean response to multiyear surface thermodynamic variability is expected to be confined to small convective regions, boundary regions of baroclinic deformation radius width, and near the equator.

Figure 3 displays the temperature and salinity census in logarithmic units at the start of the state estimate. The ocean is dominated by the very cold, intermediate salinity values of the vast abyssal interior and a calculation of net heat content change requires measurements of this cold-water sphere with volume average precisions consistent with Table 1.

Fig. 3.

Volumetric census—cubic meters of water lying in fixed intervals of temperature and salinity—of the state estimate in 1993 in logarithmic units. Total volume is about 1.3 × 1018 m3. The mean value is shown by “o” at 3.5°C and 34.8. Worthington (1981) reported a mean of 3.5°C and 34.7 ‰ on the basis of an ocean he optimistically regarded as 46% sampled. A 1 W m−2 net oceanic heating would shift the mean temperature by approximately 0.04°C in 20 yr, showing the necessity of observation of the massive cold abyssal water masses.

Fig. 3.

Volumetric census—cubic meters of water lying in fixed intervals of temperature and salinity—of the state estimate in 1993 in logarithmic units. Total volume is about 1.3 × 1018 m3. The mean value is shown by “o” at 3.5°C and 34.8. Worthington (1981) reported a mean of 3.5°C and 34.7 ‰ on the basis of an ocean he optimistically regarded as 46% sampled. A 1 W m−2 net oceanic heating would shift the mean temperature by approximately 0.04°C in 20 yr, showing the necessity of observation of the massive cold abyssal water masses.

a. Misfits

The most basic test of any least squares state estimate is the extent to which the diverse datasets have been fit to the model trajectory. A full discussion of the misfits to the approximately 2 × 109 data constraints in the estimate requires far more space than is available here. As a representative of the complete discussion, the misfits between the CTD and the state estimate in different depth ranges are shown in Figs. 4 and 5. Apart from the outliers expected in the distribution characterizing least squares residuals of Gaussian data, almost all values are close to zero and obvious basin-scale systematic offsets do not appear. Regional misfits do remain, particularly in the northern North Atlantic and parts of the Southern Ocean.

Fig. 4.

Mean differences in °C to the CTD data in the state estimate as a function of depth interval. (top to bottom) 0 to 700 m, 700 to 2000 m, 2000 to 3600 m, and 3600 m to bottom. The temperature range from top to bottom is ±1.5°, ±0.75°, ±0.5°, and ±0.25°C.

Fig. 4.

Mean differences in °C to the CTD data in the state estimate as a function of depth interval. (top to bottom) 0 to 700 m, 700 to 2000 m, 2000 to 3600 m, and 3600 m to bottom. The temperature range from top to bottom is ±1.5°, ±0.75°, ±0.5°, and ±0.25°C.

Fig. 5.

As in Fig. 4, but for the squared deviations normalized by the expected errors to produce costs. These should, on average, be consistent with a χ2 distribution with mean near unity.

Fig. 5.

As in Fig. 4, but for the squared deviations normalized by the expected errors to produce costs. These should, on average, be consistent with a χ2 distribution with mean near unity.

In a formal sense, the regions showing patterns of misfits larger than two standard deviations mean that the estimate should be rejected as inconsistent with the a priori assumptions. As with any nonlinear least squares fit, inconsistencies can arise for a number of reasons including premature termination of the optimization algorithm; inadequate model resolution or missing physics; and systematic errors in one or more of the constraining datasets. Table 2 lists the fractions of the abyssal volumes that were sampled where the misfits exceed two standard deviations. All are less than 1% of the oceanic volume. As small fractions of the total, the misfits are unlikely to have a qualitative impact on the global integrals. Interior geostrophic and hydrostatic balance tends to control the transport properties of the poorly resolved boundary currents and other special regions in calculations dominated by the very much larger interior oceanic volumes. Nonetheless, as with any least squares fit, it is a current “best estimate,” is not claimed to be “correct” in any absolute sense, and is obviously subject to quantitative change. The present solution, in terms of misfits to all of the data (whose numbers are dominated by the meteorological values, altimetry, and Argo), is deemed adequate for analysis.

b. Time scales

One of the fundamental characteristics of the ocean as it influences climate on decadal and longer time scales is its long memory—the main reason why the brevity of the instrumental record is so frustrating. Simple calculations show that the ocean responds, and thus remembers, on time scales of seconds out to thousands of years. When interpreting measurements of changes, any assumption that they have been generated by disturbances from the recent past has to be justified. The question arises specifically in the determination of the initial conditions in a calculation of change. Note that the control vector of the state estimate explicitly contains the system initial conditions—hydrography and flow.

A large number of physical mechanisms operate in the ocean as it responds dynamically and kinematically to external disturbances. Many of these adjustments will occur on time scales that are brief compared to a two-decadal time span, including such diverse mechanisms as Kelvin-like coastal and low-latitude Rossby waves, Ekman pumping changes, convective responses to changing ice cover, and changes in eddy bolus transports. Spatial scales will range from deformation radii motions and property shifts to those extending to entire ocean basins—depending directly on the physical mechanisms. On the other hand, many such processes will be present with time scales extending from multiple decades out to thousands of years. From the point of view of basin-scale heat content changes measured on a bidecadal time scale, responses will also be generated from the initial conditions in 1992. Initial values reflect any disequilibrium between modern meteorological forcing and the memory embedded in the deep ocean of fluctuations from long-ago disturbances.

From the dual (adjoint) model of the MITgcm used to obtain the state estimates, Heimbach et al. (2011) showed that changes in North Atlantic Ocean meridional heat transport exhibited a noticeable response to advected temperature changes from preceding decades and extending to great distances globally. The many mechanisms known to operate in oceanic temporal adjustment are present in the model and state estimate, and they depend strongly upon region. In temporal contrast, in another calculation employing the state estimate, Wunsch and Heimbach (2008) calculated the time for a passive tracer to reach equilibrium values over ocean basin scales, an example of which is reproduced in Fig. 6, with time scales depending upon the region, ranging from order 100 yr to nearly 10 000 yr (in the abyssal North Pacific Ocean). These long time scales are easily rationalized in a number of ways, including the diffusion times required to ultimately erase spatial gradients. For example, the diffusion e-folding times are of order L2/K, where L is a characteristic length, and K is a diagonal element of the diffusion tensor. If L ≈ 104 km (the width of the Pacific Ocean) and a horizontal diffusion coefficient is in the range of 500–1000 m2 s−1 (e.g., Ferreira et al. 2005), the characteristic time is of the order 3000 yr. Vertical distances and diffusion produce similar values. Additional long time scales can be derived, for example, from ocean volumes and their advective renewal rates.

Fig. 6.

Time in years for a passive tracer to reach 90% of its equilibrium value at 2000 m when a globally uniform concentration is imposed at the sea surface and held there (from Wunsch and Heimbach 2008). These values can be interpreted as a measure of the ocean memory and which ranges, at this depth, from several hundreds to several thousands of years. Extrapolation of the 1900-yr computation was used to estimate the much longer North Pacific equilibrium times. Note that an active tracer—one such as temperature modifying density—would have a different history, generating faster baroclinic disturbances, but final equilibrium will again likely be diffusively controlled.

Fig. 6.

Time in years for a passive tracer to reach 90% of its equilibrium value at 2000 m when a globally uniform concentration is imposed at the sea surface and held there (from Wunsch and Heimbach 2008). These values can be interpreted as a measure of the ocean memory and which ranges, at this depth, from several hundreds to several thousands of years. Extrapolation of the 1900-yr computation was used to estimate the much longer North Pacific equilibrium times. Note that an active tracer—one such as temperature modifying density—would have a different history, generating faster baroclinic disturbances, but final equilibrium will again likely be diffusively controlled.

The purpose of this paper is not the regional physics of thermal change. It is the summary estimation of the large regional changes in heat content, particularly in the abyss, as perceivable both as regional trendlike behavior on time scales exceeding the 20-yr estimate and the superimposed higher-frequency changes. These latter are interesting in their own right, but also act as a noise in attempts to determine multidecadal shifts. Regional dynamical interpretations as part of the generic problem of ocean “spinup” is left for other studies.

Depending upon geographical region, depth range, and spatial scale, changes are expected ranging from weeks, months, and years, out to those appearing as regional trends. The latter in practice may be just the expected oceanic response to past forcing—still “remembered” in the form of the continuing adjustment to the initial conditions.

To make this assertion more concrete, Fig. 7 shows one example of estimated Northern Hemisphere surface temperatures over the last 2000 yr. Translating such a curve, even if taken at face value, into a rate of atmosphere–ocean heat exchange is a major challenge. Nonetheless, for scaling purposes, suppose the approximately 0.2°C change over the last about 20 yr corresponds to an exchange between ocean and atmosphere of 1 W m−2. Then for example, the long decline from the year 1000 CE to about 1700 CE, if it too should correspond to 1 W m−2, would imply a temperature reduction of about 35 times that estimated above for a 20-yr interval. That reduction would then be overlain by the previous warming and then the rewarming over the past 300 yr. Unless existing circulation rates have been grossly underestimated, the signature of the past state must be present in any measure of basin-scale and larger heat content or temperature shifts of the past few decades. No details are available, but discovering that parts of the system are still changing in ways unconnected to the recent increase in global average temperatures would not be a surprise.

Fig. 7.

An estimate used here only for scaling purposes (Ljungqvist 2010) of Northern Hemisphere surface temperatures (ocean and land) dating to 1 CE (AD in the figure), showing multidecadal and much longer intervals of warmer and colder temperatures. The medieval warm period and the Little Ice Age are conspicuous. The gray band is the estimated two standard deviation uncertainty (likely optimistic). If translatable into air–sea heat transfers (by no means clear) then the ocean should today retain a memory of these past states as the time scales in Fig. 6 exceed this duration.

Fig. 7.

An estimate used here only for scaling purposes (Ljungqvist 2010) of Northern Hemisphere surface temperatures (ocean and land) dating to 1 CE (AD in the figure), showing multidecadal and much longer intervals of warmer and colder temperatures. The medieval warm period and the Little Ice Age are conspicuous. The gray band is the estimated two standard deviation uncertainty (likely optimistic). If translatable into air–sea heat transfers (by no means clear) then the ocean should today retain a memory of these past states as the time scales in Fig. 6 exceed this duration.

Long memory of earlier states presents a particularly difficult sampling problem in the computation of global changes. By way of example, consider a small volume of ocean with a large thermal anomaly—an anomaly resulting from long-past meteorological conditions. Suppose, in addition, that subsequently no further net heat exchange between ocean and atmosphere occurred. Unless that anomalous region has been adequately observed at the outset, as the thermal anomaly moves by advection and diffusion into other regions—that are better observed—an apparent change in heat content will be computed—only because of the spatial and temporal structure of the observing system. This possibility is one of the reasons for confining the present calculations to the comparatively recent well-observed period and for the use of all relevant data types. The state estimation approach mitigates, but does not entirely remove, the sensitivity to evolving observing systems.

3. Abyssal signals

The eddy field in the ocean appears to be rich in the lowest baroclinic mode (Wunsch 1997), which implies a major eddy noise in the deep ocean. Study of the present eddy-free motions on time scales of less than about 2 yr shows a strong coupling in both temperature and velocity between the upper and lower oceans, consistent with a primarily wind-forced response. Ponte (2012), using a different ECCO eddy-permitting state estimate (Menemenlis et al. 2005), showed that abyssal noise could seriously compromise the interpretation of sea level variability and hence degrade the heat content estimates at the levels of accuracy needed here.

a. Variances

The standard deviation of temperature variability at 2000 m is shown in Fig. 8—the central result here. For context, Fig. 9 is taken from the ECCO2, eddy-permitting state estimate of Menemenlis et al. (2005), used by Ponte (2012), which shows that the eddy-noise variance (which is likely still underestimated owing to the 18-km horizontal resolution) is about 6 times larger than the background standard deviation. Also shown (Fig. 10) is the logarithm of the ratio of the eddy-permitting variances to that of the present state vector. The considerable eddy noise is obvious although indications exist of regions in which the eddy noise is smaller than the variance of the lower-frequency shifts. Figure 11 shows the standard deviation (without eddy noise) at 3600 m. Values at both 2000 and 3600 m are small as compared to those in the thermocline. To move forward, the present analysis relies heavily on the assumption that the combination of constraints to observations and of the robust nature of the thermal wind relationships over long distances means that the state estimate faithfully tracks the large-scale thermal structures. The eddy field then represents a background noise primarily of concern in the noise-representing weights assigned to individual data points.

Fig. 8.

Base 10 logarithm of the standard deviation from monthly averages of temperature at 2000 m in the state estimate. The Atlantic and Southern Oceans carry most of the variability.

Fig. 8.

Base 10 logarithm of the standard deviation from monthly averages of temperature at 2000 m in the state estimate. The Atlantic and Southern Oceans carry most of the variability.

Fig. 9.

Estimated base 10 log of the standard deviation of temperature at 2000 m in the eddy-resolving/permitting ECCO2 state estimate. Variability on scales larger than about 3° of the latitude and longitude was suppressed to approximately isolate the synoptic eddy-scale contribution.

Fig. 9.

Estimated base 10 log of the standard deviation of temperature at 2000 m in the eddy-resolving/permitting ECCO2 state estimate. Variability on scales larger than about 3° of the latitude and longitude was suppressed to approximately isolate the synoptic eddy-scale contribution.

Fig. 10.

Base 10 logarithm of the ratio of the variance of the ECCO2 temperature variations near 2000 m (scales shorter than about 3° of latitude and longitude) to that in the version 4 state estimate without eddies. The spatial average eddy contribution is approximately 6 times the ECCO estimated variance. Regions with a logarithm below 0 are places where the interannual variability appears to exceed the eddy noise and would be of great observational significance in the unlikely event that ocean warming was expected to be globally uniform.

Fig. 10.

Base 10 logarithm of the ratio of the variance of the ECCO2 temperature variations near 2000 m (scales shorter than about 3° of latitude and longitude) to that in the version 4 state estimate without eddies. The spatial average eddy contribution is approximately 6 times the ECCO estimated variance. Regions with a logarithm below 0 are places where the interannual variability appears to exceed the eddy noise and would be of great observational significance in the unlikely event that ocean warming was expected to be globally uniform.

Fig. 11.

As in Fig. 8, but at 3600 m. Note that the high southern latitudes have become the most active regions at this depth, in contrast to the behavior nearer the sea surface.

Fig. 11.

As in Fig. 8, but at 3600 m. Note that the high southern latitudes have become the most active regions at this depth, in contrast to the behavior nearer the sea surface.

The most important result is that the standard deviation, or variance, pattern qualitatively replicates the tracer equilibrium time structure of Fig. 6. This structure is physically reasonable as regions functionally remote from atmospheric disturbances should show a muted response to short time-scale fluctuations as short wavelength and high frequencies are lost in propagation. Because the globally uniform boundary condition used for the passive tracer experiment is so different from those of atmospheric thermal disturbances, detailed resemblance is not expected.

b. Heat content regional patterns

Heat content (J m−2) between two depths z1, z2 at each horizontal location (θ, λ) is computed as

 
formula

where the heat capacity, cp = 3.8 × 103 J kg−1 °C−1, is taken as a constant. Calculation with a spatially varying cp changes nothing of significance. Figure 12 shows where heat is stored in the ocean, displaying the time-mean heat content, top to bottom, and Fig. 13 does so for the portion of the water column below 2000 m. The relatively warm North Atlantic and cold Southern and Pacific Oceans are apparent in both integrals. These patterns are important because spatial gradients, both horizontal and vertical, are determinants of the future changes in these distributions. Regions of very small horizontal gradient cannot undergo future large temporal changes from lateral advection or mixing except on very much longer time scales than the available 20 yr.

Fig. 12.

Heat content, H(0, −h) in the time mean, top to bottom using °C. Notice the strong meridional gradients at high latitudes. White contour is the boundary of mean negative temperatures and thus apparent negative heat content using a Celsius temperature scale. The relatively large heat content of the Atlantic Ocean could, if redistributed, produce large changes elsewhere in the system and which, if not uniformly observed, show artificial changes in the global average.

Fig. 12.

Heat content, H(0, −h) in the time mean, top to bottom using °C. Notice the strong meridional gradients at high latitudes. White contour is the boundary of mean negative temperatures and thus apparent negative heat content using a Celsius temperature scale. The relatively large heat content of the Atlantic Ocean could, if redistributed, produce large changes elsewhere in the system and which, if not uniformly observed, show artificial changes in the global average.

Fig. 13.

Time-mean heat content below 2000 m, H(2000, −h). The warmer Atlantic remains visible at these depths. Weak gradients in the Pacific would minimize any observed time changes owing to lateral motions or diffusion. White contour is again the boundary of zero-mean temperatures.

Fig. 13.

Time-mean heat content below 2000 m, H(2000, −h). The warmer Atlantic remains visible at these depths. Weak gradients in the Pacific would minimize any observed time changes owing to lateral motions or diffusion. White contour is again the boundary of zero-mean temperatures.

To avoid discussion of the physical accuracy of a linear or other trend, Figs. 1417 show the difference of the annual-mean values in 2011 minus those from 1993. The year 1992 is dropped as possibly showing signs of a starting transient. In the abyss, resemblances and differences to Fig. 6 can be seen. The western Atlantic and sectors of the deep Southern Ocean display a warming, with the remainder of the ocean either cooling (northwestern Indian Ocean, eastern basin of the Atlantic) or little or no change (the great bulk of the Pacific). Of most significance is the very strong regionality of the changes—expected from the numerous existing estimates of regional sea level variations.

Fig. 14.

Difference in heat content of the annual average of 2011 minus that of 1993, H(0, −h, 2011) − H(0, −h, 1993). The strong spatial structure represents a major observational challenge to determining an accurate mean change.

Fig. 14.

Difference in heat content of the annual average of 2011 minus that of 1993, H(0, −h, 2011) − H(0, −h, 1993). The strong spatial structure represents a major observational challenge to determining an accurate mean change.

Fig. 15.

As in Fig. 14, but for the top 700 m alone, H(0, −700, 2011) − H(0, −700, 1993). Annual cycle and harmonics removed. Regions of loss as well as gain depict some of the sampling difficulty. [cf. Kosaka and Xie (2013)].

Fig. 15.

As in Fig. 14, but for the top 700 m alone, H(0, −700, 2011) − H(0, −700, 1993). Annual cycle and harmonics removed. Regions of loss as well as gain depict some of the sampling difficulty. [cf. Kosaka and Xie (2013)].

Fig. 16.

As in Fig. 14, but for 2000 m to the bottom.

Fig. 16.

As in Fig. 14, but for 2000 m to the bottom.

Fig. 17.

As in Fig. 14, but for 3600 m to the bottom. Note the cooling in the deepest parts of the western North Atlantic, the entire eastern basin, and Pacific and Indian Oceans. Warming of the Antarctic Bottom Water has been discussed recently by Purkey and Johnson (2013), among others. In the present context, it is a comparatively small water mass. Warming in the Atlantic sector of the Southern Ocean is particularly conspicuous.

Fig. 17.

As in Fig. 14, but for 3600 m to the bottom. Note the cooling in the deepest parts of the western North Atlantic, the entire eastern basin, and Pacific and Indian Oceans. Warming of the Antarctic Bottom Water has been discussed recently by Purkey and Johnson (2013), among others. In the present context, it is a comparatively small water mass. Warming in the Atlantic sector of the Southern Ocean is particularly conspicuous.

At all depths, but particularly in the upper ocean, regions of warming are at least partially compensated in the global integrals by extended regions of cooling (especially the tropical Pacific and North Atlantic subtropical gyre). These patterns emphasize the problem of having adequate spatial sampling to generate mean values consistent with the accuracies in Table 1. Major unobserved regions of the abyss, as seen in Fig. 1 and Table 2, are a particular concern.

c. Time variations of global heat content

The time variations of the spatially integrated values of H,

 
formula

are shown in Fig. 18 for the integrals over varying depth ranges. The global integrals, reflecting the total ocean heat content and its changes, are problematic relative to the regional changes, representing comparatively small residuals of much larger numbers. Nonetheless, with the continuing intense interest in determining net ocean heat uptake as a confirmation of estimates of radiative forcing changes, they are calculated here because they raise, in a concrete fashion, a number of measurement issues.

Fig. 18.

Time variability of the globally integrated H(z = 0, zj, t) and denoted IH(z1, z2, t) as labeled IH(0, −100, t), IH(0, −700, t), IH(−2000, −h, t), and IH(−3600, −h, t) and the top-to-bottom integral IH(0, −h, t) in Yotta Joules (YJ; 1 YJ = 1024 J). A change of 0.1 YJ over the mean water depth of 3700 m corresponds to a temperature change of about 0.02°C.

Fig. 18.

Time variability of the globally integrated H(z = 0, zj, t) and denoted IH(z1, z2, t) as labeled IH(0, −100, t), IH(0, −700, t), IH(−2000, −h, t), and IH(−3600, −h, t) and the top-to-bottom integral IH(0, −h, t) in Yotta Joules (YJ; 1 YJ = 1024 J). A change of 0.1 YJ over the mean water depth of 3700 m corresponds to a temperature change of about 0.02°C.

The time-scale problem in models is greatly exacerbated by their known numerical drifts. ECCO state estimates have some immunity to this problem induced by the use of constraints forcing the model to those abyssal hydrographic data that do exist over the entire time interval and by constraints preventing it from moving very far from the available crude climatologies. In addition, permitting slight adjustments in the model mixing parameters served to further reduce any tendency for the model to drift.

Near-surface and total values are dominated by the annual cycle. Although the annual cycle, and its sometimes important harmonics, is comparatively well known, its large magnitude is important for the error budget of upper-ocean measurements, as even small aliases, temporal or spatial, can mask lower-frequency signals.

With the state estimate, removing the annual cycle, its first three harmonics, and the time mean of the IH is simple, with the results shown in Fig. 19. A fit of a linear trend to the global integrals with a suppressed annual cycle in the IH is also shown. In a formal sense, the apparent trends show a warming in the upper ocean and a net cooling below 2000 m. For IH(−3600, −h, t), the cooling is about 0.01°C over 19 yr. As with many climate-related records, the unanswerable question here is whether these changes are: 1) truly secular, and/or 2) a response to anthropogenic forcing, 3) whether they are instead fragments of a general red-noise behavior seen over durations much too short to depict the long time scales of Figs. 6 and 7, 4) the result of sampling and measurement biases or, 5) the result of the major changes in the temporal data density.

Fig. 19.

As in Fig. 18, but the annual cycle has been removed. Dashed–dotted lines are the best linear fits, and dashed lines are the residual. The 1997/98 ENSO event is visible primarily in IH(0, −100, t), but can also be detected below where the thermal anomaly is largely compensated. Because of the very long time scales embedded in the oceans, no particular significance is attached here to the apparent linear trends where visible, as they may well be fragments of much longer red-noise trends or systematic errors.

Fig. 19.

As in Fig. 18, but the annual cycle has been removed. Dashed–dotted lines are the best linear fits, and dashed lines are the residual. The 1997/98 ENSO event is visible primarily in IH(0, −100, t), but can also be detected below where the thermal anomaly is largely compensated. Because of the very long time scales embedded in the oceans, no particular significance is attached here to the apparent linear trends where visible, as they may well be fragments of much longer red-noise trends or systematic errors.

Time changes can sometimes be better estimated than the absolute accuracy. In the present cases, the temporal standard deviations σH(0, z) from monthly values over 20 yr are displayed in Table 3 (including the annual cycles). A rough estimate of the formal accuracy with which a temporal change can be computed between any 5-yr interval, for example, 1992–96 versus 2006–11, can be made by assuming that the 5-yr average has a standard error of , independent in the two intervals. (Because of the strong annual cycle, the monthly values are being assumed to be strongly correlated.) The difference between two estimates would have a formal standard error of , or for the total H(0, −h, t) of 1.5 × 1022 J, a heating equivalent over 20 yr of 0.07 W m−2 with a similar value for H(0, −700). Note that the apparent “pause” in global ocean heat uptake since about 2004, documented, for example, by Lyman et al. (2010, their Fig. 2), amounts to about 4 × 1022 J in about 7 yr. They show yearly 90% confidence intervals of 2–4 × 1022 J, roughly a heating error of 0.1 W m−2, consistent with those found here. Abyssal noise would contribute another 10% uncertainty to a water-column total. In any case, the small changes including the pause are at best at the very edge of what is practical precision today.

Table 3.

Standard deviation of the total heat content between the depths indicated from the 20-yr state estimate and of the equivalent heating rate over 20 yr.

Standard deviation of the total heat content between the depths indicated from the 20-yr state estimate and of the equivalent heating rate over 20 yr.
Standard deviation of the total heat content between the depths indicated from the 20-yr state estimate and of the equivalent heating rate over 20 yr.

The very important regional heterogeneity of change in heat content is obvious in the mapped figures. Temporal inhomogeneity is also considerable: Fig. 20 displays the detrended values of H(−2000, −h, t) for points in the eastern North Pacific, western North Atlantic, and the Atlantic sector of the Southern Ocean at the locations listed. Detrending was done to avoid the question of the physical nature of the lowest-frequency band. The Atlantic and Southern Oceans exhibit a great deal of excess high-frequency energy relative to the eastern Pacific Ocean, confirmed in the spectra also shown in the figure. The eastern Pacific spectral estimate shows a “redder” structure, with lower energy at all frequencies. Return time requirements for repeated sampling will evidently be different in different places.

Fig. 20.

(top) Time series of H(−2000, −h, t) from the eastern Pacific (latitude 29.5°N, longitude 155.5°W, solid curve) the western Atlantic (latitude 29.5°N, longitude 64.5°W, dashed) and the Southern Ocean (latitude 56.5°S, longitude 37.5°W, dotted). Note the visually stronger low-frequency variability from the Atlantic. (bottom) Power density spectral estimates for the three records shown above. All records approach white noise at low frequencies beyond about 10-yr period, with an order of magnitude less variance in the Pacific Ocean. A small power excess near the annual period is visible in the Atlantic values. The power laws at high frequencies s lie between about −2.2 and −3, although that characterization is oversimplified. Note that multitaper spectral methods are biased low at the longest periods. The vertical bar is an approximate 95% confidence interval.

Fig. 20.

(top) Time series of H(−2000, −h, t) from the eastern Pacific (latitude 29.5°N, longitude 155.5°W, solid curve) the western Atlantic (latitude 29.5°N, longitude 64.5°W, dashed) and the Southern Ocean (latitude 56.5°S, longitude 37.5°W, dotted). Note the visually stronger low-frequency variability from the Atlantic. (bottom) Power density spectral estimates for the three records shown above. All records approach white noise at low frequencies beyond about 10-yr period, with an order of magnitude less variance in the Pacific Ocean. A small power excess near the annual period is visible in the Atlantic values. The power laws at high frequencies s lie between about −2.2 and −3, although that characterization is oversimplified. Note that multitaper spectral methods are biased low at the longest periods. The vertical bar is an approximate 95% confidence interval.

In principle, a goal of an accuracy of 0.1 W m−2 is within reach on a decadal basis from the state estimate without having to assume anything about the form of a trend. The reader is strongly cautioned, however, that this optimistic error estimate does not include any systematic errors that are likely present in the data and the model, nor the eddy-noise contribution. Meteorological forcing errors, mainly influencing the upper ocean on a 20-yr time scale, geothermal effects in the abyss, and initial condition errors representing ongoing changes are only three of the many possibilities.

With all of the data available, the system is consistent with these comparatively small values of estimated heat content or equivalent volume-averaged temperature change. Of that total amount, approximately 10% is the contribution from below 2000 m—a value in accord with the global-mean sea level contribution portion calculated by Ponte (2012) and consistent with the estimate of K2011. It sets a limit to the precision to which an upper-ocean estimate alone can be used to calculate the change in oceanic heat storage—on this bidecadal time interval. In the convectively active regions, the abyssal contribution is much larger than 10%, and what the future holds is unknown.

d. Comparison to other studies

A large number of studies from hydrographic data of abyssal changes with some overlap with this same period have been published, usually with reference to changes in a particular region. Representative among them are Bryden et al. (1996) for the North Atlantic at 24°N, Joyce et al. (1999) for the western North Atlantic, and Purkey and Johnson (2010) for the global ocean. In these three studies, at least some of the data used are part of the ECCO state estimate (data obtained in 1992 or later), but include observations preceding that period, typically the 1980s or the 1950s [from the Atlantic survey of the International Geophysical Year (IGY)]. With the caveat that abyssal changes from the 1980s and earlier need not be the same as those occurring later, and that temporally separated hydrographic sections are contaminated by aliasing, it is still useful to briefly compare the inferred changes with those in the state estimate.

Bryden et al. (1996) inferred a weak cooling of both basins of the North Atlantic at 24°N below 2000 m in the interval 1981–92, just preceding the ECCO estimate time period. Their estimate for the longer interval 1957–92 indicated a warming to about 3000 m with cooling below.

Joyce et al. (1999), working with two 1997 meridional sections in the western Atlantic at 52° and 66°W, compared them to nearly identical measurements in the mid-1980s and to the IGY. Although a lot of detail appears, and the changes in the two available time periods are different, they found a weak indication of warming between 2000 and 3000 m at most latitudes, more pronounced in the interval 1997–IGY. The patterns are very noisy, as the ECCO estimate shows, and a major change in measurement technology took place in the interim, but again no contradiction exists with the present results.

The Purkey and Johnson (2010) study is most directly relevant, as the bulk of their data are common to the ECCO estimate—coming from the WOCE period after 1992 and later. As with most such studies, one robust inference is that noise levels are high everywhere (e.g., their Fig. 6). A strong resemblance exists between their Fig. 8a (rendered as the heating at 4000 m in 24 regional basins, constituting about 10% of the ocean volume) and Fig. 17 here. Both depict warming in the abyss at high southern latitudes, in the western basin of the Atlantic, and with cooling elsewhere. The consistency is at least reassuring, given that both studies used the same hydrographic data, but were carried out by completely different methods and with the state estimate employing a much larger and diverse dataset. The latter is dominated by altimetry and upper-ocean hydrography, but nonetheless tracks the abyssal hydrographic changes. Very different datasets are evidently qualitatively consistent.

The K2011 estimation methodology is similar to ours: a GCM at 1° resolution and 46 vertical layers [version 3 of the Geophysical Fluid Dynamics Laboratory (GFDL)/NOAA Modular Ocean Model (MOM3); Pacanowski and Griffies 2000] was used in combination with temperature data to estimate abyssal warming. Among the numerous differences, apart from the model itself, are that they combined the Green function technique of Menemenlis et al. (2005) with the Lagrange multiplier method: only temperature and salinity data were used, but the denser global observing system observations, including Argo, satellite altimetry, and scatterometry, were omitted; the datasets extended back to 1985; and the computation was run over the 40 yr beginning in 1957. In comparing their results to those of Purkey and Johnson (2010), K2011 used a much finer breakdown into 73 abyssal regions, presumably leaving a larger average residual noise level in each.

Given the numerous differences ranging from the model change to the very different database (although the 1992–present hydrography would be common to both), it is unsurprising that the K2011 results differ in some ways from the present ones, but the similarities are significant. They find regions below 3000 m of decadal-scale cooling, confined primarily to the Indian Ocean and eastern North Atlantic. On the other hand, although parts of the Pacific Ocean between 3000 and 4000 m are estimated to have been cooling, in contrast with the present results, they showed a general warming below that, albeit rather weak between 4000 and 5000 m of roughly 2–3 × 10−3 °C decade−1, and with the region below 5000 m (which we have not separated out) showing considerable warming along with the general Southern Ocean. These numbers are sufficiently small that omission of the geothermal heating is a serious concern.

Distinguishing the differences between the various estimates becomes a complex problem in defining the systematic errors, which include the details of datasets used in each study, the assumed data and representation errors, and the residual misfits of the solutions. As noted repeatedly, the available database is extremely limited, especially before 1992. The state of the art does not permit resolving these differences. Hence, the main issue facing the oceanographic community is to obtain future data so that such ambiguities do not persist into the next several decades of change.

A number of papers have appeared recently (e.g., Purkey and Johnson 2013), focusing on changes in the Antarctic Bottom Water mass, and many discussions of other regional water mass property changes have also been published. If Antarctic Bottom Water is defined as the abyssal volume with temperatures below 1°C, it constitutes less than 2% of the oceanic water volume (Worthington 1981, his Table 2.5). In the present estimate, even should significant misfits to temperature changes exist there, the global average values will not be measurably modified. A review of changes in individual water masses and their geography is beyond our present scope.

4. Sampling without the model

Most published estimates of oceanic heat content change have not employed a state estimate, but are generally described as being based upon the data alone and necessarily are commonly focused on the upper ocean. As already noticed above, the heating of the upper 700 m of the ocean by 1 W m−2 for 20 yr implies a temperature change of about 0.2°C as a water-column total. Although the upper ocean is not the focus here, an interesting and complex question is whether the observational network is capable of producing estimates of small changes with a useful accuracy. Abraham et al. (2013) have described many of these calculations in detail and provide a list of references.

Some calculations have employed so-called empirical orthogonal functions (EOFs), or singular vectors, from models not unlike the one underlying the ECCO state estimate. These are the eigenvectors of the space–time correlation matrix of the model output used as an expansion basis. For example, the 240 monthly estimates from the present ECCO state estimates define 240 orthonormal vectors whose sum can perfectly reproduce either the global temperature at any depth or the heat content. Only 240 accurate measurements of the corresponding field would be adequate. As the number of data in each month tends to greatly exceed that value (see Fig. 21), obtaining high accuracy appears easy.

Fig. 21.

Number of observations extending below 2000 m for each year (solid curve) and below 3600 m (dashed). Upper-ocean observations (not shown) greatly increase with the Argo array from the mid-2000s, introducing an important inhomogeneity with time in the estimates.

Fig. 21.

Number of observations extending below 2000 m for each year (solid curve) and below 3600 m (dashed). Upper-ocean observations (not shown) greatly increase with the Argo array from the mid-2000s, introducing an important inhomogeneity with time in the estimates.

This description of the procedure is, however, too facile. The correlation matrix eigenvectors are dependent upon the accuracy and stability of the matrix, as well as the differences in numerical values of the corresponding eigenvalues. Calculation of the resulting accuracy and stability from a finite time duration involves the underlying spatially inhomogeneous, four-dimensional space–time statistics of the state estimate. An additional problem arises when those same eigenvectors are employed for time spans much exceeding that of the model or state estimate duration—as the long oceanic memory implies evermore physical regimes will come into play with longer times.

5. Discussion, with comments on the observation problem

a. Bidecadal abyssal change

Over the 20 yr of the present ECCO state estimate, changes in the deep ocean on multiyear time scales are dominated by the western Atlantic basin and Southern Oceans. These are qualitatively consistent with expectations there of the comparatively rapid response to surface forcing. Defining the physics of those changes in terms of boundary currents, wave propagation, eddy diffusion, and the myriad of other oceanic physical processes, region by region, remains a major unfinished piece of business. In those same regions, a longer-term general warming pattern occurs below 2000 m. A very weak long-term cooling is seen over the bulk of the rest of the ocean below that depth, including the entirety of the Pacific and Indian Oceans, along with the eastern Atlantic basin. Some of this change is interpreted here as owing to a disequilibrium of the abyssal ocean to the present atmosphere, with a superimposed multiyear noise. The pattern below 3600 m is similar, with much smaller amplitude. These results differ in detail and in numerical values from other estimates, but determining whether any are correct is probably not possible with the existing datasets.

The globally integrated heat content changes involve small differences of the much larger regional changes. As existing estimates of the anthropogenic forcing are now about 0.5 W m−2, the equivalent global ocean average temperature changes over 20 yr are mostly slight compared to the shorter-term temporal variations from numerous physical sources. Small errors in data calibration, and space–time sampling and model biases, are important. Direct determination of changes in oceanic heat content over the last 20 yr are not in conflict with estimates of the radiative forcing, but the uncertainties in all the fields remain too large to rationalize, for example, the apparent pause in warming. The challenge is to develop observations so that future changes can be made with accuracies and precisions consistent with the conventional rule of thumb that they should be better than 10% of the expected signal.

b. Comments on future observations

No observing system can be designed and deployed that is capable of addressing all possible goals; specification of the particular purposes and the related accuracies and precisions is essential. Here the context of the discussion is (i) the global distribution by basin and (ii) directed at the problem of the determination of full-water-column changes in temperature (and salinity) over multiple decades. Although these choices are arbitrary to a degree, they address the important problems of sea level change and of ocean heat uptake and are basic to classical scientific understanding of how the ocean varies through time and space. Absent a full instrument deployment optimization, a plausible strategy for moving forward is to concentrate abyssal samples where both the largest short-term signals are appearing (western basin of the Atlantic, the Southern Ocean) and with the highest noise levels, with only sporadic checks in the Pacific.

Ponte (2012) has summarized the abyssal measurement problem and its possibilities. In situ abyssal measurements by Argo profilers will likely become available in the next few years (D. Roemmich 2013, personal communication). Acoustic tomographic measurements are another method for direct abyssal measurements. Satellite gravity data, such as are now available from the Gravity Recovery and Climate Experiment (GRACE) mission (Tapley et al. 2004), produce estimates of the bottom pressure fluctuations. In discussions of how to ultimately construct a feasible and useful global-scale observing system by any or all means, it is essential to define the magnitude of the signals sought, and the structure in space and time of the noise field that tends to obscure those signals. Conceivably, a continuation of the existing hydrographic sampling is adequate for some purposes.

(Although not yet analyzed, calculated salinity changes are expected to display some resemblance to those for temperature, but not to be identical, as the relevant observing technology differs considerably, as do the boundary and initial conditions. With some additional effort, the ECCO state estimate can be used to calculate the structure of changes in other properties such as oxygen, carbon, silica, and so on, which are likely to be undergoing very different space and time evolution.)

That the noise level is also greatest (Figs. 8, 9) where the largest changes appear is a challenge to any observing system. With a fuller understanding of the noise level, particularly of the abyssal eddy field, various strategies can be developed for basin-scale and global measurements of changing heat and, mutatis mutandis, the salinity and other fields. With growing confidence in the ECCO estimates, one practical strategy is to maintain a modestly augmented version of the existing observing system (“modest” in the sense of cost and ease of effort in sustaining it) and to focus on observational tests of the state estimate structures in crucial regions.

Acknowledgments

Supported in part by the National Aeronautics and Space Administration and the National Ocean Partnership Program. We had useful comments from G. Forget, M. Tingley, P. Huybers, and two anonymous reviewers. The results rest on the efforts of members of the ECCO Consortium as a whole.

REFERENCES

REFERENCES
Abraham
,
J. P.
, and Coauthors
,
2013
:
A review of global ocean temperature observations: Implications for ocean heat content estimates and climate change
.
Rev. Geophys.
,
51
,
450
483
, doi:.
Adcroft
,
A.
,
C.
Hill
,
J.-M.
Campin
,
J.
Marshall
, and
P.
Heimbach
,
2004
: Overview of the formulation and numerics of the MIT GCM. Proc. ECMWF Seminar on Recent Developments in Numerical Methods for Atmospheric and Ocean Modelling, Shinfield Park, Reading, United Kingdom, ECMWF, 139–150. [Available online at http://old.ecmwf.int/publications/library/do/references/show?id=86400.]
Balmaseda
,
M. A.
,
K. E.
Trenberth
, and
E.
Kallen
,
2013
:
Distinctive climate signals in reanalysis of global ocean heat content
.
Geophys. Res. Lett.
,
40
,
1754
1759
, doi:.
Bryden
,
H. L.
,
M. J.
Griffiths
,
A. M.
Lavin
,
R. C.
Millard
,
G.
Parrilla
, and
W. M.
Smethie
,
1996
:
Decadal changes in water mass characteristics at 24°N in the subtropical North Atlantic Ocean
.
J. Climate
,
9
,
3162
3186
, doi:.
Chelton
,
D. B.
, and
F. J.
Wentz
,
2005
:
Global microwave satellite observations of sea surface temperature for numerical weather prediction and climate research
.
Bull. Amer. Meteor. Soc.
,
86
,
1097
1115
, doi:.
Church
,
J. A.
, and Coauthors
,
2011
:
Revisiting the earth's sea-level and energy budgets from 1961 to 2008
.
Geophys. Res. Lett.
,
38
, L18601, doi:.
Dee
,
D. P.
, and Coauthors
,
2011
:
The ERA-Interim reanalysis: Configuration and performance of the data assimilation system
.
Quart. J. Roy. Meteor. Soc.
,
137
,
553
597
, doi:.
de Jode
,
G.
,
1578
:
Speculum Orbis Terrarum. Theatrum Orbis Terrarum, 277 pp
.
Emile-Geay
,
J.
, and
G.
Madec
,
2009
:
Geothermal heating, diapycnal mixing and the abyssal circulation
.
Ocean Sci.
,
5
,
203
217
, doi:.
Fenty
,
I. G.
, and
P.
Heimbach
,
2013
:
Coupled sea ice-ocean-state estimation in the Labrador Sea and Baffin Bay
.
J. Phys. Oceanogr.
,
43
,
884
904
, doi:.
Ferreira
,
D.
,
J.
Marshall
, and
P.
Heimbach
,
2005
:
Estimating eddy stresses by fitting dynamics to observations using a residual-mean ocean circulation model and its adjoint
.
J. Phys. Oceanogr.
,
35
,
1891
1910
, doi:.
Hansen
,
J.
, and Coauthors
,
2005
:
Earth’s energy imbalance: Confirmation and implications
.
Science
,
308
,
1431
1435
, doi:.
Heimbach
,
P.
,
C.
Wunsch
,
R. M.
Ponte
,
G.
Forget
,
C.
Hill
, and
J.
Utke
,
2011
:
Timescales and regions of the sensitivity of Atlantic meridional volume and heat transport magnitudes: Toward observing system design
.
Deep-Sea Res. II
,
58
,
1858
1879
, doi:.
Huybers
,
P.
, and
C.
Wunsch
,
2010
:
Paleophysical oceanography with an emphasis on transport rates
.
Annu. Rev. Mar. Sci.
,
2
,
1
34
, doi:.
Joyce
,
T. M.
,
R. S.
Pickart
, and
R. C.
Millard
,
1999
:
Long-term hydrographic changes at 52 and 66°W in the North Atlantic Subtropical Gyre & Caribbean
.
Deep-Sea Res. II
,
46
,
245
278
, doi:.
Kosaka
,
Y.
, and
S.-P.
Xie
,
2013
:
Recent global-warming hiatus tied to equatorial Pacific surface cooling
.
Nature
,
501
,
403
407
, doi:.
Kouketsu
,
S.
, and Coauthors
,
2011
:
Deep ocean heat content changes estimated from observation and reanalysis product and their influence on sea level change
.
J. Geophys. Res.
,
116
, C03012, doi:.
Ljungqvist
,
F. C.
,
2010
: A new reconstruction of temperature variability in the extra-tropical northern hemisphere during the last two millennia. Geogr. Ann.,92, 339–351, doi:.
Losch
,
M.
,
D.
Menemenlis
,
J. M.
Campin
,
P.
Heimbach
, and
C.
Hill
,
2010
:
On the formulation of sea-ice models. Part 1: Effects of different solver implementations and parameterizations
.
Ocean Modell.
,
33
,
129
144
, doi:.
Lyman
,
J. M.
,
S. A.
Good
,
V. V.
Gouretski
,
M.
Ishii
,
G. C.
Johnson
,
M. D.
Palmer
,
D. M.
Smith
, and
J. K.
Willis
,
2010
:
Robust warming of the global upper ocean
.
Nature
,
465
,
334
337
, doi:.
Marshall
,
J.
,
A.
Adcroft
,
C.
Hill
,
L.
Perelman
, and
C.
Heisey
,
1997
:
A finite-volume, incompressible Navier Stokes model for studies of the ocean on parallel computers
.
J. Geophys. Res.
,
102
,
5753
5766
, doi:.
Menemenlis
,
D.
, and Coauthors
,
2005
: NASA supercomputer improves prospects for ocean climate research. Eos, Trans. Amer. Geophys. Union,86, 89–96, doi:.
Myhre
,
G.
, and Coauthors
,
2014
: Anthropogenic and natural radiative forcing. Climate Change 2013: The Physical Science Basis, T. F. Stocker et al., Eds., Cambridge University Press, 659–740.
Pacanowski
,
R. C.
, and
S. M.
Griffies
,
2000
: The Modular Ocean Model (MOM) 3 manual. Geophysical Fluid Dynamics Laboratory (GFDL) Tech. Doc., 680 pp.
Pollack
,
H. N.
,
S. J.
Hurtrer
, and
J. R.
Johnson
,
1993
:
Heat flow from the earth’ s interior: Analysis of the global data set
.
Rev. Geophys.
,
31
,
267
280
, doi:.
Ponte
,
R. M.
,
2012
:
An assessment of deep steric height variability over the global ocean
.
Geophys. Res. Lett.
,
39
, L04601, doi:.
Purkey
,
S. G.
, and
G. C.
Johnson
,
2010
:
Warming of global abyssal and deep Southern Ocean waters between the 1990s and 2000s: Contributions to global heat and sea level rise budgets
.
J. Climate
,
23
,
6336
6351
, doi:.
Purkey
,
S. G.
, and
G. C.
Johnson
,
2013
:
Antarctic Bottom Water warming and freshening: Contributions to sea level rise, ocean freshwater budgets, and global heat gain
.
J. Climate
,
26
,
6105
6122
, doi:.
Reynolds
,
R. W.
,
T. M.
Smith
,
C.
Liu
,
D. B.
Chelton
,
K. S.
Casey
, and
M. G.
Schlax
,
2007
:
Daily high-resolution-blended analyses for sea surface temperature
.
J. Climate
,
20
,
5473
5496
, doi:.
Roemmich
,
D.
, and
C.
Wunsch
,
1984
:
Apparent change in the climatic state of the deep North Atlantic Ocean
.
Nature
,
307
,
447
450
, doi:.
Roemmich
,
D.
, and Coauthors
,
2009
:
The Argo program observing the global ocean with profiling floats
.
Oceanography
,
22
,
34
43
, doi:.
Roemmich
,
D.
,
W. J.
Gould
, and
J.
Gilson
,
2012
:
135 years of global ocean warming between the Challenger expedition and the Argo programme
.
Nat. Climate Change
,
2
,
425
428
, doi:.
Roquet
,
F.
, and Coauthors
,
2013
:
Hydrographic data collected by seals significantly reduce the observational gap in the Southern Ocean
.
Geophys. Res. Lett.
,
40
,
6176
6180
, doi:.
Talley
,
L. D.
,
2007
: Hydrographic Atlas of the World Ocean Circulation Experiment (WOCE) Volume 2: Pacific Ocean. M. Sparrow, P. Chapman, and J. Gould, Eds., WOCE International Project Office, Southampton, United Kingdom, 20 pp., plus plates. [Available online at http://www-pord.ucsd.edu/whp_atlas/pacific_index.html.]
Tapley
,
B. D.
,
S.
Bettadpur
,
J. C.
Ries
,
P. F.
Thompson
, and
M. M.
Watkins
,
2004
:
GRACE measurements of mass variability in the earth system
.
Science
,
305
,
503
505
, doi:.
von Schuckmann
,
K.
, and
P.-Y.
Le Traon
,
2011
:
How well can we derive global ocean indicators from Argo data?
Ocean Sci.
,
7
,
783
791
, doi:.
Worthington
,
L. V.
,
1981
: The water masses of the world ocean: Some results of a fine-scale census. Evolution of Physical Oceanography: Scientific Surveys in Honor of Henry Stommel, B. A. Warren and C. Wunsch, Eds., MIT Press, 42–69.
Wunsch
,
C.
,
1997
:
The vertical partition of oceanic horizontal kinetic energy
.
J. Phys. Oceanogr.
,
27
,
1770
1794
, doi:.
Wunsch
,
C.
, and
P.
Heimbach
,
2007
:
Practical global oceanic state estimation
.
Physica D
,
230
,
197
208
, doi:.
Wunsch
,
C.
, and
P.
Heimbach
,
2008
:
How long to ocean tracer and proxy equilibrium?
Quat. Sci. Rev.
,
27
,
637
651
, doi:.
Wunsch
,
C.
, and
P.
Heimbach
,
2013a
: Dynamically and kinematically consistent global ocean circulation state estimates with land and sea ice. Ocean Circulation and Climate, 2nd ed. J. C. G. Siedler, W. J. Gould, and S. M. Griffies, Eds., Elsevier, 553–579.
Wunsch
,
C.
, and
P.
Heimbach
,
2013b
:
Two decades of the Atlantic meridional overturning circulation: Anatomy, variations, extremes, prediction, and overcoming its limitations
.
J. Climate
,
26
,
7167
7186
, doi:.

Footnotes

1

The Fifth Assessment report of the Intergovernmental Panel on Climate Change (Myhre et al. 2014) estimated (their Table 8.6) a total, ocean and land, anthropogenic net radiative forcing in the wide range of 1.1 to 3.3 W m−2 for the years 1750–2011. The oceanic portion is some uncertain fraction of the total. 1 W/m2 corresponds to a total of about 360 TW.

2

A nice example can be seen in de Jode (1578).

3

Here using the Massachusetts Institute of Technology and Atmospheric and Environmental Research, Inc. (MIT–AER), version.

4

The final state estimate is obtained from the free-running forward model, using the adjusted control parameters. In this particular case, the inference of a calibration discrepancy between the infrared estimate of sea surface temperature and that of the Advanced Microwave Scanning Radiometer, whose data became available in 2002, led to a small ad hoc adjustment of the imposed surface air–temperature field in the final calculation. Both sea surface temperature products are discussed by Reynolds et al. (2007); see also Chelton and Wentz (2005).

5

A vertical temperature gradient of 1°C (1000 m)−1 and a (low) eddy diffusion coefficient of 10−5 m2 s−1 produce a diffusive heat transport of about 0.04 W m−2. A value of 0.1 W m−2 is about 34TW energy input globally.