Soil properties after forest rehabilitation by planting teak and mahogany in Java, Indonesia

Abstract We studied how rehabilitation of forests in Indonesia by planting teak, Tectona grandis (L.f.), and mahogany, Swietenia macrophylla (King), was associated with soil pH, organic matter, total nitrogen, total phosphorus, and total potassium. We also analyzed how soil properties and the environment (i.e. soil order, altitude, stand age) were associated with succession and compared rehabilitated stands with native forests. We found higher pH in teak compared to mahogany stands. The soil pH was lowest in the oldest stands (>70 years). Herb density was positively related to pH and to phosphorus, while density of seedlings and woody plants was positively related to nitrogen, potassium and phosphorus. Tree and herb species richness and tree density were positively associated with Oxisols, but negatively related to the proportion of native herbs. Species richness of herbs and density of seedlings decreased with time since rehabilitation, whereas species richness of woody plants increased. The proportion of native herbs and seedlings increased with stand age. We found few differences in soils between the planted stands and native forest. Our results demonstrated that successional vegetation of rehabilitated forests may play an important role in maintaining soil properties associated with soil order.


Introduction
Supporting ecosystem services, such as nutrient cycling and soil formation, are essential for self-perpetuating ecosystems. We can describe these services by soil properties, such as the concentrations of nutrients, organic matter, minerals, physical properties and living organisms (Schoenholtz et al. 2000). Such properties can differ between soil types and ecosystems, and with species diversity (Huston 1980). Restoration of ecosystems by planting trees after deforestation may affect soil properties through various biological processes during succession (Jia et al. 2005).
Secondary vegetation succession is the development of vegetation after disturbance, including interactive changes in soil chemical, physical, and biological properties in plant communities, from pioneer species towards more mature "climax" species (Connell and Slatyer 1977). When planted trees grow, shading will increase and herbaceous vegetation might be shaded out, while non-pioneer woody species will increase, if they can compete with the planted species.
The interaction between plants and soil can affect both biotic and abiotic elements of soil, which influence the plant community (Wardle et al. 2004;van de Voorde et al. 2011). In natural ecosystems, trees and soils interact through facilitation of soil biota (i.e. trees affect the environment for soil organisms, stimulate microbial activity, and contribute to nutrient inputs) and synergistic symbiosis (i.e. nitrogen fixing bacteria and leguminous trees) (Barrios 2007; Barrios et al. 2012). Soil fertility seems to improve with increasing tree species diversity (Huston 1980;Long et al. 2012). Soil organisms decompose organic materials, resulting in the release of CO 2 and synthesis of soil organic matter (Barrios 2007). Essential soil nutrients (i.e. nitrogen, phosphorus) are derived from the mineralization of soil organic matter through activity of soil microorganisms, which may increase plant growth (Dijkstra et al. 2006).
Massive deforestation in Java, Indonesia, took place during colonial times in the 1700s and onwards (Whitten et al. 1996). More recently, around 59 million ha of Indonesia's forests were lost between 1950 and 1997 (Tsujino et al. 2016). Loss of forests continued by about 1.5 million ha per year during (Forest Watch Indonesia 2011. The government of Indonesia has been rehabilitating degraded forestland since independence after the Second World War, mainly by planting teak Tectona grandis (L.f.) or mahogany Swietenia macrophylla (King) (Nawir et al. 2007). Here, we examined soil properties within stands of teak and mahogany plantations of various ages. The soil chemical properties we investigated were soil pH, soil organic matter (SOM), total nitrogen (N), total phosphorus (P), and total potassium (K). We tested how soil properties related to quantitative descriptors of plant communities, i.e. species richness, density of plants, proportion of native species (non-exotic), stand type (i.e. tree species planted) and stand age.
We hypothesized that the concentration of SOM, N, and K in the topsoil would increase as a result of succession, while pH and P concentration would decrease due to biological processes during succession (Aerts and Chapin 1999;Jia et al. 2005;Long et al. 2012). In addition, we compared soil properties in rehabilitated stands (mahogany and teak) with those in undisturbed native forests. We expected that soils under the native forests would have higher nutrient concentrations than the rehabilitated stands (Amponsah and Meyer 2000).

Study area
The fieldwork was conducted in May-June 2015 in the Yogyakarta province, Java Island, Indonesia, (between 110 24 0 19 00 -110 28 0 53 00 E and 7 15 0 24 00 -7 49 0 26 00 S; Figure 1), in the state forests of Gunungkidul, Bantul, and Kulonprogo regency. The Yogyakarta area has a humid tropical climate with an average humidity of 86%, and average temperature of 27 C. The mean monthly minimum and maximum temperatures are 23 C and 33 C. Average monthly rainfall is 255 mm; the rainy season lasts from October to April, and the dry season from May to September. Large parts of the area experience water shortages during the dry season. The Yogyakarta province is 3186 km 2 , the population is 3.6 million, and the forest area is 187 km 2 (6% land cover). The topography is flat to undulating, with an altitude of 100-500 m above sea level. Limestone and barren karst dominate in the Gunungkidul area (Statistics of Yogyakarta Province 2017). Modern volcanic and alluvial deposits scattered in middle miocene reworked volcanic deposits (Smyth et al. 2008) characterize the Yogyakarta region. Merapi is the nearest modern volcanic center with the last large eruption in 2010.
After deforestation in Java during colonial times, some areas were left barren. Various cultivation attempts by the colonial powers, such as coffee plantations, failed. Since Indonesia's independence in 1945, the government has gradually rehabilitated the land by planting primarily teak or mahogany ( Santoso 2012;Balai KPH Yogyakarta 2014 ). Grazing in these area is not allowed, with livestock usually being stall-fed. Intercropping systems were employed to enhance income of the local communities, but this system was limited by forest canopy closure. Food crops were normally cultivated for about 4-5 years after tree planting, after which the shade from the trees impeded further cultivation. Generally, chemical fertilizer was applied when cultivating food crops and at the time of planting of teak or mahogany (C. Udayana, pers. obs.).
Stands sampled in the present study were planted between 1941 and 2003, here referred to as teak stands and mahogany stands. Mahogany is native to Central and South America (Orwa et al. 2009) but is now widespread in tropical forests globally, including in Java, Indonesia. Teak occurs naturally in peninsular India, Myanmar, Thailand and Laos (Verhaegen et al. 2010). It is now more or less naturalized in Java (Pandey and Brown 2000).

Stand selection and characteristics
Stands were selected based on forest rehabilitation history (e.g. year of planting, tree species planted) as suggested by the Governmental Forestry Service and by the Faculty of Forestry, University of Gadjah Mada. We surveyed 24 stands in total planted with teak and mahogany. In addition, we found and surveyed three areas of remaining natural forest fragments, which had never been rehabilitated, here referred to as native forests (Supplementary Table 1; Figure 1), for comparison with planted forest stands. The small sample size of native forest stands was because there is very little native forest remaining in the area, due to the history of deforestation.
The planted stands are managed by the Governmental Forestry Service, the Faculty of Forestry, University of Gadjah Mada, and the Natural Resources Conservation Center (Yogyakarta). Mahogany and teak were typically planted at a 2 Â 4 m spacing and thinned at the age of 10 and 15 years. The trees are usually logged at the age of 35 years, depending on the forestry ministry's decision. The diameter at breast height (DBH; 1.3 m) of old mahogany stands of around 70 years was in the 30-50 cm range, and old teak stands were 30-40 cm DBH. Tree density in mahogany sites was on average about 20 trees per 100 m 2 , with an average ± SE (5 ± 0.6) tree species per 100 m 2 (including the planted species), in both in old and young stands. Teak stands had 15-20 trees per 100 m 2 in the young stands, decreasing to 5-10 in old stands, and an average of 3 ± 0.3 tree species. This can be compared to native forests that had, on average, 9 ± 0.8 tree species and a tree density of 27 ± 5.9 trees per 100 m 2 .
In the study area, soil orders were typically mediteran, latosol, rendzina, and inceptisol (Wanagama 1988;Yogyakarta 2014; Dinas and Perkebunan 2015; BPDASHL Serayu Opak Progo 2018). These soil orders were transferred to the USDA soil taxonomy (Soil Survey Staff 2014), and correspond to Alfisols, Oxisols, Mollisols, and Inceptisols, respectively. The different stand types were fairly uniform with respect to soil order, but with regional differences. In Gunungkidul, the soil is mainly Alfisols or Mollisols, while Oxisols are dominant in Bantul and Kulonprogo (Supplementary Table 1). Alfisols are moderately weathered with clear horizons, typically found under forest vegetation. Oxisols are old soils with a low natural fertility, dominated by iron oxides, quarts and weathered clay minerals, and they are common on sloping lands in the tropics/subtropics. Mollisols are fertile soils with a thick surface layer rich in organic matter and with a high base saturation, typically found under long-term grasslands. Inceptisols are young deposits that are slightly developed (Eswaran and Reich 2005).

Field procedures
We surveyed vegetation and took soil samples in three plot replicates randomly selected in each stand, approximately 100-300 m from each other, giving a total of 81 plots within 27 stands (Supplementary Table 1). We surveyed vegetation in the following categories; herbs (i.e. forbs, grasses, and ferns), seedlings (seedlings of woody species >1 cm and 50 cm high), woody plants (>0.5 m 2.5 m high), and trees (>2.5 m high). Each plot was 10 Â 10 m and this area was used for recording tree species. In the NW quarter of each plot, a 5 Â 5 m plot was used for sampling woody plants, and in the NW corner of that plot, a 1 Â 1 m plot was used for sampling seedlings and herbs. In each sampling plot we counted the number of plants of each species in the respective plant category. The plant species were identified with help of plant taxonomists of the Silviculture laboratory at the Faculty of Forestry, University of Gadjah Mada, Yogyakarta, and plant conservation service of Purwodadi Botanic Garden.
Soil samples were taken for chemical analysis from 0 to 15 cm depth, from each of the 81 plots, across all 27 forest stands. Each sample was made up of four sub-samples taken from about a meter inside each corner of the 10 Â 10 m tree plot. Before sampling, surface litter was removed. After sampling, the composite sample was air-dried, mixed, roughly sieved, and put into labeled plastic bags before being sent to the laboratory of Indonesian Agency for Agricultural Research and Development (Yogyakarta) for analyses.

Soil laboratory analyses
Soil samples were rolled and passed through a 2-mm sieve for analyses. Soil pH was determined in water with a 1:2.5 soil:water mixture (Van Reeuwijk 2002). Soil organic carbon was measured by the Walkley and Black (1934) method and soil organic matter (SOM) was obtained by multiplying percentage soil organic carbon by the Van Bemmelen factor of 1.724 (USDA 2004). Total N was analysed by the Kjeldahl method (Van Reeuwijk 2002). Total phosphorus (P) and total potassium (K) were determined by HCl 25% extraction (USDA 2004). Phosphorus concentration was measured by GENESYS TM 20 spectrophotometer. The 240 FS AS atomic absorption spectrophotometer was used to determine potassium concentration.

Data analyses
Our data were balanced with respect to the number of teak and mahogany stands and with respect to soil order found in the respective stand categories (Supplementary Table 1). We only had one site with Inceptisols and we therefore removed this site from further analyses. As region was highly confounded with soil order, we chose to use soil order in the statistical models instead of region, as this could give more biological information.
We used the five soil chemical properties (i.e. pH, SOM, N, P, and K), as response variables in linear mixed models. Log-transformation was applied for P and K to achieve a normal distribution. We first compared soil properties between rehabilitated and native forests. Then we compared soil properties between mahogany (n ¼ 11) and teak (n ¼ 12) stands among the rehabilitated sites. We used age as a predictor in the model, as age was only defined for rehabilitated stands. Predictor variables added to the model were stand type (mahogany/teak), age since planting, altitude, and soil order. We added site as a random intercept to the model. We performed a backwards selection procedure using the "drop1 command" (Zuur et al. 2009), by removing the least significant predictor until we had only significant (p < 0.050) components in the model.
To study the effect of soil properties, soil order, and stand age on the vegetation we used the following response variables: species richness, defined as the total number of species present; vegetation density, as the number of individual plants present; proportion of native species of all species; and DBH of trees. We used a Poisson error distribution and log link function when analyzing species richness and vegetation density. For the analysis of the proportion of native species, we used a binomial error distribution and logit link function. For analyzing the DBH, we used normal error distribution and identity link function (Crawley 2011). All the data analysis were carried out using R (version 3.4.2; R Development Core Team 2015).

Results
Comparing the rehabilitated stands (teak and mahogany) showed that mahogany stands had a lower average soil pH than teak stands (relative difference between teak and mahogany (for the following values, estimate ± SE is reported) 0.38 ± 0.19, Table 1; Figure 2(a)). In rehabilitated stands, pH decreased with age since rehabilitation (À0.01 ± 0.004, Table 1; Figure 2(b)). Soil pH showed no relationship with soil order or altitude (  Figure 2(c)). P also decreased with increasing altitude (À0.002 ± 0.001) but showed no relationship with stand type or age (Table 1). Total K responded in the opposite way, increasing with altitude (0.002 ± 9 Â 10 À4 , Table 1) but showed no relationship with age or soil order. The effect of stand type on total K was near significant (p ¼ 0.058; Table 1).
We found no clear differences in any of the soil properties between the native forests and the rehabilitated forests (Table 2), although the model showed a tendency towards a difference in soil pH and N due to soil orders (v 2 ¼ 5.52, df ¼ 2, p ¼ 0.063) and (v 2 ¼ 4.89, df ¼ 2, p ¼ 0.087), respectively. The effects of altitude and soil order on other soil properties were similar to those of rehabilitated stands and are therefore not repeated here.
Species richness of herbs was negatively related to stand age (À0.02 ± 0.004), whereas species richness of woody plants (0.01 ± 4 Â 10 À3 ) and trees (0.02 ± 5 Â 10 À3 ) were positively related to age (Table 3). Species richness of herbs and trees were also related to soil order, with the highest richness of herbs in Oxisols (relative difference in Oxisols vs Alfisols: 0.58 ± 0.19) and lowest in Mollisols (relative difference in Mollisols vs Alfisols: À0.34 ± 0.45, Table 3). For tree species, the highest richness was in Mollisols (relative difference in Mollisols versus Alfisols: 0.91 ± 0.33) and lowest in Alfisols (relative difference in Oxisols versus Alfisols: 0.58 ± 0.19, Table 3). There was no relationship between species richness of any plant groups and soil properties (Table 3). Species richness of seedlings did not show any relationship with any of the explanatory variables (Table 3).
Density of herbs was positively related to soil pH (1.36 ± 0.12) and to total P (0.02 ± 4 Â 10 À3 ) and negatively to SOM (À1.21 ± 0.04) and total K (À0.01 ± 5 Â 10 À3 ), but showed no relationship with N, stand age or soil order (Table 3). Density of seedlings was negatively related to the age of the stand (À0.02 ± 9 Â 10 À3 ) and soil pH (À0.58 ± 0.16) and positively to N (4.10 ± 0.82) and K (0.02 ± 8 Â 10 À3 ), but showed no relationship with SOM, P or soil order (Table 3). Density of woody plants was negatively related to soil pH (À0.29 ± 0.07) and positively to N (1.17 ± 0.30) and P (0.02 ± 2 Â 10 À3 ), but there was no relationship with SOM, K, age, or soil order (Table 3). Density of trees depended on soil order, with the highest density in Mollisols (relative difference in Mollisols versus Alfisols: 1.03 ± 0.30) and lowest in Alfisols (relative difference in Oxisols versus Alfisols: 0.28 ± 0.18). Density of trees was also negatively related to K (À0.01 ± 0.01), but showed no relationship with pH, N, SOM, P, or stand age (Table 3). The proportion of native herbs increased with increasing stand age (0.03 ± 8.6 Â 10 À3 ), and varied with soil order, being highest in Mollisols (relative difference in Mollisols versus Alfisols: 34.62 ± 17.73 Â 10 6 ) and lowest in Oxisols (relative difference in Oxisols versus Alfisols: À0.96 ± 0.40, Table 3). The proportion of native seedlings was positively related to stand age (0.02 ± 0.008), as was the DBH of the planted trees (0.55 ± 0.08, Table 3). The proportion of native plants and DBH showed no relationship with any soil properties, while the proportion of native woody plants and trees showed no relationship with any of the explanatory variables (Table 3).

Discussion
As expected, we found soil pH decreased with age in both rehabilitated stand types. Mahogany stands, which were somewhat older than teak stands, had a lower soil pH than teak stands. However, we did not find the expected effect of stand type or age on other soil properties, which is surprising in relation to earlier studies (Huston 1980;Long et al. 2012). Li et al. (2013) examined available nutrients and showed an effect of succession, but we examined the total concentration of nutrients and we could not see any effect of succession over a timespan of >70 years. We found no clear differences in any of the soil properties between native forests and the planted stands. This could indicate that soil parameters of rehabilitated stands approached those of native forests with time. We should however keep in mind the low number of native forest sites sampled.
We saw that soil pH decreased with increasing age of the planted stands, supporting Aweto (1981), Perumal et al. (2017), and Li et al. (2013) who observed a decline in soil pH over 10, 18, and 30 years, respectively, as result of a succession. We showed that this decline continued with age. Accumulation of litter on the soil surface may indirectly contribute to the soil acidity due to a high carbon concentration. Comparing the two tree species used for rehabilitation, mahogany stands had a lower pH than teak stands, after controlling for age (on average mahogany stands were older than teak stands, 50 ± 3 years versus 39 ± 4, respectively). In humid tropical climates, high temperatures and precipitation can decrease soil pH due to leaching Figure 2. Soil properties in rehabilitated stands estimated from linear mixed models showing (a) differences in pH between stand types, (b) relationship between soil pH and age since rehabilitation, (c) total P in relation to soil order. Shaded area in (b) and bars in (a) and (c) shows 95% confidence intervals.  (Shamshuddin and Daud 2011). Acidity of soil is also caused by mineralization of organic matter and N fixation by legumes (Fageria and Nascente 2014). Further, pH tended to be related to soil order, being highest in the fertile Mollisols and lowest in Oxisols. This could be because Oxisols, which is poor in K, Ca, and Mg, is highly weathered (Von Uexkull 1986).
The results showed that total P and total K responded differently to altitude, with P decreasing with increasing altitude and K increasing with increasing altitude, in both the analyses including native forest and only rehabilitated stands. The soil orders were to some extent sorted according to altitude, with Oxisols somewhat more common in the higher areas and Alfisols in the lower, possibly explaining part of this difference. Given that the land was sloping, P could also have been washed away in soil particles by surface runoff, thus being lower in the topsoil at higher altitude. Regarding total K, old and weathered soils often have lower K concentration than young soils, and especially young volcanic soils (Graham and Fox 1971). Presence of young or medium young volcanic rocks in higher areas in the region may explain our increase of K with altitude.
Species richness of herbs and tree seedlings decreased with increasing age of the stand, while the density of herbs was unaffected. Species richness of woody plants and trees increased, but density decreased with increasing age. From succession theory, we know that after disturbance of the soil surface by deforestation and restoration activities, herbs colonize the bare ground, and the establishment of non-native species can be promoted (Stoddard et al. 2011). In our study, the proportion of exotic herbs decreased with increasing stand age, while the proportion of native species increased. Martin et al. (2005) used the proportion of native species as an attribute to measure restoration success. Species richness of woody plants and trees increased as the stands matured, as many of these species are adapted to growing in a dense forest. As their size and the size of the planted trees increased, their density decreased due to competition for space, light, water and nutrients.
Our study showed that herb density was positively related to soil pH, whereas seedling density and woody species density were negatively correlated to pH, which agrees with their responses to the age of the stands. Herbs responded in the opposite way to seedlings and woody plants as a result of succession. We also saw that the density of herbs and of woody plants were positively related to total P but negatively to total K, while density of seedlings and of woody plants were positively related to total N. Our previous analyses showed that Leguminosae was the largest plant family found in the study area (Udayana et al. 2019). It shows that leguminous herbs, seedlings, and woody plants, are important components for biological nitrogen fixation, which may benefit the build-up of soil N (Barrios 2007).
We know that a high total concentration of nutrients does not necessarily mean that these are available for plants, but there is an equilibrium between the total and the plant available pools. Nutrients stimulate plant growth, but plants also maintain nutrients within the soil/plant system. As an example: a poor plant density would make the soil more prone to erosion, thus would lead to loss of topsoil and nutrient loss. Dense vegetation, with deep roots, would help nutrient circulation as fine roots decompose and new roots take up the nutrients (Schoenholtz et al. 2000;Healey and Gara 2003).
The negative association between plant density and K is surprising. There might be additional elements in the K rich soils masking the result of K. Volcanic deposits may contain high K, but also too high levels of certain trace elements. The negative effect of K could be caused by imbalances in the uptake of other elements (Mengel and Kirkby 1980), such as trace elements on which we have no data.

Conclusion
Ecosystem services, related to vegetation and soil may increase with age since rehabilitation, when trees, woody plants, and native herbs increase. We would have expected to see soil organic matter and nutrient reservoirs increase with the growing vegetation, but we got few such responses. Although the oldest stands were >70 years, this is still a short time for the reservoirs to change. On the other hand, changes in soil pH were detected, which decreased through time. The soil pH under mahogany stands was lower than teak stands. Soil nutrient elements were positively associated with density of herbs, seedlings, and woody plants showing that litter abundance of the understory may supply the topsoil with a number of nutrients, such as nitrogen and phosphorus through decomposition. The results also showed that species richness of trees and density of trees were positively associated with Oxisols which had the highest phosphorus concentration. It seems that trees are a key factor for improving nutrient-poor Oxisols.