SIMULATING CROP YIELD, SOIL NITROGEN, AND ORGANIC CARBON IN NO-TILLAGE CROP SEQUENCES IN A SUBTROPICAL CLIMATE IN BRAZIL

Brazil stands out worldwide for its high grain production in areas of direct sowing. The objective of this study was to simulate and assess the relationship of soil organic carbon content and nitrogen, crop yield, and biomass of two crop sequences under the no-tillage system in a subtropical region of São Paulo, Brazil, using CSM-CROPGRO-Soybean and CSM-CERES-Maize models. The modeling was carried out considering the meteorological conditions of Jaboticabal, SP, Brazil. The treatments consisted of combining two summer crops (maize and soybean) with maize cultivation as a winter crop. The average biomass and productivity for corn were 15594 kg ha −1 and 5996 kg ha −1 , respectively, and for soybeans they were 5905 kg ha −1 and 3441 kg ha respectively. For soil organic carbon and nitrogen, a small variation was observed between years, and in addition there was a decline in their levels after a year with low biomass production. In our study, the RMSE and MAPE values between the observed and simulated productivity by the model were 2.21 kg ha −1 and 44.24%, respectively. The analysis of main components for the cultivation of corn explained 83.9% of the variability, and for the cultivation of soy, 93.5%. Among the tested models, the CROPGRO was the one with the best accuracy.


INTRODUCTION
The no-tillage system (NTS) has been presented as an alternative to mitigate the emission of greenhouse gases (GHG) arising from agricultural practices (Lal, 2015;Bayer et al., 2016;Paustian et al., 2016;De Araújo Santos et al., 2019). It is considered a low carbon agriculture system, resulting in increases in soil carbon stocks after some years of its implementation (Lal, 2015;De Araújo Santos et al., 2019;Silva et al., 2019). Besides, the NTS is also characterized by improving the physical, chemical, and biological structure of the soil (Raphael et al., 2016;Rosolem et al., 2016;Calonego et al., 2017).
However, some studies have reported that in the NTS, the increase in organic matter will only be effective when a species that is efficient in biological nitrogen fixation process is incorporated into the crop rotation . It is estimated that for every 10 units of carbon sequestered in the soil, there is a need to immobilize one unit of nitrogen (Six et al., 2002).
Aiming at the feasibility of employing decision support systems in different regions and climates, computational models seem to be a valuable alternative to estimate the amount of carbon and organic nitrogen in the soil, depending on climatic conditions and cultural and soil management practices, in a given period (Weber et al., 2016). The decision support system for agrotechnology transfer (DSSAT) is widely used to simulate crop yield, development, and income. The residual component of soil organic matter (SOM) of the CENTURY model was incorporated into the DSSAT, allowing for simulations and for conducting long-term sustainability analyses (Liu et al., 2017).
The objective of this study was to simulate and assess the relationship of soil organic carbon content and nitrogen, crop yield and crop biomass of two sequencing crops (soybeans and maize) under the NTS in a subtropical region of São Paulo, Brazil, using the CSM-CROPGRO-Soybean and CSM-CERES-Maize models.

Location and characterization of the experimental area
This experiment was carried out in the municipality of Jaboticabal, SP, at the coordinates 21°15′22′′ S and 48°18′58′′ O at 550 m altitude, during July and August 2016 for the 2015/2016 agricultural year. The climate of the region, according to the classification of Thornthwaite (1948), is of the B1rB'4a', humid mesothermal, with little water deficiency, presenting an average annual temperature of 22.2°C with the average of the hottest month being over 22ºC and the average coldest month being above 18ºC. There is an average annual precipitation of 1.425 mm, with higher volumes from October to March (Rolim & Aparecido, 2016).
The soil was classified as eutrophic Red Latosol, clay texture (Santos et al., 2013). Since 2001, the area has been under the NTS and the results presented in this study correspond to the years 2003 to 2016. Before the implementation of the system, it was used for the production of soybeans and corn in a conventional soil tillage system for 25 years.
The treatments consisted of a combination of two sequences of summer crops with one winter crops. For the summer, there was either monoculture of maize (Zea mays L.) (MM) or soybean (Glycine max L.) (SS), and maize was the winter crop. For more information on the crop treatments and experimental design, consult Marcelo et al., (2009), De Araújo Santos et al., (2019 and Xavier et al., (2020).

Modeling procedures and input variables
The CROPGRO (Jones et al., 2001), CERES (Jones & Kiniry, 1986) and CENTURY (Parton et al., 1994) models were applied to simulate crop yield, crop biomass, and soil organic carbon and nitrogen from 2004 to 2016 for the cultivation of maize and soybeans. Site-specific input variables, such as soil texture (sand, silt, and clay contents), soil density, and SOC and SNT contents, and average annual productivity for soybeans and corn were extracted from published works (Marcelo, 2007;Marcelo et al., 2009;Marcelo, 2011;Martins et al., 2012) with results for the years 2003,2005,2006,2007,2009, and 2010 from the experimental area described in 2.1. The data from meteorological elements for Jaboticabal, SP, Brazil were used in the DSSAT program. The daily inputs were maximum temperature, minimum temperature, wind speed, relative humidity, precipitation, and global solar radiation. These variables were obtained from the Meteorological Station of the Faculdade de Ciências Agrárias e Veterinárias da Universidade Estadual Paulista "Júlio de Mesquita Filho" (FCAV-UNESP). The FAO 56 method was used to estimate crop evapotranspiration (ETc).

Data analysis and model evaluation metrics
The data were initially analyzed using descriptive statistics (mean, standard error of the mean, standard deviation, minimum, maximum, coefficient of variation, asymmetry, and kurtosis) of the simulated data. Subsequently, the percentage deviation (%) (Equation 1) of the estimated values from the observed values of productivity were established (kg ha −1 ): The performances of the models were evaluated using linear regression analysis, in which the independent variables were those of the observed data and those dependent on the results extracted from the DSSAT. The adjusted coefficient of determination (R 2 adj.) (Equation 2) was calculated according to Cornell & Berger (1987) and the square root of the mean error (RMSE) (Equation 3) and absolute percentage of the error (MAPE) (Equation 4) (Willmott, 1981): Where: N is the number of points in the data sample, K is the number of independent regressors, that is, the number of variables in the model, excluding the constant. Where: N is the number of points in the data sample; Yobsi is the observed value of Y, and Yesti is the estimated value of Y.
Principal component analysis is an exploratory multivariate technique that condenses the information contained in a set of original variables into a set of smaller dimensions, composed of new latent variables, preserving a relevant amount of the original information. The new variables are the eigenvectors (main components) generated by linear combinations of the original variables, constructed with the eigenvalues of the covariance matrix.
The main components whose eigenvalues were higher than the unit were considered according to the criterion established by Kaiser (1958). The coefficients of the linear functions, which define the main components, were used to interpret their meaning, using the sign and the relative size of the coefficients as an indication of the weight to be assigned to each variable. Only coefficients with high values were considered for interpretation, usually those greater than or equal to 0.50 in absolute value. These analyses were processed using R (R Development Core Team, 2017).

RESULTS AND DISCUSSION
The mean of maize biomass for all simulated years was 15594 ± 597 kg ha −1 , while for soybean it was 5905 ± 164 kg ha −1 ( Table 1). The crop yield was 5996 ± 275 and 3441 ± 121 kg ha −1 for maize and soybean, respectively (Table 1). According to the Companhia Nacional de Abastecimento-CONAB (2010), the forecast of maize yield for the 2010 harvest was 3906 kg ha −1 in Brazil, with a surplus of 2000 kg ha −1 compared to the DSSAT simulation. The soybean yield forecast from CONAB in 2010 was 2629 Engenharia Agrícola, Jaboticabal, v.40, n.4, p.536-544, jul./aug. 2020 kg ha −1 , with a value lower than the means found for the simulated experiment in the DSSAT. The values of soil organic carbon and nitrogen were the same for both crop sequences (Table 1). The values of the coefficient of variation of the organic carbon of the soil were low, thus demonstrating that the variation of this variable during the 13-year chronosequence was quite limited.
The use of maize in crop sequences is a practice highly recommended in the literature since maize has a high potential for biomass production when compared to soybeans (De Araújo . Furthermore, the physiological processes of maize cause the crop to produce a biomass rich in carbon (Taiz & Zeiger, 2010), resulting in a slower breakdown and hence increasing the soil carbon content (Martins et al., 2012).
Another reason for the possible relationship of the highest levels of soil organic carbon in areas with maize cultivation in the sequence of crops is that their roots are the main pathways for the release of plant exudates into the soil, which end up interfering with microbial activity that causes them to release exudates that enrich the soil with organic carbon (Austin et al., 2017;Faucon et al., 2017). TABLE 1. Descriptive statistics of crop biomass (kg ha −1 ) and crop yield (kg ha −1 ) for maize in the sequence maize-maize and soybean in the sequence soy-maize and soil surface total organic carbon at maturity (OCTAM, kg ha −1 ), soil organic carbon at maturity (OCAM, kg ha −1 ), soil surface total nitrogen at maturity (ONTAM, kg ha −1 ), and soil nitrogen at maturity (ONAM, kg ha −1 ). The results of the temporal variation of crop yield and weight showed similar behaviors for both soybean ( Figure 1A) and maize ( Figure 1B). Olibone et al., (2010) and Caires et al., (2015) observed that crop weight production indicates, in most cases, good yield values.
The levels of organic carbon decreased slightly after 2010 (Figure 1). This fall in SOC and SON levels occurred just after the years when crop yields were low. This indicates that the production of crop weight did not exceed that of previous years when compared to years with a higher crop biomass/yield. In this way, the microorganisms may have consumed the organic matter present in the area, causing the decrease. Such behavior was observed for the two sequences studied. It is worth mentioning that although there was a decrease in the levels of total organic carbon, this variation was very low (Table 1).
The levels of organic nitrogen, as well as those of carbon, also decreased as time went on, and this was already expected since as the levels of organic matter in the soil decreased, there was also a decrease in nitrogen (Batjes, 1996).
Engenharia Agrícola, Jaboticabal, v.40, n.4, p.536-544, jul./aug. 2020 FIGURE 1. Temporal variation of crop biomass (kg ha −1 ) and crop yield (kg ha −1 ) for maize in the sequence maize-maize and soybean in the sequence soy-maize and soil surface total organic carbon at maturity (OCTAM, kg ha −1 ), soil organic carbon at maturity (OCAM, kg ha −1 ), soil surface total nitrogen at maturity (ONTAM, kg ha −1 ), and soil nitrogen at maturity (ONAM, kg ha −1 ).   The model adjustment (R²adj. 0.91, p < 0.05) was significant only for soybean productivity (Figure 2). The adjusted coefficient of determination indicated an optimal fit of the model. With R²adj. 0.91, this means that 91% of the variation in soybean yield was explained by the model. Yang et al., (2014) found an R 2 of 0.96 between simulated and observed data for soybean yield in an area without irrigation, while in irrigated areas, Dogan et al., (2007) found values of R 2 of 0.94 and 0.88 for 2003 and 2004, respectively. Thus, it can be seen that the CROPGRO model can satisfactorily simulate soybean productivity under different management conditions in Jaboticabal, São Paulo, Brazil.
Still, for the performance of the model, it is worth noting that both RMSE and MAPE are good metrics to use for calibrating the model because the RMSE has the same unit of measurement as the simulated variables and the MAPE is given as a percentage (Yang et al., 2014). In our study, the RMSE and MAPE values between the observed and simulated productivity by the model were 2.21 kg ha −1 and 44.24%, respectively. Coelho et al., (2018), simulating sugarcane productivity for the municipality of Jaboticabal, found RMSE values for the cultivars studied ranged between 1.91 and 2.58 Mg ha −1 . For the cultivation of maize, the two factors together explained 83.9% of the total variation of the original data. The first process (PC 1) represents 49.7% of the variability and the second process (PC 2) 34.2%. The process contained in PC 1 is the most important for this study, as it is derived from the highest eigenvalue and has the highest percentage of explanation (49.7%), with the variables that most contribute to this being represented by weight (−0, 83), yield (−0.86), soil total organic carbon at maturity (OCTAM, kg ha −1 ), and soil organic carbon at maturity (OCAM, kg ha −1 ), (−0.86) (Figure 3). FIGURE 3. Biplot chart showing the biomass (kg ha −1 ), yield (kg ha −1 ) soil surface total organic carbon at maturity (OCTAM, kg ha −1 ), soil organic carbon at maturity (OCAM, kg ha −1 ), soil surface total nitrogen at maturity (ONTAM, kg ha −1 ), and soil nitrogen at maturity (ONAM, kg ha −1 ) for maize.
According to the signs of the factorial loads, PC 1 is negative and strongly correlated with yield, followed by OCTAM, OCAM, and biomass (Table 3). We can understand the relationship between nitrogen, productivity, and biomass in CP1 as being related to conservationist practices, such as the direct seeding system (NTS), which disturbs the soil less and applies a greater contribution of organic waste (Li et al., 2017), and this accumulation of organic matter in the soil directly influences the productivity of the crop. For CP 2, the factorial loads were also negative for ONTAM (−0.95) and ONAM (−0.95) ( Table 3). CP 2 is only being influenced by soil nitrogen. Nitrogen is the nutrient required in the greatest quantity by plants (Obour et al., 2017), but nitrogen absorption has very complex dynamics in the soil (Silva et al., 2005), which may justify the process of the second main component containing only nitrogen.
In the main component analysis for soybean cultivation, both factors explained 93.5% of the total variation in the original data. CP 1 represents 61.2%, and PC 2 32.3% (Figure 4). The process contained in CP 1, namely the variables that contribute the most to this, are represented by OCTAM (0.93), OCAM (0.96), ONTAM (0.95), and ONAM (0.95). The factorial loads are all negative (Table 4). PC 1 is represented by C and N. They are the main components of SOM and their stocks will vary depending on the rates of addition, particularly by waste. In agricultural systems, the stocks of organic C and N are also influenced by the management system adopted (Souza et al., 2009). PC 2 is represented by biomass (−0.96) and productivity (−0.98), and is also represented by negative factor loads. Productivity is influenced by the amount of biomass. These variables are favored by the vegetation cover formed due to the use of the no-till system (Muraishi et al., 2005).

CONCLUSIONS
Among the tested models, CROPGRO had the best accuracy. Therefore, DSSAT can be used to simulate soybean yield under the NTS for the climatic conditions of Jaboticabal/SP. The levels of carbon and organic nitrogen in the soil showed little variation over the years, but there was a small decrease after years with low biomass production. The use of multivariate techniques is a useful tool to verify the relationship of organic carbon and nitrogen in the soil with crop productivity when using simulated data.

ACKNOWLEDGMENTS
This study was financed in part by the Coordenação de Aperfeicoamento de Pessoal de Nível Superior (CAPES, Brazil) under the funding code 001 and the Conselho Nacional de Desenvolvimento Científico e Tecnológico (CNPq) (Processo 142384/2017-8).