knitr::opts_chunk$set( collapse = TRUE, comment = "#>" )
The 'WetlandEconsRPackage' provides a suite of functions (5 functions) that could be used in concert to analyse the economics of wetland drainage in an agricultural landscape. In this vignette, we will explain how to use all the functions to estimate the net present value of wetland drainage. First, we will explain how to install the package and make it ready for analysis. Second, we will show how to use the first function (wetlandscape) to simulate wetland characteristics in a sub-basin. Third we will show how to use the family of net-present value functions to estimate the net-present value of wetland drainage at the quarter section level and at the higher resolutions.
The package can be installed in 2 simple steps. First, install the devtools package. Second, use the devtools to install the package which is stored at github. Finally, we can use the library function to activate the package for analysis.
# step1- install.packages("devtools") #Step2- devtools::install_github("tidyverse/WetlandEconsRPackage") library(WetlandEconsRPackage)
The functions that are in the package are wetlandscape, npv_3yrotation, npv_4yrotation, npv_5yrotation, prop_rotation. I urge you to read more about the functions using the help:"?WetlandEconsRPackage::function". For instance, read the wetlandscape function with ?WetlandEconsRPackage::wetlandscape. The help document will show you the parameters in the functions and their meanings.
The wetland function in the WetlandEconsRPackage is used to simulate wetlands in a sub-basin (a typical Canadian Prairie Pothole region landscape). The simulation is at the quarter section level, which is a farmland 64.75 ha(160 acres) in size. It produces a dataframe with the following columns: a) wetland size (ha), b) wetland groupings (tier), c) probability of harvesting (pr), d) number of wetlands (nwl), e)wetland location(wl_location).
We follow Pomeroy et al. (2002) to group the wetlands into 5 tiers. The wetland tiers and how they are connected is shown in figure 1 below:
The wetland sizes (ha) for the tiers are created with uniform distribution the wetland minimum size and maximum size parameters (which are controlled by the user). For instance, if I want the wetland sizes in tier 1 to be between 1 and 3 hectares, I will set w1_min = 1 (wetland minimum size) and w1_max = 3 (wetland maximum size).
Similarly, the probability of harvesting parameter (pr), which shows the probability that a farmer can harvest from a converted wetland, for wetlands in the tiers is controlled bv the user.The parameter is between 0 (minimum pr) and 1 (maximum pr); a parameter of 1 means that the farmer will be able to obtain 100% of crop yield from the drained wetland and 0 means the farmer will get 0 yield from the drained wetland. For instance, if we believe wetlands in tier 1 have probability for wetlands between 50% and 70% we will set the minimum and maximum pr as p1_min = 0.5 and p1_max = 0.7, respectively.
The number of wetlands per quarter section for the entire subbasin is created using a random discrete distribution, with the parameters minimum number of wetlands(nwl_min) and maximum number of wetlands (nwl_max). Again, the wetland drainage cost for the entire subbasin is assumed to be randomly distributed with the parameters minimum drainage cost (dc_min) and maximum drainage cost(dc_max) which are controlled by the user. Similarly the locations of wetlands in a the quarter sessions (center or corner) are randomly generated for each quarter session.
Use Case: lets say we want to simulate wetland distribution for 100 quarter sections. We believe wetland sizes become progressively bigger as we move from tier 1 wetlands; also, and the probability of harvesting increases as we move from tier 1 wetlands. Also, we believe the number of wetlands per quarter section range from 1 to 4 for the entire sub-basin, and the minimum and maximum drainage cost range from 200 to 600 $/ha. The dataframe corresponding to the above assumptions is created with the wetlandscape function below:
library(DT) Subbasin_wetlands <- wetlandscape( nq = 100, w1_min = 0.1, w1_max = 0.2, w3_min = 0.21, w3_max = 0.6, w6_min = 0.61, w6_max = 0.8, w12_min= 0.81, w12_max= 0.9, w24_min= 0.9, w24_max= 50, p1_min= 0.5, p3_min= 0.7, p6_min= 0.8, p12_min= 0.85, p24_min = 0.86, dc_min = 200, dc_max = 600, nwl_min = 1, nwl_max= 4 ) datatable(Subbasin_wetlands)
The NPV_5yrotation is a function that is used to estimate the net present value of wetland drainage for a 5-year annual crop rotation plan. This means 5 different crops or 4 different crops and a fallow are alternated yearly until the end of the planning horizon. The arguments in the function are production profits for the crops, wetland size, drainage cost, discount rate, and planning horizon; the user is able to specific the levels of the parameters based on his or her study assumptions.
Use Case: In this use case we will use a Canola-Spring Wheat-Flax-Barley-Fallow annual rotation to estimate the net present value for wetland drainage for a 50-year planning horizon. We will assume a 8.7% discount rate. Also, we will use the 'wetlandscape' function to simulate wetland characteristics as in 3.1. Again, for illustrative purposes, we will assume profit levels for the crops and add randomness to it so that it varies across the quarter-sections.In real studies we could obtain crop production data from provincial or federal agencies. The other crop rotation functions, namely NPV_4yrotation and NPV_3yrotation can also be implemented like the NPV_5yrotation.
npv_parameters <- data.frame(canola_profit = 300 * runif(nrow(Subbasin_wetlands), 0.6, 1), spwheat_profit = 250 * runif(nrow(Subbasin_wetlands), 0.6, 1), flax_profit = 150 * runif(nrow(Subbasin_wetlands), 0.6, 1), barley_profit = 100 * runif(nrow(Subbasin_wetlands), 0.6, 1), fallow = 0* runif(nrow(Subbasin_wetlands), 0.6, 1)) %>% cbind(Subbasin_wetlands) npv <- npv_parameters %>% dplyr::mutate( npv = unlist(purrr::pmap(list(np1= canola_profit, np2= spwheat_profit, np3= flax_profit, np4= barley_profit, np5= fallow, wl= wl, dc=dc, r=0.87, t=50),npv_5yrotation)) ) %>% dplyr::select(wl, npv) datatable(npv)
The Prop_rotation function works like the yearly rotation functions, but in year year of the rotation multiple crops are cultivated proportional to to drained wetland area. The arguments of the function are crop profits, proportion of wetland areas devoted to a specified crops, wetland size, drainage cost, discount rate and planning horizon.
Use Case: We will illustrate this function by maintaining the assumptions in 3.2. The only difference lies in the how the crop rotation is performed yearly. Instead of one crop per year as in 3.2, this function will divide the wetland area among the specified crops.We will use the same data in 3.2.
npv_prop <- npv_parameters %>% dplyr::mutate( npv_proportional = unlist(purrr::pmap(list( #crops that are proportionally cultivated on drained wetland for year 1 y1np1= canola_profit, y1np2= spwheat_profit, y1np3= flax_profit, y1np4= barley_profit, y1np5= fallow, y1pwlnp1 = 0.5, y1pwlnp2 = 0.1, y1pwlnp3 = 0.2, y1pwlnp4 = 0.2, y1pwlnp5 = 0, #crops that are proportionally cultivated on drained wetland for year 2 y2np1= canola_profit, y2np2= spwheat_profit, y2np3= flax_profit, y2np4= barley_profit, y2np5= fallow, y2pwlnp1 = 0.3, y2pwlnp2 = 0.2, y2pwlnp3 = 0.2, y2pwlnp4 = 0.3, y2pwlnp5 = 0, #crops that are proportionally cultivated on drained wetland for year 3 y3np1= canola_profit, y3np2= spwheat_profit, y3np3= flax_profit, y3np4= barley_profit, y3np5= fallow, y3pwlnp1 = 0.2, y3pwlnp2 = 0.5, y3pwlnp3 = 0.2, y3pwlnp4 = 0.1, y3pwlnp5 = 0, #crops that are proportionally cultivated on drained wetland for year 4 y4np1= canola_profit, y4np2= spwheat_profit, y4np3= flax_profit, y4np4= barley_profit, y4np5= fallow, y4pwlnp1 = 0.4, y4pwlnp2 = 0.4, y4pwlnp3 = 0.1, y4pwlnp4 = 0.1, y4pwlnp5 = 0, #crops that are proportionally cultivated on drained wetland for year 5 y5np1= canola_profit, y5np2= spwheat_profit, y5np3= flax_profit, y5np4= barley_profit, y5np5= fallow, y5pwlnp1 = 0.5, y5pwlnp2 = 0.1, y5pwlnp3 = 0.2, y5pwlnp4 = 0.2, y5pwlnp5 = 0, wl = wl, dc = dc, r = 0.87, t = 50), prop_rotation)) ) %>% dplyr::select(wl, npv_proportional) datatable(npv_prop)
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.