AquaCropR: Crop Growth Model for R

Abstract

The Food and Agriculture Organization (FAO) AquaCrop model, run either via a standalone graphical user interface (GUI) or via a matlab application programming interface (API) (AquaCrop-OS), has been successfully tested on many crop species and under multiple scenarios. However, with these current versions, it is difficult for users to adapt formulae, add functionality or incorporate the model into other applications such as decision support tools. Here, we report on the release of a version of AquaCrop written in R. Performance of the model was tested using published datasets of wheat (Triticum aestivum L.) and maize (Zea mays L.), comparing output from AquaCropR with these other versions of AquaCrop. Our goal in developing this version was to widen the use and improvement of AquaCrop through open access.

1. Introduction

The negative impact of climate change on agricultural production and trade, combined with a growing world population, are great concerns shared by farmers and decision-makers [1]. In the context of changing environments, crop models are useful tools to simulate crop yields. Crop models have been evolving to translate processes related to crop growth and development into mathematical equations. Crop models have been used over decades for agronomic purposes and have more recently incorporated advances in the modelling of environmental footprints, biotic constraints, trait and gene effects and the upscaling of global change impacts [2]. Crop models can be incorporated into tools for decision support, which can help farmers to manage crops and take corrective measures towards improving crop performance. Therefore, adoption of climate-smart agriculture practices is crucial to mitigate climate-related losses in yield and food quality, to reduce greenhouse gas emissions and to increase sustainable food production [3]. The mechanistic modelling approach is useful to improve our understanding of how crops develop, grow and produce a marketable yield, and is valuable for describing and predicting how endogenous (e.g., genetic makeup) and exogenous factors (e.g., weather and soil type) can impact crop performance.

AquaCrop, developed by the Food and Agriculture Organization (FAO) of the United Nations, is one of many published mechanistic models that can describe plant development from emergence to maturity, and how the efficiency and magnitude of physiological processes are moderated by multiple factors [4,5]. It was designed to be a generic model, based on fundamental biophysical principles and intended to be parsimonious and adaptable to specific crop species (http://www.fao.org/aquacrop). AquaCrop is now one of the most widely used models to simulate biomass and yield of different crops under a variety of climates and geographic locations.

AquaCrop v6.0 is currently available via an API provided by FAO (http://www.fao.org/aquacrop/software/en/) or by a Matlab API for the Open Source version AquaCrop-OS, developed by Foster et al. [6]. The FAO API runs as a standalone GUI, which allows the user to input climate and crop management data and returns predictions on yield performance based on the user’s inputs. Since the GUI is an executable file, users cannot modify, add or extend the model to account for more localised constraints, emerging problems and newly available data. The Matlab API was developed to fill this gap; however, since it is a script, the model can only be run using Matlab, which is a licence-based platform. To solve the problem of accessibility and extensibility, we have developed the AquaCropR API, which runs on the open source R platform [7]. AquaCropR is translated from AquaCrop-OS, which is based on an early version of AquaCrop GUI that did not simulate crop performance well in very dry environments. The latest version of AquaCrop GUI v6.0 addresses this and other issues.

To facilitate the use and extension of AquaCrop, we have deposited AquaCropR on the GitHub so that others can contribute to the development and enhancement of the model.

2. AquaCropR Overview

The AquaCropR package runs under the R environment [5] and was translated from the Matlab version of AquaCrop, AquaCrop-OS [6], which was based on equations incorporated in AquaCrop v5.0a developed by FAO [8,9]. Key features of AquaCropR that allow for rapid and easy extensibility, updating and integration with other systems are described in the following sections.

2.1. Open Source

AquaCropR is available on the GitHub repository (anyelacamargo/aquacropr) to allow for user collaboration and fast extension of the model. This is a key feature of AquaCropR, because it will facilitate the accessibly and use of the AquaCrop model to the wider community.

2.2. eXtensible Markup Language (XML)

AquaCropR uses the XML standard for data integration, data exchange and queries. Use of XML standard allows users to transfer data from other repositories to AquaCropR and to export and exchange data using the same format.

2.3. AquaCropR—New Functions

In addition to the core functionality of the current AquaCrop v6.0 and AquaCrop-OS that has been incorporated into AquaCropR, several functions have been added to the AquaCropR software to aid the user, described in the following subsections.

2.3.1. Evapotranspiration (ETo)

Weather variables, such as maximum and minimum temperature, precipitation and evapotranspiration, are needed in daily time steps in order to run the AquaCrop model. When potential evapotranspiration from a reference surface (ETo) data are not available, AquaCropR allows the user to calculate it using the Penman-Monteith equation, a generally accepted standard according to the FAO guidelines for computing crop evapotranspiration [10]. AquaCropR’s in-built function allows users to calculate ETo using available weather data such as mean temperature and solar radiation. The Penman-Monteith equation and guidelines can be found in Allen et al [9].

2.3.2. Solar Radiation

Global solar radiation is a critical input for AquaCrop, as the energy balance of the crop drives the water productivity function and therefore biomass accumulation. In the absence of evapotranspiration, AquaCropR calculates net solar radiation from incident solar radiation derived from a weather station radiometer or radiation modelled from satellite-derived data (e.g., https://power.larc.nasa.gov/ [11]). It is possible to roughly estimate global radiation from sunshine hours and site latitude, for example using adaptations of the Ångström-Prescott model (e.g., Muzathik et al. [12]), but the easy online access to satellite-derived radiation data makes this unnecessary.

2.3.3. Sowing Date Estimation

AquaCropR has a built-in function to estimate optimal sowing dates based on rainfall accumulation over a given number of successive days. Other issues that can be considered in future AquaCropR releases for an ‘optimal’ sowing date are: soil trafficability, suitable conditions for seed bed cultivations, weather forecast, probability of soil capping after a heavy rainfall event post-drilling or soil moisture content at depth of planting.

2.3.4. Optimisation Module

AquaCropR has a built-in function for parameter optimisation; currently, parameter optimisation is possible for the Canopy Growth Coefficient (CGC) and Maximum Canopy Cover as a fraction of ground area (CCx).

The procedure to use the optimiser module is included in the ‘AboutUsingParameterOptimisation’ of the AquaCropR manual available from the AquaCropR repo.

3. Installation and Running AquaCropR

3.1. Installation

To install AquaCropR, follow these steps:

• Clone the AquaCrop repository
  git clone https://bitbucket.org/anyelacamargo/aquacropr.git
• Run R
• Install the devtools package (if not available already)
  install.packages(“devtools”)
• Load the devtools package.
  library(devtools)
• Set working directory to AquaCrop’s location
  setwd(‘your location/AquaCropR’)
• Set working directory a level above AquaCropR
  setwd(‘..’)
• Install AquaCropR
  Install(‘AquaCropR’)
• Use the following command to test AquaCropR’s installation
  ?? AquaCropR

You should be able to see the package’s help.

3.2. Running a Simulation

The following is the procedure to run a validation test set to simulate wheat yields. More details on the dataset are given in Section 4.

• Clone the AquaCropExamples repository
  git clone https://github.com/anyelacamargo/AquaExamples.git
  Check that the ‘input_wheat_tunis_s1′ folder was pulled from the ‘AquaExamples’ repository.
• Load R
• Load AquaCropR package
  library(AquaCropR)
• Load crop simulation parameter files.
  FileLocation = ReadFileLocations(pathname/input_wheat_tunis_cropot/filesetup.xml’)
  Where pathname correspond to the location folder where the ‘AquaExamples’ repo was installed (e.g., ‘repo/AquaExamples’)
• Initialise the model using crop and simulation parameters
  InitialiseStruct <- Initialise(FileLocation)
• Run the simulation and view the results dataset. Check the AquaCropR tutorial for definition of the variables returned in the Simulation Results dataset.
  Outputs <- PerformSimulation(InitialiseStruct)

The output folder contains all the results from the simulation. For more details about each variable, refer to the AquaCropR tutorial (AquaCropR.pdf), also located in the ‘AquaExamples’ repository. All the wheat simulations can be run by running ‘simulate_dummy.R’ in the same folder, described in more detail in Section 4.

4. Testing and Evaluation

To ensure that AquaCropR can be used as alternative to the AquaCrop GUI counterpart (version 6.0), a number of test simulations were carried out using the same comparison exercises developed by Foster et al. [6] for a maize experiment, and by FAO for multiple wheat experiments used as training exercise for an AquaCrop training workshop [8].

4.1. Maize Exercise

The aim of this test was to validate the performance of AquaCropR against AquaCrop. The test scenario considers an irrigated maize crop in Nebraska, USA. Irrigation is triggered when soil moisture is depleted to 70% of the total available soil water (i.e., TAW, the water held between permanent wilting point and field capacity over the available rooting depth). The test is run over a 30-year period (1986–2015) using actual weather records (data and plots are available on https://github.com/anyelacamargo/AquaExamples. Inputs files (weather, crop management, soil properties, irrigation, crop calendar, water content) for the AquaCrop-OS tests are publicly available in tab separated format. We used the same files but transformed into the xml standard. Output files are also publicly available; we used them to compare the results from both software applications.

4.2. Wheat Exercises

The aim of the wheat exercises was to compare AquaCropR output with the AquaCrop GUI (v. 6.0) using four different simulation scenarios for a wheat crop in Tunis, Tunisia (36°48′ N, 10°10′ E) for 20 years (1983 to 2002). These are the same exercises that are used by FAO for AquaCrop training purposes [5]. In each scenario, the .CRO (crop file), .SOL (soil file) and .SW0 (Initial water content) files available with the AquaCrop GUI were used in the tests. The remaining input settings were set as default as no further information was available. Observed empirical data were not available either. Therefore, to compare outputs from both applications, the AquaCrop GUI was run for each test scenario and year and results were collected individually and stored in a table. Weather data were sourced from the POWER Data Access Viewer [11]), which allows automatic download and use of public data and were used in as inputs for AquaCropR and the AquaCrop GUI APIs.

To collect data from the AquaCropR, the software was run in batch through a script and result files were collected at the end on the run. Result tables from AquaCrop GUI and AquaCropR were compared using standard statistical metrics.

4.2.1. Wheat Scenario 1: Default Soil

This test assessed AquaCropR performance in predicting wheat yields under dryland growing conditions, in comparison with the AquaCrop GUI using similar inputs. The initial water content in the soil profile of a uniform sandy loam was characterised by a wet top soil and dry subsoil. Crop-specific parameters were those given in the input file ‘WheatGDD.CRO’ (for Growing Degree Day (GDD)) provided with the AquaCrop GUI. In each year, sowing day was on 15th October and soil characteristics are in

Table 1. Soil characteristics for wheat scenarios 1 and 2. Soil characteristics for scenarios 3 and 4 were the same as scenario 1. Abbreviations: TAW, total available soil water; PWP, permanent wilting point; FC, field capacity; SAT, water content at saturation; Ksat, saturated hydraulic conductivity.

Scenario	Soil Texture Class	Soil Layer (m)	Soil Water Content (Volume %)				Ksat (mm d−1)
			TAW	PWP	FC	SAT
1	sandy loam	0–0.3	12	10	22	41	500
	sandy loam	0.3–0.9	12	10	22	41	500
2	clay loam	0–0.3	16	24	40	50	155
	silt loam	0.3–1.7	22	11	33	46	500

.

4.2.2. Wheat Scenario 2: A Typical Local Soil

This test assessed yield on a local soil typical of the growing region in Tunis, with a wet topsoil and dry subsoil. Soil characteristics are indicated in Table 1. Other crop parameters, sowing date and weather inputs were the same as Scenario 1.

4.2.3. Wheat Scenario 3: A Late-Maturing Wheat Variety

This test assessed yield for a wheat variety with a longer crop cycle than the default variety used in the AquaCrop tutorials. The aim of this exercise was to test if the two AquaCrop versions showed similar yield responses to seasonal variations in growing conditions. In each year, sowing day was 15 October and soil characteristics were the same as in Scenario 1 (Table 1). The input file ‘WheatGDD.CRO’ was modified for a local variety for these parameters, as indicated in [8]:

• Time from sowing to emergence (GDD) = 289
• Time from sowing to flowering (GDD) = 1073
• Time from sowing to senescence (GDD) = 2835
• Time from sowing to physiological maturity (GDD) = 3390
• Duration of flowering (GDD) = 264
• Planting density (plants m⁻²) = 350
• Maximum canopy cover (%) = 90

4.2.4. Wheat Scenario 4: Applying an Irrigation Schedule

In this scenario, wheat yield response to deficit irrigation was run to test the response of the models to different water inputs. According to local practices, farmers sow irrigated wheat at the start of the winter period into dry soil, which had reached permanent wilting point by 15 August. The soil characteristics were the same as in Scenario 1 (Table 1). Since supplies of irrigation water are limited, only 30 mm water is applied after sowing (1 December as used in the AquaCrop GUI original exercise) to ensure germination. Subsequently, two further applications of 40 mm each were applied on 1 February and 1 March.

4.3. Results

4.3.1. Maize

Simulated yield of an irrigated maize crop was on average 13.44 t/ha for both AquaCrop-OS and AquaCropR, with no significant deviations between the model outputs observed across the tested years (Root Mean Square Error (RSME) = 0; Pearson Correlation (R²) = 1) (Figure 1). This shows that the functions coded in AquaCropR faithfully reproduced those in AquaCrop-OS.

4.3.2. Wheat

Results of the default wheat simulations are shown in Figure 2. Yield potential of rainfed wheat was on average 8.6 ton ha⁻¹ (Figure 2A) using the AquaCrop GUI and AquaCropR (RSME = 0.32; R² = 0.91). Only small deviations in outputs between the models were observed across the tested years. Further speculation on the cause of the small deviations between the two versions are discussed below.

The soil type tests showed a decline in yield when a local soil was used (Figure 2B) in comparison to a uniform sandy loam (Figure 2A), with all other parameters equal. This shows that AquaCropR is sensitive to changes in soil moisture and that the yield responses are reasonable and accurate, compared with the GUI version AquaCrop 6.0, (RSME = 0.47; R² = 0.97). In years where there was some discrepancy between the models, AquaCropR tended to output slightly larger yields than the GUI, which suggests the model was supplying slightly more water to the crop than in the GUI version.

Using the crop configuration inputs for a local variety showed greater year to year fluctuations in yields (Figure 2C) compared with tests using the default wheat variety (Figure 2A). The main difference can be attributed to a longer grain filling period (flowering to maturity; 2317 GDD for the local variety versus 1150 GDD in the default setting in WheatGDD.CRO), and shorter duration from sowing to flowering (1073 for the local variety versus 1250 GDD in the default). In years with favourable growing conditions, the local variety yields surpassed those of the default, and in poor years the local variety yielded less than the default (Figure 2C). It is surmised that the timings of phenological growth stages of the local variety made it more sensitive to prevailing environmental conditions, such as increased heat later in the crop cycle. There were also greater deviations between the outputs of the AquaCrop GUI and AquaCropR compared with the default crop settings (Figure 2C; RSME = 0.15; R² = 0.86). Compared with the AquaCrop GUI, AquaCropR tended to output smaller yields in low yielding years, and output larger yields in good years using the local variety characteristics. The AquaCrop training manual tutorial files for the Tunis scenario 3 state that the observed grain yield in trials in the 1990/1991 season was 9.15 t/ha [8]. This agrees well with the AquaCropR simulated yield of 9.39 t/ha (Figure 2C). In the absence of further observed yield data for the years in these test environments, it is not clear which model output is closer to reality. However, there are numerous reports in the literature that show validations of the AquaCrop GUI model using observed data. If yields in the poor years were limited by moisture to a greater extent than in the good years, this suggests that the AquaCrop GUI may overestimate soil water availability in severe stress, which has been suggested by others [13], or the inverse case in AquaCropR. Further work is needed to drill down into the input conditions for those years to understand how the models deviate in those circumstances, and to determine which model settings need to be adjusted.

The test comparing simulated yields across years using a fixed, deficit irrigation schedule shows that the model performed under conditions of limited soil moisture supply (weather data and plots are available on https://github.com/anyelacamargo/AquaExamples/) similarly to the AquaCrop GUI, and the agreement between models was good (RSME = 0.95; R² = 0.77). AquaCropR output slightly larger yields in several years, but showed larger deviations with small yields in six of the 22 years (Figure 2D). Since it is important for the reader and AquaCropR user to understand what factors were important for these differences between model outputs,

Table 2. Differences between AquaCrop GUI and AquaCropR.

	Simulate	AquaCrop GUI V6.0	AquaCropR
Simulation of the early development of the canopy cover under water stress	The expansion rate of seedling canopy cover (CC) is no longer limited by water stress when the crop germinates.	Updated accordingly	Remains as original. Under development for new version.
Simulation of root deepening in a dry subsoil	If the soil water depletion at the front of root zone expansion exceeds a specific threshold, the root deepening will slow down in Version 6.0, and can even become inhibited if the soil water content at the front is at permanent wilting point	Updated accordingly	Remains as original. Under development for new version.
Simulation of soil water stress	The depletion in the total root zone with the depletion in the top soil is compared at each time step. This determines which part of the soil profile is the wettest and will determine the degree of water stress.	Updated accordingly	Remains as original. Under development for new version.
Simulation of canopy cover decline	The equation to model the decline in green crop canopy has been modified to make the simulation of senescence duration less divergent for different CCx and especially for small crops.	$CC=C{C}_x\left[1-0.05\left(e^{\frac{3.33 CDC}{C{C}_t+2.29}t}-1\right)\right]$	$CC=C{C}_x\left[1-0.05\left(e^{\frac{CDC}{C{C}_t}t}-1\right)\right]$
Simulation of cold stress	Assigning the crop transpiration coefficient (KcTr) as the target parameter for the cold stress coefficient.	$T{r}_i=K{s}_i\left(K{s}_{Tr,i}K{C}_{Tr,i}\right)ET{o}_i$	$T{r}_i=K{s}_i K{c}_{Tr.i}ET{o}_i$
Calibration and simulation of salinity stress	The calibration of the crop for soil salinity stress is developed and is unlinked with the calibration for soil fertility. Simulation of salinity stress considers ECe (Electrical Conductivity of the saturated soil-paste extract), which is the indicator for soil salinity and ECsw (Electrical Conductivity of the soil water). When the soil dries out, the increase of ECsw results in the increase of osmotic effects causing a stronger closure of the stomata, and the simulation of a stronger reduction of crop transpiration.	New module	Not available. Under development for new version.
Extra soil characteristics	Possibility to specify the penetrability of a soil horizon which describes the effect on the expansion rate of the root zone. Possibility to specify gravel in a soil horizon	New options	Not available. Under development for new version.

is designed to highlight key differences between AquaCropR/OS and AquaCrop GUI, v 6.0 [14]. Looking at the table, the user can immediately realise that under no stress conditions both models are likely to simulate similar crop growth and performance. However, under high stress conditions, AquaCrop GUI is more likely to simulate actual crop development and performance—this was the reason for v 6.0 being developed. As the source code for v 6.0 is not publically available, we are currently using the v6.0 manual [14] to upgrade AquaCropR.

In addition to the grain yield plots shown here, biomass, canopy cover, soil infiltration, irrigation (when tested) and Et₀ plots are also available in the ‘AquaExamples’ repository on GitHub. The plotted data are located in the output folder of each test (e.g., input_wheat_tunis_cropot\output).

5. Summary

This study developed the open source R version of AquaCrop, which we called AquaCropR. To demonstrate AquaCropR usability and versatility we used crop management data from published sources and weather data from publicly available repositories to fit simulate five crop scenarios. The results of the test demonstrated that AquaCropR can be easily used to simulate or predict yields under different crop and experimental conditions.

AquaCropR is translated from AquaCrop-OS which in turn is based on an early version of AquaCrop GUI, which did not simulate well crop performance under very dry environments. Therefore, when comparing AquaCropR against AquaCrop-OS, we expected the same results. Yield comparisons between AquaCrop-OS and AquaCropR outputs for maize crops growing under no stress conditions showed highly similar values (RSME = 0; R² = 1).

On the contrary, when comparing AquaCropR against AquaCrop GUI, we expected some degree of disagreement because the models are slightly different (the latter has been updated to deal with simulations of crop growing under dry environments). This behaviour was also demonstrated in Figure 2C,D which corresponded to yield comparisons of wheat crops growing under some sort of stress (e.g., S3 Late maturity wheat variety and S4 Water deficit).

In the current work, we are conducting collaborative maize field trials under a range of N inputs over several seasons. These have been running for two years to collect data from maize trials. The field data will be used to tune, test and expand the AquaCropR API.

The open access to source code in R, and the flexibility to adapt this to new conditions and crops should further widen the use of AquaCrop. The danger of open access is that erroneous changes can be made to formulae that become incorporated into newer versions of the software as it evolves. It is recommended that new versions of AquaCropR use these standardised test datasets to validate performance before using the model. In future work, addition of a soil nutrient module will permit better simulations under nutrient-limited conditions, as the current versions assume non-limiting conditions, or only allow election of a decline in crop coefficient due to a nitrogen ‘deficiency’. Another future facility we are working on is the addition of mathematical optimisation to fit model parameters instead of using constant parameters. The current AquaCropR version is able to optimise Canopy growth coefficient (CGC) and Maximum Canopy Cover in Fraction Soil (CCx) to AquaCropR, but a more robust optimization module is currently under expansion.