1. Introduction
This vignette is a short tutorial to the usage and functionalists of AquaBEHER R package. It is directed to first-time package users who are familiar with the basic concepts of R. The vignette presents the use of several key functions of the package, some useful hints and guidelines. AquaBEHER computes and integrates daily reference evapotranspiration (Eto) and a soil water balance model to estimate parameters of crop and soil water balances for agricultural crops. Using the computed daily soil water balance parameters, the package can estimates the rainy season calendar (Onset, Cessation and Duration) based on agroclimatic approach for a predefined window.
2. Installation and Loading
Installing the latest development version from the online repository using the devtools package. Note that to utilize the full functionality of devtools on Windows, Rtools must be installed.
# install.packages("devtools")
# devtools::install_github("RobelTakele/AquaBEHER", dependencies = TRUE, type = "source",
# build_manual = TRUE, build_vignettes = TRUE)
library(AquaBEHER)
library(ggplot2)
3. Required Climate Data
The methods for calculating evapotranspiration from meteorological data require various physical parameters. Some of the data are measured directly in weather stations. Other parameters are related to commonly measured data and can be derived with the help of a direct or empirical relationship.
The meteorological factors determining evapotranspiration are weather parameters which provide energy for vaporization and remove water vapor from the evaporating surface. The principal weather parameters to consider are presented below.
- Maximum temperature
- Minimum temperature
- Solar radiation
- Dew point temperature or relative humidity
- Wind speed
In addition, georeferenced information on the location of the climate record is required:
- Latitude
- Longitude
- Elevation
data(AgroClimateData)
str(AgroClimateData)
#> 'data.frame': 14975 obs. of 14 variables:
#> $ GridID: chr "MOZ0007149" "MOZ0007149" "MOZ0007149" "MOZ0007149" ...
#> $ Lat : num -15.1 -15.1 -15.1 -15.1 -15.1 ...
#> $ Lon : num 39.3 39.3 39.3 39.3 39.3 ...
#> $ Elev : num 392 392 392 392 392 ...
#> $ WHC : num 97.8 97.8 97.8 97.8 97.8 ...
#> $ Year : num 1982 1982 1982 1982 1982 ...
#> $ Month : num 1 1 1 1 1 1 1 1 1 1 ...
#> $ Day : num 1 2 3 4 5 6 7 8 9 10 ...
#> $ Rain : num 0 0 0 1.91 0 ...
#> $ Tmax : num 32.2 33.1 33.5 32.8 32.7 ...
#> $ Tmin : num 23.1 23.1 23.1 23.6 22.8 ...
#> $ Rs : num 23.9 26.4 25 24.2 23.4 ...
#> $ Tdew : num 20.2 20.5 20.5 20.8 21.4 ...
#> $ Uz : num 4.72 4.28 3.62 2.54 1.48 ...
head(AgroClimateData)
#> GridID Lat Lon Elev WHC Year Month Day Rain
#> 1 MOZ0007149 -15.09238 39.2519 392.1337 97.84914 1982 1 1 0.000000
#> 2 MOZ0007149 -15.09238 39.2519 392.1337 97.84914 1982 1 2 0.000000
#> 3 MOZ0007149 -15.09238 39.2519 392.1337 97.84914 1982 1 3 0.000000
#> 4 MOZ0007149 -15.09238 39.2519 392.1337 97.84914 1982 1 4 1.907393
#> 5 MOZ0007149 -15.09238 39.2519 392.1337 97.84914 1982 1 5 0.000000
#> 6 MOZ0007149 -15.09238 39.2519 392.1337 97.84914 1982 1 6 0.000000
#> Tmax Tmin Rs Tdew Uz
#> 1 32.24396 23.11500 23.86698 20.21160 4.723783
#> 2 33.07202 23.12585 26.38375 20.48284 4.279407
#> 3 33.49679 23.12602 25.00704 20.45689 3.622179
#> 4 32.76818 23.60351 24.16475 20.83896 2.535047
#> 5 32.65872 22.79294 23.44483 21.36882 1.477617
#> 6 31.80630 22.43975 21.99277 21.29297 1.953415
4. Potential Evapotranspiration
For many agricultural applications, it is relevant to get an estimate of the potential evapotranspiration (PET). Different methods are developed for estimating Eto. Most of them use empirical equations to determine the value of PET from weather variables. The AquaBEHER package provides options for estimating reference evapotranspiration (Eto) using the FAO Penman-Monteith, Priestley Taylor and Hargreaves-Samani formulations.
Usage
calcEto(AgroClimateData, method = “PM”, crop = “short”)
- The function give a list of
- daily estimations of Eto in (mm/day)
- daily estimations of extraterrestrial radiation (MJ/m2/day)
- daily estimations of slope of vapor pressure curve (kPa/°C)
Example
The calcEto compute with inputs of data frame containing daily values of meteorological parameters:
<- calcEto(AgroClimateData, method = "PM", crop = "short")
PET
str(PET)
#> List of 7
#> $ ET.Daily : num(0)
#> $ Ra.Daily : num [1:14975] 40.9 40.9 40.9 40.9 40.9 ...
#> $ Slope.Daily : num [1:14975] 0.217 0.221 0.224 0.222 0.217 ...
#> $ Ea.Daily : num [1:14975] 2.37 2.41 2.41 2.46 2.54 ...
#> $ Es.Daily : num [1:14975] 3.82 3.94 4 3.94 3.85 ...
#> $ ET.formulation: chr "Penman-Monteith FAO56"
#> $ ET.type : chr "Reference Crop ET"
Graphical comparison of the evapotranspiration (mm/day) calculated using the FAO Penman–Monteith formulation and the Hargreaves-Samani formulation:
# Compute Eto using hargreves-samani formulation using the example data from 'AgroClimateData':
data(AgroClimateData)
<- calcEto(AgroClimateData, method = "HS")
Eto.HS
# Now compute Eto using Penman-Monteith formulation for hypothetical grass (short crop):
<- calcEto(AgroClimateData, method = "PM", Zh = 10)
Eto.PM
plot(Eto.PM$ET.Daily[1:1000], ty="l", xlab="Days since 1996", ylab="Eto (mm/day)", col="black", lwd = 1, lty = 2)
lines(Eto.HS$ET.Daily[1:1000], col="blue", lwd = 2, lty = 1)
legend("bottom",c("Eto: Penman–Monteith ","Eto: Hargreaves-Samani"),
horiz=TRUE, bty='n', cex=1,lty=c(2,1),lwd=c(2,2), inset=c(1,1),
xpd=TRUE, col=c("black","blue"))
Notice: It appears that the FAO Penman–Monteith formulation presents enhanced day-to-day variations of evapotranspiration than the the Hargreaves-Samani formulation.
5. Soil Water Balance
This function perform daily computations of soil water balance parameters for the root zone. Soil water changes daily in response to rainfall, evapotranspiration, runoff and deep drainage.
Assumptions
- Atmospheric conditions affect the rate at which crops use water.
- The soil has uniform cross-section of homogeneous volume with a measured depth and a unit area.
- A well-established, dense grass crop is growing, which completely covers the soil surface.
Example
The calcWatBal compute with inputs of data frame containing daily values of Rain, Eto and soil water holding capacity.
<- calcEto(AgroClimateData, method = "PM", Zh = 10)
PET
# Add the estimated PET 'ET.Daily' to a new column in AgroClimateData:
$Eto <- PET$ET.Daily
AgroClimateData
# Estimate daily water balance for the soil having 100mm of WHC:
<- 100
soilWHC
<- calcWatBal(data = AgroClimateData, soilWHC)
watBal
str(watBal )
#> 'data.frame': 14975 obs. of 20 variables:
#> $ GridID: chr "MOZ0007149" "MOZ0007149" "MOZ0007149" "MOZ0007149" ...
#> $ Lat : num -15.1 -15.1 -15.1 -15.1 -15.1 ...
#> $ Lon : num 39.3 39.3 39.3 39.3 39.3 ...
#> $ Elev : num 392 392 392 392 392 ...
#> $ WHC : num 97.8 97.8 97.8 97.8 97.8 ...
#> $ Year : num 1982 1982 1982 1982 1982 ...
#> $ Month : num 1 1 1 1 1 1 1 1 1 1 ...
#> $ Day : num 1 2 3 4 5 6 7 8 9 10 ...
#> $ Rain : num 0 0 0 0 0 ...
#> $ Tmax : num 32.2 33.1 33.5 32.8 32.7 ...
#> $ Tmin : num 23.1 23.1 23.1 23.6 22.8 ...
#> $ Rs : num 23.9 26.4 25 24.2 23.4 ...
#> $ Tdew : num 20.2 20.5 20.5 20.8 21.4 ...
#> $ Uz : num 4.72 4.28 3.62 2.54 1.48 ...
#> $ Eto : num 6.45 6.8 6.49 5.84 5.23 ...
#> $ R : num 0 0 0 0 0 0 0 0 0.684 1 ...
#> $ AVAIL : num 0 0 0 0 0 ...
#> $ TRAN : num 0 0 0 0 0 ...
#> $ DRAIN : num 0 0 0 0 0 0 0 0 0 0 ...
#> $ RUNOFF: num 0 0 0 0 0 ...
# Plotting the water balance output for the climatological year from 2019 to 2020 using ggplot2:
.19T20 <- watBal[watBal$Year %in% c(2019, 2020),]
watBal<- as.Date.character(paste0(watBal.19T20$Year, "-", watBal.19T20$Month, "-", watBal.19T20$Day))
date.vec
ggplot(data = watBal.19T20) +
geom_line(aes(y = AVAIL, x = date.vec, fill = "AVAIL"), size = 0.8, color = "red") +
geom_col(aes(y = Rain, x = date.vec, fill = "Rain"), size = 1) +
scale_x_date(date_breaks = "1 month", date_labels = "%b-%Y") +
scale_fill_manual(name = " ", values = c('AVAIL' = "red", 'Rain' = "blue")) +
scale_y_continuous(expand = c(0, 2)) +
labs(y="Moisture (mm)", x=NULL) +
theme_linedraw() +
theme(axis.title = element_text(size = 14, colour = "black", family = "Times New Roman"),
axis.text = element_text(size = 10, colour = "black", family = "Times New Roman"),
axis.text.x = element_text(size = 10, colour = "black", family = "Times New Roman", angle = 45, vjust = 0.5))
#> Warning: Using `size` aesthetic for lines was deprecated in ggplot2 3.4.0.
#> ℹ Please use `linewidth` instead.
#> This warning is displayed once every 8 hours.
#> Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
#> generated.
#> Warning in geom_line(aes(y = AVAIL, x = date.vec, fill = "AVAIL"), size = 0.8,
#> : Ignoring unknown aesthetics: fill
6. Rainy Season Calandar
The onset and cessation dates of the rainy season were determined for
each climatological year Figure 1
. The term climatological
year represents the period between two driest periods, which is
traditionally defined based on a calender year starting from the driest
month and has a fixed length of 12 months.
Various methods have been developed to estimate the rainy season calendar, i.e. the onset, cessation and duration of the rainy season. Common method used for crop production applications is the agroclimatic approach. As per agroclimatic approach, a normal rainy season (growing season) is defined as one when there is an excess of precipitation over potential evapotranspiration (PET). Such a period meets the evapotransiration demands of crops and recharge the moisture of the soil profile FAO 1977; 1978; 1986). Thus, the rainy season calendar defined accordingly:
Onset
The onset of the rainy season will start on the first day after
onsetWind.start
, when the actual-to-potential
evapotranspiration ratio is greater than 0.5 for 7 consecutive days,
followed by a 20-day period in which plant available water remains above
wilting over the root zone of the soil layer.
Cesation
The rainy season will end, cessation, on the first day after
onsetWind.end
, when the actual-to-potential
evapotranspiration ratio is less than 0.5 for 7 consecutive days,
followed by 12 consecutive non-growing days in which plant available
water remains below wilting over the root zone of the soil layer.
Duration
The duration of the rainy season is the total number of days from onset to cessation of the season.
Example
Using the example climate data provided by the AquaBEHER package, compute the rainy season calandar:
# seasonal calndar is estimated for the onset window ranges from 01-September to 31-January having a soil with 100mm of WHC
<- 100
soilWHC <- "1996-09-01" # earliest possible start date of the onset window
onsetWind.start <- "1997-01-31" # the latest possible date for end of the onset window
onsetWind.end <- "1997-06-30" # the latest possible date for end of the cessation window
cessaWind.end
<- calcSeasCal(watBal, onsetWind.start, onsetWind.end, cessaWind.end, soilWHC = 100)
seasCal.lst
str(seasCal.lst)
#> List of 3
#> $ :'data.frame': 40 obs. of 6 variables:
#> ..$ Year : num [1:40] 1982 1983 1984 1985 1986 ...
#> ..$ onset.Year : num [1:40] 1982 1983 1984 1985 1986 ...
#> ..$ onset.Month: num [1:40] 11 12 11 11 11 1 12 11 12 12 ...
#> ..$ onset.Day : num [1:40] 11 11 9 19 20 10 18 13 26 3 ...
#> ..$ onset.JD : num [1:40] 315 345 314 323 324 10 353 317 360 337 ...
#> ..$ onset.Value: int [1:40] 72 102 70 80 81 132 109 74 117 94 ...
#> $ :'data.frame': 40 obs. of 6 variables:
#> ..$ Year : num [1:40] 1982 1983 1984 1985 1986 ...
#> ..$ cessation.Year : num [1:40] 1983 1984 1985 1986 1987 ...
#> ..$ cessation.Month: num [1:40] 5 4 5 5 4 4 5 3 5 5 ...
#> ..$ cessation.Day : num [1:40] 12 18 9 10 22 21 6 15 6 21 ...
#> ..$ cessation.JD : num [1:40] 132 109 129 130 112 112 126 74 126 142 ...
#> ..$ cessation.Value: num [1:40] 255 232 252 253 235 235 249 197 249 265 ...
#> $ :'data.frame': 40 obs. of 4 variables:
#> ..$ Year : num [1:40] 1982 1983 1984 1985 1986 ...
#> ..$ onset.YYYYDOY : chr [1:40] "1982-315" "1983-345" "1984-314" "1985-323" ...
#> ..$ cessation.YYYYDOY: chr [1:40] "1983-132" "1984-109" "1985-129" "1986-130" ...
#> ..$ Duration : int [1:40] 183 130 182 173 154 103 140 123 132 171 ...
# plotting year to year variation of onset cessation and seasonal duration
<- data.frame(Year = seasCal.lst[[1]][,c("Year")],
seasCal.dF Onset = seasCal.lst[[1]][,c("onset.JD")],
Cessation = seasCal.lst[[2]][,c("cessation.JD")],
Duration = seasCal.lst[[3]][,c("Duration")])
ggplot(data = seasCal.dF) +
geom_line(aes(y = Onset, x = Year, color = "Onset"), size = 1) +
geom_line(aes(y = Cessation, x = Year, color = "Cessation"), size = 1) +
geom_area(aes(y = Duration, x = Year, color = "Duration"), size = 0.8, alpha = 0.4)+
scale_color_manual(name = "Calendar", values = c('Onset' = "blue", 'Cessation' = "red", 'Duration' = "grey")) +
labs(y="Day of a year (DOY)", x=NULL) +
theme_bw()
7. References
Allen, R.G.; Pereira, L.S.; Raes, D.; Smith, M. Crop Evapotranspiration: Guidelines for Computing Crop Water Requirements; FAO Irrigation and Drainage Paper no. 56; FAO: Rome, Italy, 1998; ISBN 92-5-104219-5.
Doorenbos, J. and Pruitt, W.O. 1975. Guidelines for predicting crop water requirements, Irrigation and Drainage Paper 24, Food and Agriculture Organization of the United Nations, Rome, 179 p.
Hargreaves, G.H. and Samani, Z.A. (1985) Reference Crop Evapotranspiration from Temperature. Applied Engineering in Agriculture, 1, 96-99.
\[\\[0.2in]\]
The Genetics Group at the Center of Plant Sciences is a geographically and culturally diverse research team working on data-drivem agicultural innovation combining crop genetics, climate, and participatory approaches. We are based at Scuola Superiore Sant’Anna, Pisa, Italy.
You can contact us sending an email to Matteo Dell’Acqua (mailto:m.dellacqua@santannapisa.it) or Mario Enrico Pè (mailto:m.pe@santannapisa.it). You can find out more about us visiting the group web page (http://www.capitalisegenetics.santannapisa.it/) and following us on Twitter @GenLab_SSA
We are committed to the free software and FAIR principles. This set of repositories collects our latest developments and provide reusable code.