In situ measurements of thermal diffusivity in sediments of the methane-rich zone of Cascadia Margin, NE Pacific Ocean

Thermal diffusivity (TD) is defined as the ratio of the ability of a solid to conduct thermal energy to its capacity for storing thermal energy (Incropera et al., 2013); this ratio best describes how the solid responds to time-varying external thermal forcing. The TD of marine sediments is very poorly characterized due to the challenge of acquiring quantitative in situ values in an extreme environment. Previous measurements have used laboratory fabricated sediments or sediments recovered from coring (Hurwitz et al., 2012); methods which do not reproduce the in situ environment due to either mechanical disturbance or the ephemeral nature of pore water hydrate (liquid-gas) phases. Furthermore, a time-series of measurements on the order of days or longer is required to accurately compute TD for penetration depths of adequate length to be representative of the sediment column. Acquisition of such data is possible, however, with the use of a remotely operated vehicle (ROV) or similar submersible. To our knowledge, in situ thermal diffusivity measurements have not been obtained previously from any portion of the Cascadia Margin (Canadian, Washington, Oregon, northern California), although they have been made in other environments and in other sediment types (Jackson and Richardson, 2000;Wheatcroft et al., 2007). Near shore measurements of sediment TD are also uncommon, although the acquisition procedures are far less demanding, as described by Thomson (2010) for the collection of sediment TD values in the very shallow coastal waters of Washington State.

The goal of this experiment was to determine the range of TD values for marine sediments in 500 to 1100 m of water, which includes the upslope boundary of the methane hydrate stability zone on the Cascadia Margin in the NE Pacific. The Washington portion of the Cascadia North American margin is a focus area for two US National Science Foundation programs (GeoPRISM and EarthScope). The accretionary wedge on this margin has been described at considerable length in the literature (e.g. Schmalzle et al., 2014; McCrory et al., 2014 and references within), as have the sediment types and physical properties for the alternating layers of pelagic clays and sandy silt turbidites on the Washington margin (see Davis and Hyndman, 1989; Atwater et al., 2014 and references within). Heat flow on the Washington margin varies from 120 mW m^–2 at the deformation front to less than 50 mW m^–2 near the edge of the shallow shelf (Johnson et al., 2013).

Characterizing the methane hydrate stability zone at mid-latitude continental margins has global significance, as methane hydrate deposits are the most climate-sensitive reservoirs of carbon on Earth. Globally, an estimated 99% of the hydrate reservoir is held in continental margin sediments at depths below 500 m (Collett et al., 2009; Ruppel, 2011). Slight increases in near-bottom water temperatures can potentially destabilize large quantities of hydrate, releasing bubbles of methane gas with 30 times the greenhouse potency of CO₂ (Denman et al., 2007; Riedel et al., 2010). One consequence of global warming over the past 100 years is the impact of rising temperatures of the near-bottom seawater that is in physical contact with the large deposits of methane hydrate present on continental slopes. Because the stability of methane hydrate within these massive carbon reservoirs depends on both temperature and pressure, warming seawater will cause dissociation of the hydrate phase and release of CH₄ gas on a time scale of only decades (Boswell and Collett, 2011; Ruppel, 2011; Phrampus and Hornbach, 2012). The time scale over which these thermal disturbances propagate into the sediment-hosted methane hydrate deposits depends directly on the thermal diffusivity of the uppermost seafloor. Sediment thermal diffusivity is a physical property that is only poorly known from laboratory studies; very few measurements have been made in situ.

Previous studies of contemporary methane hydrate decomposition have focused on high-latitude Arctic regions (Berndt et al., 2014; Westbrook et al., 2009), although more recent studies have shown that similar hydrate dissociation is occurring at mid-latitudes, with large volumes of methane currently being released on the Cascadia portion of the North American continent (Hautala et al., 2014). This release of oxidized carbon from hydrate decomposition can lower the pH of surrounding seawater, threatening local marine organisms sensitive to acidification. Key to understanding the dynamic upslope limit of methane hydrate stability is the rate of heat transfer from warming bottom water into the overlying sediments; accurate values of thermal diffusivity are required to advance this field of inquiry.

TD measurements were acquired during the R/V ATLANTIS cruise AT26-04, where the primary research goal was to study heat flow and fluid flux on the Washington continental margin (Johnson et al., 2013). During this cruise, an instrumented fiberglass probe was inserted into the seafloor for approximately 6 days at each of two sites on the continental slope; one deployment was at 552 m depth and the other at 1049 m. Site depths were selected to include the upper limit of methane hydrate stability and an additional region well below this depth to provide a control for future hydrate stabilities. Prior to deployment, high-resolution swath bathymetry and water column data were inspected to select sites with uniform low-slope bathymetry and no methane bubbles detected in the water-column data. In addition, a visual inspection of the selected deployment locations were performed via the high-resolution camera of the ROV Jason II; the probe was deployed only in flat, heavily sedimented regions with no visible evidence of hard surface material exposures (carbonate or hydrate) or fluid discharge.

The TD probe was approximately one meter in length and instrumented with four temperature sensors; it was inserted using the manipulator of the ROV Jason II. The probe was designed to hold four Onset TidbiT v2 temperature loggers at three sub-surface spacings: 5 cm, 10 cm, and 15 cm. To reduce vertical resistance during insertion and provide support during deployment, a 5.1 × 5.1 cm fiberglass L-channel with low thermal conductivity served as a rigid backbone for the probe. Mounting configuration ensured that no two consecutive sensors would lie directly above another to minimize insertion disturbance to the sediment. A schematic of the probe is shown in Figure 1.

The solid epoxy-encased TidbiT Version 2 sensors had a rated absolute accuracy of ± 0.21 °C, resolution of 0.02 °C, and a factory-rated depth limit of only 305 m. As deployment depths would exceed this limit, a pre-cruise pressure test was necessary to identify sensors that might fail. All four sensors were pressure-tested to a simulated depth of 700 m in the Pressure Test Vessel of the School of Oceanography at the University of Washington. All sensors passed the pressure test and were certified for use in the field.

All temperature sensors were calibrated prior to the cruise, as required for marine temperature sensors (Johnson et al., 2010). This laboratory calibration procedure included measurement of the offset of each thermistor with respect to a primary sensor (an Antares thermistor) by long-term soaking in a stirred ice-bath in a sediment core cold room (4°C). This calibration was done after the sensors had been placed in the high pressure vessel to accurately determine their sensitivity when deployed at water depths deeper than their factory certification. These offset values are reported in Table S1. For the calculation of thermal diffusivity, the need for highly accurate relative temperatures is much greater than the absolute temperatures measured, as it is the difference between sensor temperatures over the different measurement depths that determines the thermal diffusivity.

The first probe deployment, at a water depth of 552 m, began on August 4, 2013 and lasted 147 h; the second deployment, at a depth of 1049 m, began on August 15 and lasted 140 h (Figure 2). The specifics of both deployments are listed in Table 1. All four thermistors logged temperature at 1 minute intervals from their respective positions on the probe. This configuration ensured that one sensor remained 5 cm above the sediment-water interface (gray bar in Figure 1), that the first buried sensor was 5 cm below the interface, and the last two sensors were at consecutive 5 cm depth intervals deeper within the sediment. Upon recovery of the probe, the calibration offsets were applied to the data acquired from each individual sensor.

The maximum observed temperature variation for a single sensor during the first deployment was ± 0.4°C, more than 25 times the resolution of the sensor. This bottom water temperature variation was sufficient to cause a fluctuation of ± 0.15°C at 15 cm, the maximum depth sampled. This amplitude variation is 7.5 times the resolution of the calibrated sensor, resulting in a dataset with statistically significant values for calculating TD at this site. A series of oscillations in bottom water temperature with a clear tidal period of approximately 24 h were observed in these data. There was also a systematic decrease in bottom water temperature observed over the full duration of the deployment (Figure 3). This observed long-term trend could be due to seasonal variability or other long time scale phenomena not fully captured in the six day window.

The maximum temperature range for the diurnal variation of ± 0.25°C observed during the second, deeper deployment was similar to the first dataset. This 0.15°C loss of variation between sites only reduced the ratio between sensor resolution and magnitude of water temperature variation from a factor of 25 to 20. Although smaller in amplitude, the frequency of the bottom water temperature oscillations at the deeper site was unexpectedly higher than at the shallower site, with an average period of 7 h (Figure 4). The smaller amplitude temperature variation at the second site, combined with the much higher frequency of temperature oscillations, resulted in a considerably attenuated signal at greater depth in the sediment, with the 15 cm-deep sensor only observing two temperature steps through the entire deployment.

To estimate temperatures within the sediment column at depth z and time t due to bottom water thermal forcing, a range of TD values was modeled in Matlab Version R2011a using the differences between temperature at a given time from the bottom water sensor and the sensor at the depth of interest; see Text S1 for details of method. The program was custom-written and is also included in the Text S1. The analysis was run with 100 input TD values; the standard deviation (SD) between the model and true temperature values at each depth was calculated for each of the input TDs. These SDs were then plotted for the two deployments (Figure 5).

For each depth, the TD value corresponding to the lowest standard deviation between the modeled and observed temperature variations was chosen as the best approximation for the true TD. This process was repeated for the data from the second deployment, resulting in a total of six TD values; two for each sensor depth within the sediment column (Table 2). The model fit, calculated using the chosen TD value, is plotted for each sediment depth on Figures 3 and 4. The discretization of input TD values introduced a maximum error of 3.26% for each final TD value, as the nature of an input matrix cannot resolve an exact output value.

In order to compare the overall variation in thermal diffusivity between the two deployment depths, a second analysis was performed. This computation followed the same procedure as for the previous analysis, but the datasets from all three sensor depths at each site were analyzed using the same TD value. For each artificial TD value produced by the thermal model, the standard deviations were summed for all three depth intervals. Thus, for each artificial TD value, a measure of its fit across the full dataset was obtained. The TD value with the best fit (lowest summed deviation) for the full dataset was then chosen as the best representative value for the depth-integrated profile for each site. The purpose of this final technique was to reduce the noise for a series of calculations in what is assumed to be a 15-cm thick sediment layer of homogenous thermal diffusivity. The results yielded a TD of 4.33 × 10^–7 m² s^–1 for the first shallower location and 1.15 × 10^–7 m² s^–1 for the second, deeper location.

The magnitude of the temperature change recorded by the sensors has a direct effect on the accuracy of our final result. Thus, for the model runs where the corresponding dataset had greater variability, the range of calculated standard deviation values is likewise greater, making possible the selection of a minimum value with a high degree of confidence (see Figure 5). However, for the runs where the variability of the corresponding dataset was low (specifically the 15-cm depth of the second deployment), there was almost no difference in the standard deviation values, producing an unreliable TD estimate. For the two deployments conducted during this project, only the 15-cm sensor from the second deployment fell into this unreliable category, as all other sensors recorded significant variability and produced distinct standard deviation curves. Here, our significance criterion was a variation in temperature at least three times the resolution of the sensor over a temperature oscillation cycle. For our experiment, this translates to a change in temperature of at least 0.06°C in 14 h or less. The accuracy of our results were also affected by the insertion of the probe, as differences between the true and model sensor depths have the potential to introduce errors. A full sensitivity analysis was performed to address this and is included in the Text S1.

Turcotte and Schubert (1982) described a simple method for estimating penetration depth (also called diffusion length) of external thermal forcing:

Here L is the diffusion length, α is the overall thermall diffusivity found at each deployment locations, and Δt is the time interval of interest. Calculations using overall diffusivity values from each deployment were performed for six different time intervals of interest (listed in Table 3). Penetration depths of 0.10, 3.69, and 11.68 m for bottom water forcing periods of Δt = 7 h, 1 y, and 10 y using the overall diffusivity value were calculated from the first deployment (552 m water depth). The estimates for the second deployment (1049 m water depth) found diffusion lengths of 0.05, 1.90, and 6.02 m for time intervals of 7 h, 1 y, and 10 y. Though the diffusion length solutions to this equation are low in accuracy and do not account for expected variations in physical properties with depth due to compaction, they provide a reasonable first-order approximation based on actual in situ thermal diffusivity measurements and highlight the potential for deeper penetration of warming bottom waters at the upper boundary of this methane hydrate zone.

Recent studies (Pohlman et al., 2009; Phrampus and Hornbach, 2012; Brothers et al., 2014; Berndt et al., 2014; Hautala et al., 2014) have shown that increasing seawater temperatures at mid-water depths on continental margins globally are releasing methane hydrate-derived carbon at flux rates that can provide positive feedback to the present level of global warming. In some cases, hydrate-induced slope failures on the continental margin slopes can produce tsunamis capable of inundating coastal communities (Rao et al., 2002; Lopez et al., 2010). The Cascadia Margin in the NE Pacific has both over-steepened slopes and methane-rich sediments containing abundant solid hydrate that is vulnerable to warming-induced slope failure (Riedel et al., 2002; Booth-Rea et al., 2008; Torres et al., 2009). There is considerable societal need to examine the tsunami hazards associated with the present and on-going phase dissociation associated with ocean warming: our current experimental values of in situ thermal diffusivity can contribute to the boundary conditions of these models.

This study aimed to provide estimates of the thermal diffusivity of marine sediments on the Cascadia Margin through in-situ measurements acquired using a simple, yet highly functional device. The resulting, new in situ TD data for sediments located on a methane-rich continental margin span the relevant ocean depth interval where gas hydrates are known to be dissociating (Torres et al., 2009; Johnson et al., 2013). Eight TD values were obtained: three depth-specific values for each of two deployment locations, and two depth-independent values that characterize the surface sediments of the two deployment sites. The first deployment, at a water depth of 552 m, had an overall TD of 4.33 × 10^–7 m² s^–1, while the second deployment, at a water depth of 1049 m, had an overall TD value of 1.15 × 10^–7 m² s^–1. Using a simple analytical model for heat transfer, these diffusivity values were used to predict 1/e thermal penetration depths. These penetration depths can be improved using more complex, realistic numerical models, but are the first estimates based on data obtained in situ from a hydrate and methane-rich continental margin. These new data can be used to improve the ability to model and predict how rising seawater temperatures in the future may impact the relatively unstable reservoir of methane–derived carbon, which has the potential to provide a strong positive feedback to the already substantial anthropogenic greenhouse gas emissions.

Prior to publication, all data used in this study, including site location and thermistor logs, will be posted on the GeoPRISM data archives and the National Geophysical Data Center database, as part of the NSF Open Access policy for all data from this cruise.

© 2015 Homola, Paul Johnson and Hearn. This is an open-access article distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited.

The crew of the R/V Atlantis and operating personnel for the ROV Jason II were essential to the success of this study. Special thanks are given to Tor Bjorklund for assistance with design and construction of the temperature probe, and to Robert Harris, OSU, for his assistance and advice during the acquisition and processing phases of this project.