Can C-Band SAR be used to estimate soil organic carbon storage in tundra?

A new approach for the estimation of Soil Organic Carbon (SOC) pools North of the tree line has been developed based on synthetic aperture radar data (SAR). SOC values are directly determined from backscatter values instead of upscaling using land cover or soil classes. The multi-mode capability of SAR allows application across scales. It can be shown that measurements in C-band under frozen conditions represent vegetation and surface structure properties which relate to soil properties, specifically SOC. It is estimated that at least 29 PgC are stored in the upper 30 cm of soils North of the tree 5 line. This is approximately 25% less than stocks derived from the soil map based Northern Circumpolar Soil Carbon Database (NCSCD). The total stored carbon is underestimated since the established empirical relationship is not valid for peatlands as well as strongly cryoturbated soils. The approach does however provide the first spatially consistent account of soil organic carbon across the Arctic. Furthermore, it could be shown that values obtained from 1 km resolution SAR correspond to accounts based on a high spatial resolution (2 m) land cover map over a study area of about 7 x 7 km in NE Siberia. The approach can 10 be also potentially transferred to medium resolution C-band SAR data such as ENVISAT ASAR Wide Swath with 120 m resolution but it is in general limited to regions without woody vegetation. Comparisons to the length of unfrozen period indicates the suitability of this parameter for modelling of the spatial distribution of soil organic carbon storage. 1 Biogeosciences Discuss., doi:10.5194/bg-2016-157, 2016 Manuscript under review for journal Biogeosciences Published: 3 May 2016 c © Author(s) 2016. CC-BY 3.0 License.

1 Introduction Figure 1. Location of field sites with high resolution land cover based SOC maps and area covered by the CAVM map as well as ENVISAT ASAR GM data (medium grey) with~120 m (Closa et al., 2003). Data availability of this mode is lower since it was acquired on request only. They are commonly gridded to 75 m x 75 m (e.g. Bartsch et al., 2007;Santoro et al., 2011;Reschke et al., 2012).

Soil organic carbon data
In-situ measurements of SOC used in this study have been collected from five different sites across the arctic: Shalaurovo, Kytalyk and Arymas in Siberia, Zackenberg on Greenland and Tulemalu in Canada (Fig. 1). All are located within the tundra 5 biome and are characterized by continuous permafrost terrain. For these sites, the investigated SOC depth increments are 0-5cm, 0-30 cm and 0-100cm. Two types of input data are used. Soil pedon point data and maps of SOC derived from thematic up-scaling of the soil pedon data using high to very high resolution optical satellite and airborne data (Palmtag et al., 2015;Siewert et al., 2015;Hugelius et al., 2010). In order to obtain these maps collected soil pedons were grouped according to the thematic classes in these schemes. Simple arithmetic means and standard deviations were then calculated for each SOC 10 storage depth increment per thematic class (and for the calculation of 0-30 cm and 0-100 cm SOC stocks). These means were subsequently weighed by the proportional representation of each thematic class in the study area to arrive to a weighed 4 Biogeosciences Discuss., doi:10.5194/bg-2016-157, 2016 Manuscript under review for journal Biogeosciences Published: 3 May 2016 c Author(s) 2016. CC-BY 3.0 License. almost 0 kg m -2 at alpine and barren ground locations to more than 80 kg m -2 for peat bogs. The maximum of non peat sites is approximately 35 kg m -2 for 100 cm and 15 kg m -2 for 30 cm. Table 1 provides further details and the data sources and land cover map thematic content.
The Northern Circumpolar Soil Carbon Database (NCSCD) by Hugelius et al. (2013Hugelius et al. ( , 2014 provides SOC stocks in the circumpolar permafrost region. The NCSCD is a polygon-based digital database compiled from harmonized regional soil 10 classification maps in which data on soils have been linked to pedon data from the northern permafrost regions to calculate SOC content and mass. It includes SOC values for 0-30 cm, 0-100 cm, 0 -200 cm and 0-300 cm. For this study, only the NCSCD area North of the arctic treeline as defined in the CAVM (Walker et al., 2002) is considered.

15
Radar backscatter is dependent on sensor parameters such as incidence angle, polarisation and wavelength as well as target parameters like surface roughness and vegetation structure as well as dielectric properties (Ulaby et al., 1982). Roughness and permittivity are the governing factors in case of bare soil (Oh et al., 1992). The dielectric constant highly depends on moisture content, leading to higher backscatter values in the microwave range under wet soil conditions (Woodhouse, 2006). Regions with soils conditions close to saturation near the surface can be therefore identified using C-band SAR data. This has been 20 demonstrated applicable for peatland detection at high latitudes (Reschke et al., 2012). The wet and at the same time high SOC areas have a low bulk density over several tens of cm and are water/ice rich (more than 60% at e.g. Kytalyk, Weiss et al.
( 2015)). The dielectric constant is significantly lower under frozen conditions. Frozen soils cause therefore similar backscatter like dry soils which has been specifically exploited for C-band applications (e.g. Wagner et al., 1999;Park et al., 2011).  backscatter is thus determined by the above surface remains of vascular plants, surface roughness, near surface soil texture and, if present, also snow cover. The latter has little influence in early winter (Naeimi et al., 2012). Interaction of C-band signals with snow are lower than for shorter wavelengths (Ulaby and Stiles, 1981). The signal may also penetrate a few cm into the soil. In undisturbed environments (no buildings or agriculture) it can be assumed that scattering is governed by soil type and vegetation cover. The influence of vascular plants on signal interaction is however limited at C-Band (approximately 5.6 cm 5 wavelength, Waring et al. (1995)). Surface roughness thus plays an important role for spatial differences in backscatter during frozen conditions in tundra regions. Specifically data acquired in HH (horizontally sent and received) polarization are expected to represent soil conditions better than VV (vertically sent and received) polarization (Brown et al., 2003). Vertically (with respect to the Earth surface) polarized waves interact more with vertically structured vegetation parts (stems) than horizontally polarized waves. HH as well as HV polarizations are thus more sensitive to roughness than VV polarizations (Holah et al., It has been shown for C-band (Jagdhuber et al., 2014) that volume scattering (at anisotropic particles) dominates for peatland soil during unfrozen conditions and it changes to surface scattering when frozen. The dielectric contrasts between scattering components decrease and surface roughness indeed determines the magnitude of backscatter. Wetter, carbon rich soils have a smoother surface (with respect to the 5.6 cm wavelength) than drier low carbon sites in the high arctic what leads to the 15 hypothesis that C-band backscatter can be used as proxy for SOC content (Fig. 2). In order to exclude the influence of spatial variations in soil moisture as well as deep snow cover only acquisitions from December (frozen soil, limited snow cover) are used for this study. Radiopositioning Integrated by Satellite) orbit files) and digital elevation data (Shuttle Radar Topography Mission-improved U.S. Geological Survey GTOPO30 digital elevation model), geocoded images were produced with sub pixel accuracy (Park et al., 2011;Pathe et al., 2009). The data are resampled into a fixed 15 arc-seconds grid (datum WGS-84), within 0.5°by 0.5°t iles, to allow efficient spatial and temporal analysis. The data (> 8000 scenes north of 60°N) was normalised to a reference angle of 30°by fitting a linear model to the backscatter data (Pathe et al., 2009;Sabel et al., 2012) in order to remove the 10 influence of local incidence angle on radar backscatter. The model provides an estimate of the slope in units of decibel per degree of incidence angle which characterizes the decrease of the radar backscatter from near range to far range. The model is calibrated for each pixel separately using the acquisitions from overlapping orbits (Wagner et al., 2008). No data can be processed with the used orthorectification tools (SGRT) for scenes which cross the dateline with this toolbox. This leads to a data gap in far Eastern Russia. The dataset has been eventually resampled to a grid with polar stereographic projection with 15 500 x 500 m pixel size.
On average 45 acquisitions have been available per pixel. Since GM data exhibit a comparably high noise (Park et al., 2011) temporal and/or spatial statistical measures (averaging, filtering etc.) need to be applied. The minimum of the entire record for each pixel was therefore used in this study instead of single values representing a certain date. Data are derived as sigma nought and converted to dB. The dataset has been masked for lakes and glaciers based on the map classes of the Global Lakes and 20 Wetlands Database (Lehner and Döll, 2004) and GlobCover (Bicheron et al., 2008) as well as tree line (Walker et al., 2002).

Determination of relationship between backscatter and SOC
The C-band backscatter is directly compared to locally up-scaled SOC maps and underlying pedon (point) data (table 1, Fig. 1).
Not all classes and the full range of SOC values can be found at single sites. A region with lower SOC (Zackenberg) and a site with medium to high values (Kytalyk) are therefore used in combination to obtain a representative range for the establishment 25 of the empirical relationship for up-scaling. The maps of the remaining sites have been used for validation.
Zonal mean values (a zone refers to a land cover class) have been extracted for the SOC classes available for Kytalyk and Zackenberg for model calibration. The advantage of the zonal mean opposed to the pedon data is the scale comparability to the GM data.
The Pearson correlation has been derived for the zonal means and eventually a function determined by least-squares re-30 gression. The obtained function has been subsequently applied to the circumpolar dataset. The land cover based SOC maps available from Tulemalu, Arymas and Shalaurovo have been used for validation. Regional differences have been assessed using the soil map based NCSCD (v2.2) by Hugelius et al. (2013). It has been converted to a 500 x 500 m gridded dataset with sep- arate layers for each SOC class and percentages of the soil types turbel, histel and histosols. This dataset has not been applied for training since its based on different types of data sources around the Arctic. The impact of soil type on the SOC retrieval is in addition investigated using the information available from the original pedon (point) information from all study sites since this information is not preserved in the land cover classifications. 3.5 Transfer of the approach to WS data

15
The higher spatial resolution (but lower sampling rate and inconsistent coverage) data from ENVISAT ASAR Wide Swath (also HH polarization) have been used in order to test the transferability of the approach across scales for the Kytalyk study site in NE Siberia. Due to the limited data availability, normalization cannot be applied as for the GM data (approach by Wagner et al., 2008;Sabel et al., 2012). Sabel et al. (2012) andWagner et al. (2008) exploit the availability of a representative range of incidence angles for a certain location by using acquisitions from several overlapping orbits. A conventional method which 20 corrects for local terrain related effects only (as available with the free NEST toolbox by the European Space Agency) has been used instead. The radiometric normalization available with NEST only accounts for terrain effects. This leads to a location specific bias with respect to the circumpolar GM based dataset. It was therefore required to adjust the WS data to the value range of the GM dataset. The incidence angle differs however by less than 0.1 • across the Kytalyk and Zackenberg sites. A single offset value per site can be therefore used to adjust the WS backscatter to GM. It has been derived from the average 25 regional backscatter of both datasets.

SOC determination
The range of dB of the GM data which corresponds to the SOC values represented in the reference datasets is about 8 dB for land cover class averages and almost 10 dB for pedon data. An R² of 0.86 was determined for the linear relationship between   (Table 2). SOC is lower at the latter site what agrees with the in situ records.

Comparison with independent datasets
Large negative deviations of more than 10 kg m -2 from the high resolution land cover maps (Table 1) are only found for peat bogs (Fig. 8) which are located at Tulemalu. This is also consistent with the pedon derived information for soil types (Fig. 9).  The mean difference between the NCSCD and the GM result is 3.8 kg m -2 and 5.8 kg m -2 for 30 and 100 cm respectively (standard deviation 6.3 kg m -2 and 15.1 kg m -2 ). SOC totals within the CAVM domain are listed in Table 3.
The differences increase with increasing SOC in the NCSCD (Fig. 6). SOC values from GM for both 30 and 100 cm are mostly higher across Northern America and lower across Siberia (Fig. 4). Transitions between areas of positive and negative value regions are sharp reflecting boundaries of maps which underlie the NCSCD (Fig. 7). SOC values change at country  The residual plots (Fig. 10) for the depth of organic layer and cryoturbated carbon also confirm that the SAR method is biased low in sites with substantial cryoturbation and deep O-horizons. The SAR method is biased high for sites with limited cryoturbation and/or less than 10 cm organic layer thickness. The differences do not relate to cryoturbated carbon in case of 30 cm estimates (R²=0.14), but to some extent for 100 cm values (R²=0.5).
The results obtained from the NCSCD for the NDVI classes suggest a bi-modal behaviour with the first maximum for NDVI  GM and NCSCD averages for the length of unfrozen period classes (Fig. 9) differ from each other. Maximum SOC in GM corresponds to about 110-120 days of unfrozen period length. A local maximum can be also found for NCSCD over that period but SOC is higher for more than 150 days.

5
In case of the NCSCD as well as the GM records, an increase of SOC with increasing length of the unfrozen period can be shown (Fig. 9). The variability increases for unfrozen period lengths over 120 days due to substantially low number of samples.   Samples correspond to pixels of 12.5 x 12.5 km Quickbird based results (Fig. 11, Table 3). The river flood plain shows lower SOC than the thermokarst landscape to the north.
Drained lake basin patterns and associated gradients are still captured with WS but not with GM. These differences are not captured in the NCSCD database. NCSCD SOC 0-30 cm over the complete area is 17 kg m -2 what adds up to 0.81 TgC.
The same satellite data source as in this study has been used for discrimination of wetness levels by Widhalm et al. (2015).
Areas with backscatter below -16.5 dB have been shown to correspond to wet areas with potentially higher methane emissions.

5
These areas would correspond to SOC 0-100 cm larger than 25 kg m -2 . This is confirmed by the fen and moist tundra records from the reference datasets (Table 1, Fig. 8). Medium or mixed wetness corresponds to about 20-25 kg m -2 . Several willow, grass, fen as well as dry tundra samples fall into this category. Dry areas as defined in Widhalm et al. (2015) correspond to mostly fell, heath and boulder classes. The majority of pedon records from histels and non-permafrost mineral soils can be also found in this category (see Fig. 12). Turbels and orthels can be associated with the mixed and wet classes with about 50% of 10 the turbels and 25% of the orthels in the wet class based on the pedon data. Orthels are also represented with about 50% in the wet class using the GM quantification.