Soil texture derived from topography in North-eastern Amazonia

ABSTRACT We present a 1:100,000 scale soil texture map of Paragominas county (Pará, Brazil), covering 19,330 km2. The method allows rapid production of a soil texture map of a large area where the strength of a duricrust controls the relief. It is based on an easily accessible explanatory variable, topography, which is represented using a Digital Elevation Model. The method makes it possible to map the spatial distribution of the texture of the topsoil layer. Modelling was complemented by field observations to identify the laws governing the spatial organisation of soil textures. The spatial variability of the elevation above sea-level of the duricrust was obtained by Kriging. The error rate of the resulting map is 26%, and the observations of the four soil texture units were respectively 78%, 90%, 41% and 60% accurately located.


Introduction
Soil fertility is an important driver explaining such land use dynamics as deforestation and crop expansion. Soil fertility in the humid tropics depends strongly on their texture (Quesada et al., 2009), because they determine erosion and compaction, cation exchange capacity and water holding capacity. For example in eastern Amazonia, as the cation exchange capacity is better in clayey soils than in sandy soils, the clayey soils are more suitable for agriculture and have been preferred areas for soybean and maize expansion since the 2000s (Piketty et al., 2015). Knowledge of soil properties such as soil texture over large areas is then indispensable when planning territorial development and environmental conservation. However, there is a general lack of soil data concerning tropical areas, particularly when they are forested.
Effective time saving approaches are thus needed to map soil properties. In tropical areas, it is necessary to find methodologies that rapidly and effectively capture information about the spatial variability of soils and reduce the need for intensive and expensive sampling (Mora-Vallejo, Claessens, Stoorvogel, & Heuvelink, 2008).
The current research concerns a region in northeastern Amazonia where the relief is controlled by a duricrust surface. Duricrust (a hard crust or nodular layer formed at or near the land surface by weathering processes (e.g. ferricrete, calcrete, silcrete)), is widespread in the tropics. Iron-rich duricrusts (also called petroplinthites, ferricretes or lateritic crusts) are made up of a hardened accumulation of alumina and/or iron oxide forming an indurated layer (International Union of Soil Science Working Group [IUSS WRB], 2006;Quesada et al., 2011). In a laterite profile of a soil known as plinthosol, duricrust overlies a clayey mottled zone and weathered zone (saprolite). Because the duricrust resists erosion, it is frequently exposed by erosion. Particularly in old low relief regions, the duricrust caps plateaus that are dissected by rivers. Erosion is slower on the plateaus because of the indurated layer (Eze, Udeigwe, & Meadows, 2014;Twidale & Bourne, 1998). The hypothesis tested here is that the duricrust can be mapped by analysing elevation data. The distribution of the duricrust can then be used to predict other layers inherited from the same pedogenic processes and located at either higher or lower elevations, especially when their texture is different.
The objective of this work is to map soil texture units using soil-landscape modelling in a county of Amazonia. A suitable methodology is proposed for the conditions that prevail in north-eastern Amazonia: (i) the geological structure is not too complex, consisting of monoclinal sedimentary units, (ii) geomorphic surfaces are inherited from lateritic pedogenesis, (iii) the soil texture is closely linked to the landscape morphology and (iv) the soil texture explains the differences in resistance to erosion. landscape morphology of this region is highly influenced by an indurated lateritic layer dating from the Paleogene, which forms tabular relief separated by wide valleys. The study is focused on the County of Paragominas which covers 19,330 km 2 .

Geology
The County of Paragominas is in the Grajaú geological basin (Figure 1). The Grajaú Basin covers 100,000 s km 2 extending from the Tocantins River to the southeast of São Luis do Maranhão. The Grajaú Basin is filled to a depth of 800 m by clastic sediments of Cretaceous to Paleogene age (Soares Junior, Costa, & Hasui, 2008). In Pará, the basin is represented by kaolinitic sandstone of the Itapecurú series (Albian-Cenomanian) and the Ipixuna formation (Upper Cretaceous-Lower Tertiary). These sets are capped by a duricrust which developed over a wide area of flattening during the Paleogene, called the 'Sulamericana Surface' (Rossetti, 2004). A clay layer, called 'Belterra clay', covers the duricrust to a thickness of around 10 m (King, 1962;Kotschoubey, Truckenbrodt, & Calaf, 2005b). During the late Oligocene, erosion led to major dissection of the relief and to the formation of wide valleys. Finally, during the Quaternary, new climatic and/or tectonic activity led to the formation of the existing drainage network (Kotschoubey, Truckenbrodt, & Hieronymous, 1996).
Around Paragominas, the common profile observed along the escarpments of the plateau, from bottom to top, consists of Kotschoubey, Calaf, Costa Lobato, Saba Leite, and Duarte Azevedo (2005a) ( Figure 2): (1) bed rock of friable sandstone and clayey sandstone; (2) mottled clays of irregular thickness varying from 10 to 30 m, the mottled black, red, orange and white colour is due to iron oxidation and reduction in the plinthic horizon; (3) duricrust (frequently 5 m thick) composed mainly of bauxite is sometimes massive, but usually with a nodular or columnar structure, and consists of gibbsite and hematite; (4) Belterra clay which is generally 10 m thick, but can reach 15 m or be absent; these yellow clays are composed of kaolinite.

Landforms
North-eastern Pará and western Maranhão are made up of dissected plateaus inclined towards the northwest. This morphological entity is named 'Planalto Setentrional Pará -Maranhão' (Barbosa & Pinto, 1974), and extends 400 km from north to south. The plateaus near the town of Paragominas have an elevation of 160-190 m above sea-level and separated by valleys up to several kilometres wide ( Figure 3).

Soils
The most prevalent soils are Ferralsols, called latossolos according to the Brazilian Agricultural Research Corporation (Embrapa) classification (Embrapa, 1981) or Oxisols according to the USDA classification (Soil Survey Staff, 2014). Ferrasols are deeply weathered, red or yellow soils of the humid tropics. They have diffuse horizon boundaries, a clay assemblage dominated by low-activity clays (mainly kaolinite) and a high sesquioxide content (International Union of Soil Science Working Group [IUSS WRB], 2014). In the study area, their textures are highly dependent on the nature of the bed rock and on the topography (Sombroek, 1966). On plateaus covered with Belterra clay, they are clayey Ferralsols (70-80% kaolinite) (Kotschoubey et al., 2005a). On the upper slopes of the valleys, they are gravel soils, formed by the dismantling of the duricrust. In the valley bottoms and lower slopes, they are loamy sand Ferralsols formed by weathering of sandstone (Rodrigues, Silva, Oliveira Junior, Gama, & Valente, 2003).

Climate
The climate in the study area is humid tropical (type Aw, according to the Köppen classification). A mean annual temperature of 26.3°C and annual precipitation of 1693 mm were recorded at Paragominas station between 1992 and 2010 (Andrade, 2011). Ninety per cent of annual precipitation falls in the wet season, which lasts from December to the end of May (Pinto et al., 2009).

Vegetation
The natural vegetation is a dense tropical moist forest, but 45% of the study area has been deforested since the 1960s, mainly for cattle farming, and since the 2000s for soybean and maize cropping (Piketty et al., 2015;Pinto et al., 2009).

Methods
Toposequences were studied by combining field observations and spatial analysis. The resulting soil-landscape units were projected onto a Digital Elevation Model (DEM) and then sampled in the field. The aim of the field work was to understand the overall relief and soil relationships and to conduct local surveys to identify the texture of the topsoil layer. Spatial analysis with a geographical information system (GIS) was used to characterise the spatial relationships between soil texture and topography, which enabled us to map out the Paleogene duricrust surface, the associated layers above and below the duricrust, and then to produce a soil texture map.

Field observations
Preliminary field work consisted mainly of a systematic description of the links between landform and soil texture. Four morphopedological units were defined: plateaus with clayey soils, plateau convex edges covered by gravel soils, steep slopes covered by mottled clay on the upper part and by loamy sand on the lower part, and valleys and plains with loamy sand. Observation points were selected along toposequences to represent the four morphopedological units in various locations of the study area, including outside the county border in order to perform the interpolation. Visual and tactile identification of the texture were performed, following the method described in Juo and Franzluebbers (2003): a small sample of soil is moistened, while kneading, the sample's stickiness (clay), smoothness (silt) and grittiness (sand) are evaluated. Then, the sample is assigned to a class in accordance with the USDA/FAO soil texture triangle. The soil colour was also used to distinguish Belterra clay and mottled clay: Belterra clay is yellow while mottled clay consists of black, red, orange and white spots, due to iron/aluminium oxidation and reduction. As Ferralsols have diffuse horizon boundaries, identification focused only on the topsoil layer which is representative of the deeper texture (Quesada et al., 2011). The time saved in comparison with digging auger holes made it possible to increase the number of points to obtain better spatial coverage of the different morphopedological units in the study area. A total of 417 points were selected, which represents a mean density of 1 point per 46 km 2 for the whole county. There are 107 points of clay, 158 of duricrust, 37 of mottled clay and 115 of sand. The points were located using a global positioning system (GPS) receiver and incorporated in a GIS.

Digital Elevation Model
The topography was represented by the Topodata DEM (Valeriano & Rossetti, 2012). Topodata was produced by INPE (Instituto Nacional de Pesquisas Espaciais) by resampling the Shuttle Radar Topography Mission (SRTM) 3 ′′ resolution DEM data, to obtain a DEM with a resolution of 1 ′′ , that is, about 30 m.

Soil map
A soil map of the county was compiled by Embrapa -Amazônia oriental (Rodrigues et al., 2003) at a scale of 1:250,000. The method used was based on a field survey plus photointerpretation of radar images (RADAM project). The soils were classified as Plintossolos, Gleissolos, Argissolos, Neossolos and Latossolos respectively Plinthosols, Gleysols, Acrisols, Arenosols and Ferralsols of the World Reference Base (International Union of Soil Science Working Group [IUSS WRB], 2014).
The method proposed in this article is distinguished from the Embrapa method by a low cost approach based on morphopedology and GIS analysis, without a soil classification objective but oriented specifically on the textural properties, a determining factor for agricultural expansion and the land planning.

Spatial analysis
The sub-horizontal structure and the generally regular thickness of the layers at a regional scale (Kotschoubey et al., 2005a) led us to assume that the relative soillandscape relationships are universal for this area, and that the vertical position and horizontal location of each layer could be derived from the elevation above sea-level of the Paleogene duricrust using the DEM.
The elevation asl of the observation points was deduced from overlay with the Topodata DEM ( Figure  4). Textural classes of clay (Belterra), duricrust, mottled clay and loamy sand were extracted to form four geographical data layers, with their elevation above sealevel as an item value. The surface of the Paleogene duricrust is a continuation of the southeast-northwest inclination, as previously revealed by the geological studies noted above. It occurs on the different plateaus and is interrupted by the dissecting valleys (Figure 3). The surface cannot be obtained directly from plateau elevation because it is partly covered by clay. However, it can be represented by a theoretical surface after interpolation between the observation points of the duricrust. The surface was built using a Kriging model.
Kriging is suitable for data with a spatial autocorrelation, which is the case of the elevation asl of a duricrust resulting from a flattening area. Kriging estimates the statistical relationships among the sample points. It assumes that the distance between sample points reflects a spatial correlation that can be used to explain data variation. The algorithm implemented in Esri ArcGIS ® first requires the creation of an empirical semivariogram representing the variance for each pair of observation points in the study area and secondly in the adjustment of a model to fit the semivariogram. Among the five models available, the spherical function was chosen using ordinary kriging. A variable radius was used to calculate each cell value from 12 neighbouring observation points. The modelled surface of the Paleogene duricrust was then used to infer the  location of the other texture layers by subtracting the elevations from the elevations of the DEM Topodata. A classification was applied to the results to compile the map of soil texture, where the values of classes representing the mean elevation differences between the texture units were derived from the profile reported in the literature and the field observations presented in Table 1.
The validity of the results was checked first against the vertical position of the observation points on each geological layer in relation to the modelled duricrust surface by spatial subtraction and used to test agreement between the model results and the regional geology. The validity of the results was also checked by the overlay of the soil texture map on the observation points and on the Embrapa soil map to calculate the respective true positive rates for each soil texture unit (Figure 4).

Results
The soil texture map (Main Map) shows plateaus covered by Belterra clay. The lateritic duricrust exposures at the convex edges of the plateaus are explained by the resistance of the duricrust to mechanical erosion. The upper slopes of the valleys consist of mottled clay (the lower part of a weathering profile), while the lower slopes and valley bottoms consist of loamy sand derived from Itapecurú and Ipixuna sandstone.

Vertical differences between observations and the interpolation surface
According to the geological relationships reported in Section 2.1.1., Belterra clay should be above the Paleogene duricrust, and have a thickness of up to 15 m, mottled clay should be beneath the duricrust and have a thickness ranging from 10 to 30 m, and the loamy sand should be below the mottled clay.
The modelled surface elevation asl was very close to the elevation of the duricrust at the observation points. 77% of the observed elevations of the Belterra clay were above the interpolation surface, which is in agreement with the regional profile in which the Belterra clay overlies the Paleogene duricrust; the median is +13 m ( Figure 5). Some negative values could be due to a mismatch between the surface interpolation at these points, but also to either errors in tactile and visual identification (existence of pure clay layers in the clayey sandstones or lighter coloured mottled clay can being confused with Belterra clay) or by errors in the calculated duricrust elevation resulting from Topodata DEM uncertainty. Eighty-five per cent of mottled clay observations lay beneath the interpolation surface of the duricrust with a −12 m median ( Figure 5), which corresponds to the actual depth of the mottled clay  Figure 4. Flow chart of field observations and spatial analysis used to compile the map of soil texture.
below the duricrust (mottled clay thickness is of the order of 10-30 m). Eighty-four per cent of loamy sand was located below the interpolation surface of the duricrust with a greater difference in elevation (−27 m median) than for mottled clay. This is in agreement with observations of the clayey sandstone of the Itapecurú or Ipixuna Formations ( Figure 5).

Spatial differences between observations and the soil texture map
The observation points were combined with the modelled soil texture map to validate the spatial dimensions of the data (Table 2). Assuming the observation points are accurate (some errors could have been made in tactile and visual identification) the error rate is 26%. For each true positive rate: 78% of Belterra clay and 90% of duricrust observations were accurately located ( Table 2). The good agreement of the duricrust map is linked to the modelling method, in which the interpolation surface was adjusted to fit the observation points. The results for mottled clay and loamy sand were less accurate, with a 41% and 60% true positive rate respectively.

Comparison with the Embrapa soil map
Overlaying the Embrapa map with the resulting soil texture map revealed a relatively good match of the texture type between the two maps ( Table 3): 89% of the Belterra clay area corresponded to clayey Ferralsols, 72% of the mottled clay area to clayey Ferralsols, 58% of loamy sand area to medium texture Ferralsols, but only 7% of the duricrust area overlaid pisolitic Ferralsols. The topographic position of the duricrust and its spatially narrow pattern could explain the lack of agreement between the maps.

Conclusions
In regions where the relief is strongly controlled by a duricrust, the method presented here makes it possible to cheaply and rapidly create a soil texture map. The method is especially suitable for large areas partly covered by forest where difficult access is a problem. The method was shown to be effective in estimating the spatial variability of soil textures in the study area, due to the regular vertical organisation of the different soil textures. However, the method does have certain limits linked to the uncertainty of the DEM, errors in the identification of textures at the observation points, uncertainties concerning the duricrust interpolation surface, the irregular thickness of the layers of mottled clay, and more generally, the difficulty involved in representing natural phenomena with complex patterns.
Software ESRI ArcGIS 10.1. was used for all analysis and map production.