HOW MUCH CARBON IS STORED IN DESERTS ? AN APPROACH FOR THE CHILEAN ATACAMA DESERT USING LANDSAT-8 PRODUCTS

The Atacama Desert in northern Chile is known as the driest place on Earth, with an average rainfall of about 15 mm per year. Despite these conditions, it contains a rich variety of flora with hundreds of species characterised by their extraordinary ability to adapt to this extreme environment. These biotic components have a direct link to important ecosystem services, especially those related to carbon storage and sequestration. No quantitative assessment is currently available for these services and the role of the desert in this matter remains unclear. We propose an approach to estimate above-ground biomass (AGB) using Landsat-8 data, which we tested in the Taparacá region, located in the northern section of the desert. To calibrate and validate the models, we used field data from 86 plots and several spectral indexes (NDVI, EVI and SAVI) obtained from the provisional Landsat-8 Surfacereflectance products. We applied randomised branch sampling and allometry principles (non-destructive methods) to collect biomass samples for all plant biological types: wetlands, steppes, shrubs and trees. All samples were dried in an oven until they reached constant weight and the final values were used to extrapolate dry matter content (AGB) to each plot in terms of kg m-2. We used all available scenes from September 2014 to August 2015 to calculate the maximum, minimum and average value for each index in each pixel within this period. For modeling, we used the method based on classification and regression trees called random forest (RF), available in the statistical software R-Project. The explained variance obtained by the RF algorithm was around 80-85%, and it improved when a wetland vector layer was used as the predictive factor in the model to reach the range 85-90%. The mean error was 1.45 kg m-2 of dry matter. The best model was obtained using the maximum and mean values of SAVI and EVI indexes. We were able to estimate total biomass storage of around 8 million tons (~ 4 million tons of C) for the whole region. * Corresponding author


INTRODUCTION
The Atacama Desert in northern Chile has been described as the driest place on earth (Lettau, 1978).Despite these conditions, it contains nearly 550 species of vascular plants representing 225 genera and 80 families (Anderson et al., 1990).The most diverse families are the Asteraceae, Nolanaceae, Cactaceae, Boraginaceae, and Apiaceae (Aronson, 1990;Benoit-C, 1989).These species of flora are characterised by their extraordinary ability to adapt to this extreme environment.These biotic components have a direct link to important ecosystem services, particularly those related to carbon storage and sequestration.Deserts can provide benefits that can meet the local inhabitants' demands (De Groot et al., 2012).These benefits include water, food supply, medicine from plants, carbon storage and sequestration and some raw materials.However, there are few studies about these services and how can they be quantified (Bidak et al, 2015).No quantitative assessment is currently available for these services and the role of the Atacama desert in this matter remains unclear.Taking into consideration the described satiation, we propose an approach to estimate aboveground biomass (AGB) using Landsat-8 data in the northern section of the desert.

Study Area
The Atacama Desert is an elongated strip of land along the northwest coast of Chile and far southwest of Peru, essentially bounded by the Pacific Ocean on the west.It extends nearly 1600 km and reaches a maximum width of 180 km (Hogan & Fund, 2014).Mean annual precipitation is about 15 mm and values below 5 mm are not uncommon.Most of this amount is due to drizzle from coastal stratus and very unusual rainfall episodes (Garreaud and Rutllant, 1996;Vargas et al., 2006).The study area is a section of the desert named the Tarapacá Region (Figure 1), an administrative division of Chilean territory of 42,226 km² (20°17′00″S 69°20′00″O).

Field data collection and allometry
Sampling was conducted in the field to collect biomass data covering main vegetation and land cover types in the region: grasslands ("pajonales"), shrublands ("tolares"), wetlands, bare soil and tree formations ("queñoa" and "tamarugo").A total of 86 samples were registered using concentric circular plots.Trees, shrubs and herbs were counted and measured in rings with 10 m, 5 m and 1 m radii for each plot, respectively.Individual biomass values were obtained by using allometry principles and randomised branch sampling (non-destructive methods) to collect biomass samples for all plant biological types using individuals of mean size: wetlands, steppes, shrubs and trees (Gregoire et al., 1995).All samples were dried in an oven until they reached constant weight and the final values were used to extrapolate AGB to each plot in terms of kg m-2 (Figure 2).Biomass values in terms of kg m-1 were obtained by extrapolating biomass branch samples to full individuals and from individuals to plots.

Landsat Data and Pre-processing
We considered all available scenes, path 001-002 and row 073-075 (WRS), from September 2014 to August 2015 of provisional Landsat-8 Surface-reflectance products and downloaded NDVI, EVI and SAVI spectral indexes (USGS, 2015).Additionally, we used MOD13 NDVI and EVI products (MODIS-TERRA), 250 m resolution, for the same period of time.
For all indexes we calculated the maximum, minimum and average value for each index in each pixel within this period.Using coordinates from each plot, we extracted all indexed values from both sensors.

Ancillary data
We also gathered some available ancillary data that provided additional explanation power for the biomass (Table 1).All variables were re-sampled to 30 m resolution.

Modelling approach
We used the method based on classification and regression trees called random forest (RF), available in the statistical software R-Project.This algorithm, proposed by Breiman (2001), constructs a classification or regression tree using a process called bagging.RF calculates the relevance of predictor variables by assigning a score that depends on the changes in the error when the values of this variable are permuted (%IncMSE) and the increase in error when each variable is removed (IncNodePurity).If the effect is larger, the variable is assigned greater importance (Reif et al., 2006).
Considering this approach, AGB of each plot was modelled using all spectral indexes (Landsat and MODIS) and ancillary data (Table 1) as predictors.

RESULTS
Figure 3 shows a graph of the variance of AGB per vegetation type.It is noted that, on average, deserts and salt flats have minimum biomass values (~ 0).By contrast, the wetland formation "Bofedal-Pajonal" type (BOPA) has the highest value (0.65 kg m-2), followed by Queñoal, a tree-dominated vegetation type that presents the greatest variance between samples.
Figure 4 shows the importance of each variable obtained by running the initial RF algorithm.Letter "L" is used to indicate Landsat products, while "Mo" is reserved for MODIS ones."Max", "Min" and "Mean" are used to indicate maximum, minimum and mean values in the period.In a second RF model we only used the first five variables according to the IncNodePrurity index: Vegas : Binary raster identifying wetlands LMaxSavi : Landsat-8 SAVI maximum for the period.LMaxEvi : Landsat-8 EVI maximum for the period.LMeanSavi : Landsat-8 SAVI average for the period.
And the final selected RF model is:

AGB=1+Vegas+LMaxSavi+LMaxEvi+LMeanEvi+LMeanSavi+ϵ
This model explains 87.09% of the variance of the data and exhibited a RMSE of 1.4529 kg m-2 (bias 0.0191).Figure 5 shows how the residuals (difference between observed and predicted) are concentrated in the centre of the histogram, indicating a non-biased model fit.
By using the final RF model we produced a predictive AGB for the complete Tarapacá Region (Figure 6).This product was used to estimate a total AGB accumulation of nearly 5.5 million tons.
Finally, by using the TESSA approach (Toolkit for Ecosystem Service Site-Based Assessment -http://tessa.tools/)we were able to estimate underground biomass using AGB values as predictors for each vegetation class.Final estimations are presented in

CONCLUSIONS
We developed a method to quantify the accumulated AGB in a section of Atacama Desert, in the Tarapacá Region.The absence of reference data made it necessary to collect new field data.
The results indicate that the greatest accumulation of biomass is, as expected, in the upper area of the region.In particular, it draws attention to the elevated accumulation in the highland marshes (wetlands) and surrounding vegetation.These systems have an important role in storing and filtering water in high altitude environments, mainly by the type of rhizomatous vegetation that coexists with the aquatic grasslands that form layers of organic matter in constant aerobic decomposition, forming a mattress of layers that enable water purification.
Overall, it may be said that this region, with its 4,222,580 ha, has accumulated a total biomass of around 8 million tons, i.e., about four million tons of C.These results can be incorporated into further assessments of ecological services of the desert.

Figure 3 :
Figure 3: Boxplots of AGB in each vegetation type or land cover.Bofedal, BOPA and Salar are different kind of wetlands; Desierto means bear soil.Matorral, DETO, QUPA and Tolar are types of shrubland.Pajonal is a type of grassland; Queñoal and Tamarugo are vegetation dominated by trees.

Figure 4 .
Figure 4. RF IncNodePurity indicates the amount the error increase following the removal of each variable.

Figure 5 .
Figure 5. Histogram of the selected RF model residuals.

Figure 6 .
Figure 6.Predicted ABG for the whole Tarapacá Region.

Table 2 .
AGB and BGB estimation for Tarapacá Region.