# Making a Model

This tutorial will cover the topic of creating a Cropbox model.

## Growing Degree-Day

You might have heard the terms like growing degree days (GDD), thermal units, heat units, heat sums, temperature sums, and thermal-time that are used to relate the rate of plant or insect development to temperature. They are all synonymous. The concept of thermal-time or thermal-units derives from the long-standing observation and assumption that timing of development is primarily driven by temperature in plants and the relationship is largely linear. The linear relationship is generally held true over normal growing temperatures that are bracketed by the base temperature (Tb) and optimal temperature (Topt). Many existing crop models and tree growth models use thermal-unit approaches (e.g., GDD) for modeling phenology with some modifications to account for other factors like photoperiod, vernalization, dormancy, and stress. The growing degree days (GDD) is defined as the difference between the average daily air temperature (T) and the base temperature below which the developmental process stops. The bigger the difference in a day, the faster the development takes place up to a certain optimal temperature (Topt). The Cumulative GDD (cGDD) since the growth initiation (e.g., sowing, imbibition for germination) is then calculated by:

\begin{align} \mathrm{GDD}(T) &= \max \{ 0, \min \{ T, T_{opt} \} - T_b \} \\ \mathrm{cGDD} &= \sum_i^n \mathrm{GDD}(T_i) \\ \end{align}

In this section, we will create a model that simulates GDD and cGDD.

### System

Let us start by making a system called GrowingDegreeDay. This can be done using a simple Cropbox macro @system.

@system GrowingDegreeDay

#### Variables

From the equation, let's identify the variables we need to declare in our system. In the equation for GDD, we have two parameters Topt and Tb. Since they are fixed values, we will declare them as preserve variables, which are variables that remain constant throughout a simulation.

@system GrowingDegreeDay begin
Tb ~ preserve
To ~ preserve
end

Tb and To are parameters that we may want to change depending on the simulation. To make this possible, we will assign them the parameter tag, which allows the tagged variables to be altered through a configuration for each simulation. Note that we will not assign values at declaration because we will configure them when we run the simulation.

@system GrowingDegreeDay begin
Tb ~ preserve(parameter)
To ~ preserve(parameter)
end

Lastly, we will tag the variables with units. Tagging units is the recommended practice for many reasons, one of which is to catch mismatching units during calculations.

@system GrowingDegreeDay begin
Tb ~ preserve(parameter, u"°C")
To ~ preserve(parameter, u"°C")
end

In the GDD equation, T represents the average daily temperature value necessary to calculate the GDD. Likewise, the variable in our system will represent a series of daily average temperatures. The series of temperature values will be driven from an external data source, for which we will create a separate system later on for data extraction. For the GrowingDegreeDay system, we will declare T as a hold variable, which represents a placeholder that will be replaced by a T from another system.

@system GrowingDegreeDay begin
T ~ hold
Tb ~ preserve(parameter, u"°C")
To ~ preserve(parameter, u"°C")
end

We declared all the necessary variables required to calculate GDD. Now it is time to declare GDD as a variable in the system. Because GDD is a variable that we want to evaluate and store in each update, we will declare it as a track variable with T, Tb, and To as its depending variables.

@system GrowingDegreeDay begin
T ~ hold
Tb ~ preserve(parameter, u"°C")
To ~ preserve(parameter, u"°C")

GDD(T, Tb, To) => begin
min(T, To) - Tb
end ~ track(min = 0, u"K")
end

Note that we have tagged the unit for GDD as u"K". This is to avoid incompatibilities that u"°C" has with certain operations.

Now that GDD is declared in the system, we will declare cGDD as an accumulate variable with GDD as its depending variable. Recall that accumulate variables perform the Euler method of integration.

@system GrowingDegreeDay begin
T ~ hold
Tb ~ preserve(parameter, u"°C")
To ~ preserve(parameter, u"°C")

GDD(T, Tb, To) => begin
min(T, To) - Tb
end ~ track(min = 0, u"K")

cGDD(GDD) ~ accumulate(u"K*d")
end

We have declared all the necessary variables for GrowingDegreeDay.

#### Mixins

Now let's address the issue of the missing temperature values. We will make a new system that will provide the missing temperature data we need for simulating GrowingDegreeDay. We will call this system Temperature. The purpose of Temperature will be to obtain a time series of daily average temperature values from an external data source.

@system Temperature

For this tutorial, we will be using weather data from Beltsville, Maryland in 2002. The data is available here.

If you are unsure how to read CSV files, check the following documentation. You can also follow the template below (make sure path corresponds to the path to your data file).

using CSV, DataFrames

weather = CSV.read(path, DataFrame)

The data should look something like the following:

first(weather, 3)
3×10 DataFrame
Rowyeardoyrad (W/m^2)Tavg (°C)Tmax (°C)Tmin (°C)rainfall (mm)date (:Date)GDD (K)cGDD (K)
Int64Int64Float64Float64Float64Float64Int64DateFloat64Float64
12002135295.814.922.18.602002-05-156.96.9
22002136297.918.027.74.902002-05-1610.010.0
32002137224.221.327.314.332002-05-1713.313.3

Notice that the column names have units in parentheses. The unitfy() function in Cropbox automatically assigns units to values based on names of the columns (if the unit is specified).

weather = unitfy(weather)
first(weather, 3)
3×10 DataFrame
Int64Int64Quantity…Quantity…Quantity…Quantity…Quantity…DateQuantity…Quantity…
12002135295.8 W m^-214.9 °C22.1 °C8.6 °C0 mm2002-05-156.9 K6.9 K
22002136297.9 W m^-218.0 °C27.7 °C4.9 °C0 mm2002-05-1610.0 K10.0 K
32002137224.2 W m^-221.3 °C27.3 °C14.3 °C3 mm2002-05-1713.3 K13.3 K

In the Temperature system, there is one variable that we will declare before declaring any other variable. We will name this variable calendar.

@system Temperature begin
calendar(context) ~ ::Calendar
end

The purpose of calendar is to have access to variables inside the Calendar system such as init, last, and date, which represent initial, last, and current date, respectively.

Note

calendar is a variable reference to the Calendar system (one of the built-in systems of Cropbox), which has a number of time-related variables in date format. Declaring calendar as a variable of type Calendar allows us to use the variables inside the Calendar system as variables for our current system. Recall that context is a reference to the Context system and is included in every Cropbox system by default. Inside the Context system there is the config variable which references a Config object. By having context as a depending variable for calendar, we can change the values of the variables in calendar with a configuration.

The next variable we will add is a variable storing the weather data as a DataFrame. This variable will be a provide variable named data.

@system Temperature begin
calendar(context) ~ ::Calendar
data ~ provide(parameter, index=:date, init=calendar.date)
end

Note that we have tagged the variable with a parameter tag so that we can assign a DataFrame during the configuration. We will set the index of the extracted DataFrame as the "date" column of the data source. The init tag is used to specify the starting row of the data that we want to store. calendar.date refers to the date variable in the Calendar system, and is a track variable that keeps track of the dates of simulation. The initial value of date is dependent on calendar.init which we will assign during configuration. By setting init to calendar.date, we are making sure that the provide variable extracts data from the correct starting row corresponding to the desired initial date of simulation.

Now we can finally declare the temperature variable using one of the columns of the DataFrame represented by data. Because this variable is driven from a source, we will be declaring a drive variable named T. The from tag specifies the DataFrame source and the by tag specifies which column to take the values from.

@system Temperature begin
calendar(context) ~ ::Calendar
data ~ provide(parameter, index=:date, init=calendar.date)
T ~ drive(from=data, by=:Tavg, u"°C")
end
Main.Temperature

We finally have all the components to define our model. Because GrowingDegreeDay requires values for T from Temperature, let's redeclare GrowingDegreeDay with Temperature as a mixin. Because we want to run a simulation of GrowingDegreeDay, we also want to include Controller as a mixin. Recall that Controller must be included as a mixin for any system that you want to simulate.

@system GrowingDegreeDay(Temperature, Controller) begin
T ~ hold
Tb ~ preserve(parameter, u"°C")
To ~ preserve(parameter, u"°C")

GDD(T, Tb, To) => begin
min(T, To) - Tb
end ~ track(min = 0, u"K")

cGDD(GDD) ~ accumulate(u"K*d")
end
Main.GrowingDegreeDay

### Configuration

The next step is to create a configuration object to assign the values of parameters. Recall that data, T, Tb, and To are empty variables at the moment.

As covered in the Configuration section, we can make a single Config object with all the configurations we need for our systems.

Given the nature of GDD, this model is a daily model. To run a daily simulation, we need to configure the step variable in the Clock system from 1u"hr" to 1u"d". This will change the time interval of the simulation from hourly (default) to daily.

c = @config :Clock => :step => 1u"d"

Next we will add the configurations for GrowingDegreeDay. The only parameters we have to configure are Tb and To.

c = @config (
:Clock => (
:step => 1u"d"
),
:GrowingDegreeDay => (
:Tb => 8.0u"°C",
:To => 32.0u"°C"
)
)

Next we will pair the aforementioned DataFrame weather to data in Temperature

c = @config(
:Clock => (
:step => 1u"d"
),
:GrowingDegreeDay => (
:Tb => 8.0u"°C",
:To => 32.0u"°C"
),
:Temperature => (
:data => weather
)
)

Lastly, we will configure the init and last parameters of the Calendar system, which will define the time range of our simulation.

c = @config(
:Clock => (
:step => 1u"d"
),
:GrowingDegreeDay  => (
:Tb => 8.0u"°C",
:To => 32.0u"°C"
),
:Temperature => (
:data => weather
),
:Calendar => (
:init => ZonedDateTime(2002, 5, 15, tz"America/New_York"),
:last => ZonedDateTime(2002, 9, 30, tz"America/New_York")
)
)

Config for 4 systems:

 Clock step = 1 d
 GrowingDegreeDay Tb = 8.0 °C To = 32.0 °C
 Temperature data = 139×10 DataFrame…
 Calendar init = ZonedDateTime(2002, 5, 15, tz"America/New_York") last = ZonedDateTime(2002, 9, 30, tz"America/New_York")

### Simulation

Now that we have fully defined GrowingDegreeDay and created a configuration for it, we can finally simulate the model.

s = simulate(GrowingDegreeDay;
config = c,
stop = "calendar.stop",
index = "calendar.date",
target = [:GDD, :cGDD]
)

first(s, 10)
10×3 DataFrame
Rowcalendar.dateGDDcGDD
DateQuantity…Quantity…
12002-05-156.9 K0.0 d K
22002-05-1610.0 K6.9 d K
32002-05-1713.3 K16.9 d K
42002-05-184.5 K30.2 d K
52002-05-191.6 K34.7 d K
62002-05-202.1 K36.3 d K
72002-05-210.8 K38.4 d K
82002-05-223.6 K39.2 d K
92002-05-236.7 K42.8 d K
102002-05-2412.1 K49.5 d K

### Visualization

To end the tutorial, let's visualize the simulation using the plot() and visualize() functions.

We can input the DataFrame from our simulation in the plot() function to create a plot.

Here is a plot of GDD over time.

plot(s, "calendar.date", :GDD; kind=:line)

We can also simultaneously run a new simulation and plot its result using the visualize() function.

Here is a plot of cGDD over time.

visualize(GrowingDegreeDay, "calendar.date", :cGDD; config=c, stop="calendar.stop", kind=:line)

## Lotka-Volterra Equations

In this tutorial, we will create a model that simulates population dynamics between prey and predator using the Lotka-Volterra equations. The Lotka-Volterra equations are as follows:

\begin{align} \frac{dN}{dt} &= bN - aNP \\ \frac{dP}{dt} &= caNP - mP \\ \end{align}

Here is a list of variables used in the system:

SymbolValueUnitsDescription
t-$\mathrm{yr}$Time unit used in the model
N--Prey population as number of individuals (state variable)
P--Predator population as number of individuals (state variable)
b-$\mathrm{yr^{-1}}$Per capital birth rate that defines the intrinsic growth rate of prey population
a-$\mathrm{yr^{-1}}$Attack rate or predation rate
c--Conversion efficiency of an eaten prey into new predator; predator's reproduction efficiency per prey consumed)
m-$\mathrm{yr^{-1}}$Mortality rate of predator population

### System

Let's begin by creating a system called LotkaVolterra. Since this is a system that we want to simulate later on, we must include Controller as a mixin.

@system LotkaVolterra(Controller)

We will first declare a time variable with a yearly unit, which we will use for plotting the model simulations later on. Recall that context.clock.time is a variable that keeps track of the progression of time. We are simply declaring a variable to keep track of the time in years.

@system LotkaVolterra(Controller) begin
t(context.clock.time) ~ track(u"yr")
end

Next, we will declare the parameters in the equations as preserve variables. preserve variables are variables that remain constant throughout a simulation.

@system LotkaVolterra(Controller) begin
t(context.clock.time) ~ track(u"yr")

b: prey_birth_rate            ~ preserve(parameter, u"yr^-1")
a: predation_rate             ~ preserve(parameter, u"yr^-1")
c: predator_reproduction_rate ~ preserve(parameter)
m: predator_mortality_rate    ~ preserve(parameter, u"yr^-1")
end

Now let's declare the prey and predator populations as variables. The Lotka-Volterra equations describe the rates of change for the two populations. As we want to track the actual number of the two populations, we will declare the two populations as accumulate variables, which are simply Euler integrations of the two population rates. Note that a variable can be used as its own depending variable.

@system LotkaVolterra(Controller) begin
t(context.clock.time) ~ track(u"yr")

b: prey_birth_rate            ~ preserve(parameter, u"yr^-1")
a: predation_rate             ~ preserve(parameter, u"yr^-1")
c: predator_reproduction_rate ~ preserve(parameter)
m: predator_mortality_rate    ~ preserve(parameter, u"yr^-1")

N(N, P, b, a):    prey_population     =>     b*N - a*N*P ~ accumulate
P(N, P, c, a, m): predator_population => c*a*N*P -   m*P ~ accumulate
end

By default, accumulate variables initialize at a value of zero. In our current model, that would result in two populations remaining at zero indefinitely. To address this, we will define the initial values for the two accumulate variables using the init tag. We can specify a particular value, or we can also create and reference new parameters representing the two initial populations. We will go with the latter option as it allows us to flexibly change the initial populations with a configuration.

@system LotkaVolterra(Controller) begin
t(context.clock.time) ~ track(u"yr")

b: prey_birth_rate            ~ preserve(parameter, u"yr^-1")
a: predation_rate             ~ preserve(parameter, u"yr^-1")
c: predator_reproduction_rate ~ preserve(parameter)
m: predator_mortality_rate    ~ preserve(parameter, u"yr^-1")

N0: prey_initial_population     ~ preserve(parameter)
P0: predator_initial_population ~ preserve(parameter)

N(N, P, b, a):    prey_population     =>     b*N - a*N*P ~ accumulate(init=N0)
P(N, P, c, a, m): predator_population => c*a*N*P -   m*P ~ accumulate(init=P0)
end
Main.LotkaVolterra

### Configuration

With the system now defined, we will create a Config object to fill or adjust the parameters.

First, we will change the step variable in the Clock system to 1u"d", which will make the system update at a daily interval. Recall that Clock is a system that is referenced in all systems by default. You can technically run the model with any timestep.

lvc = @config (:Clock => :step => 1u"d")

Config for 1 system:

 Clock step = 1 d

Next, we will configure the parameters in the LotkaVolterra system that we defined. Note that we can easily combine configurations by providing multiple elements.

lvc = @config (lvc,
:LotkaVolterra => (
b = 0.6,
a = 0.02,
c = 0.5,
m = 0.5,
N0 = 20,
P0 = 30
)
)

Config for 2 systems:

 Clock step = 1 d
 LotkaVolterra b = 0.6 a = 0.02 c = 0.5 m = 0.5 N0 = 20 P0 = 30

### Visualization

Let's visualize the LotkaVolterra system with the configuration that we just created, using the visualize() function. The visualize() function both runs a simulation and plots the resulting DataFrame.

visualize(LotkaVolterra, :t, [:N, :P]; config = lvc, stop = 100u"yr", kind = :line)

### Density-Dependent Lotka-Volterra Equations

Now let's try to make a density-dependent version of the original Lotka-Volterra model which incorporates a new term in the prey population rate. The new variable K represents the carrying capacity of the prey population.

\begin{align} \frac{dN}{dt} &= bN-\frac{b}{K}N^2-aNP \\ \frac{dP}{dt} &= caNP-mP \\ \end{align}

#### System

We will call this new system LotkaVolterraDD.

@system LotkaVolterraDD(Controller)

Since we already defined the LotkaVolterra system, which already has most of the variables we require, we can use LotkaVolterra as a mixin for LotkaVolterraDD. This makes our task a lot simpler, as all that remains is to declare the variable K for carrying capacity and redeclare the variable N for prey population. The variable N in the new system will automatically overwrite the N from LotkaVolterra.

@system LotkaVolterraDD(LotkaVolterra, Controller) begin
N(N, P, K, b, a): prey_population => begin
b*N - b/K*N^2 - a*N*P
end ~ accumulate(init = N0)

K: carrying_capacity ~ preserve(parameter)
end
Main.LotkaVolterraDD

#### Configuration

Much like the new system, the new configuration can be created by reusing the old configuration. All we need to configure is the new variable K.

lvddc = @config(lvc, (:LotkaVolterraDD => :K => 1000))

Config for 3 systems:

 Clock step = 1 d
 LotkaVolterra b = 0.6 a = 0.02 c = 0.5 m = 0.5 N0 = 20 P0 = 30

#### Visualization

Once again, let's visualize the system using the visualize() function.

visualize(LotkaVolterraDD, :t, [:N, :P]; config = lvddc, stop = 100u"yr", kind = :line)

#### Calibration

If you want to calibrate the parameters according to a particular dataset, Cropbox provides the calibrate() function, which relies on BlackBoxOptim.jl for global optimization methods. If you are interested in local optimization methods, refer to Optim.jl package for more information.

For this tutorial, we will use a dataset containing the number of pelts (in thousands) of Canadian lynx and snowshoe hare traded by the Hudson Bay Trading Company in Canada from 1845 to 1935. The data is available here.

If you are unsure how to read CSV files, check the following documentation. You can also follow the template below (make sure path corresponds to the path to your data file).

using CSV, DataFrames

pelts = CSV.read(path, DataFrame)

The data should look something like the this:

first(pelts, 3)
3×3 DataFrame
RowYear (yr)HareLynx
Int64Float64Float64
1184519.5830.09
2184619.645.15
3184719.6149.15

Recall that we can use the unitfy() function in Cropbox to automatically assign units when they are specified in the column headers.

pelts = unitfy(pelts)
first(pelts, 3)
3×3 DataFrame
RowYearHareLynx
Quantity…Float64Float64
11845 yr19.5830.09
21846 yr19.645.15
31847 yr19.6149.15

Let's plot the data and see what it looks like.

visualize(pelts, :Year, [:Hare, :Lynx], kind = :scatterline)

For our calibration, we will use a subset of the data covering years 1900 to 1920.

pelts_subset = @subset(pelts, 1900u"yr" .<= :Year .<= 1920u"yr")
21×3 DataFrame
RowYearHareLynx
Quantity…Float64Float64
11900 yr12.827.13
21901 yr4.729.47
31902 yr4.7314.86
41903 yr37.2231.47
51904 yr69.7260.57
61905 yr57.7863.51
71906 yr28.6854.7
81907 yr23.376.3
91908 yr21.543.41
101909 yr26.345.44
111910 yr53.111.65
121911 yr68.4820.35
131912 yr75.5832.88
141913 yr57.9239.55
151914 yr40.9743.36
161915 yr24.9540.83
171916 yr12.5930.36
181917 yr4.9717.18
191918 yr4.56.82
201919 yr11.213.19
211920 yr56.63.52

Before we calibrate the parameters for LotkaVolterra, let's add one new variable to the system. We will name this variable y for year. The purpose of y is to keep track of the year in the same manner as the dataset.

@system LotkaVolterra(Controller) begin
t(context.clock.time) ~ track(u"yr")
y(t): year            ~ track::Int(u"yr", round)

b: prey_birth_rate            ~ preserve(parameter, u"yr^-1")
a: predation_rate             ~ preserve(parameter, u"yr^-1")
c: predator_reproduction_rate ~ preserve(parameter)
m: predator_mortality_rate    ~ preserve(parameter, u"yr^-1")

N0: prey_initial_population     ~ preserve(parameter)
P0: predator_initial_population ~ preserve(parameter)

N(N, P, b, a):    prey_population     =>     b*N - a*N*P ~ accumulate(init=N0)
P(N, P, c, a, m): predator_population => c*a*N*P -   m*P ~ accumulate(init=P0)
end
Main.LotkaVolterra

We will now use the calibrate() function to find parameters that fit the data. Keep in mind that the search range for each parameter will be determined by you. We will use the snap option to explicitly indicate that the output should be recorded by 365-day intervals to avoid excessive rows in the DataFrame causing unnecessary slowdown. Note that we will use 365u"d" instead of 1u"yr" which is technically equivalent to 365.25u"d" following the convention in astronomy. For information regarding syntax, please check the reference.

lvcc = calibrate(LotkaVolterra, pelts_subset;
index = :Year => :y,
target = [:Hare => :N, :Lynx => :P],
config = :Clock => (:init => 1900u"yr", :step => 1u"d"),
parameters = LotkaVolterra => (;
b = (0, 2),
a = (0, 2),
c = (0, 2),
m = (0, 2),
N0 = (0, 100),
P0 = (0, 100),
),
stop = 20u"yr",
snap = 365u"d"
)

Config for 1 system:

 LotkaVolterra b = 0.909952 a = 0.0346191 c = 0.634675 m = 0.681484 N0 = 12.7112 P0 = 15.4088

As you can see above, the calibrate() function will return a Config object for the system.

Using the new configuration, let's make a comparison plot to visualize how well the simualation with the new parameters fits the data.

p1 = visualize(pelts_subset, :Year, [:Hare, :Lynx]; kind = :scatterline)
visualize!(p1, LotkaVolterra, :t, [:N, :P];
config = (lvcc, :Clock => (:init => 1900u"yr", :step => 1u"d")),
stop = 20u"yr",
kind = :line,
colors = [1, 2],
names = [],
)

Now let's try calibrating the density-dependent version of the model. Since we made a slight change to LotkaVolterra, let's make sure to define LotkaVolterraDD again.

@system LotkaVolterraDD(LotkaVolterra, Controller) begin
N(N, P, K, b, a): prey_population => begin
b*N - b/K*N^2 - a*N*P
end ~ accumulate(init = N0)

K: carrying_capacity ~ preserve(parameter)
end
Main.LotkaVolterraDD

Don't forget to add K among the parameters that we want to calibrate.

lvddcc = calibrate(LotkaVolterraDD, pelts_subset;
index = :Year => :y,
target = [:Hare => :N, :Lynx => :P],
config = :Clock => (:init => 1900u"yr", :step => 1u"d"),
b = (0, 2),
a = (0, 2),
c = (0, 2),
m = (0, 2),
N0 = (0, 100),
P0 = (0, 100),
K = (0, 1000)
),
stop = 20u"yr",
snap = 365u"d"
)

Config for 1 system:

 LotkaVolterraDD b = 0.739399 a = 0.024501 c = 1.49285 m = 0.945061 N0 = 13.3898 P0 = 4.71093 K = 203.162

Once again, let us make a comparison plot to see how the density-dependent version of the model fares against the original dataset.

p2 = visualize(pelts_subset, :Year, [:Hare, :Lynx]; kind = :scatterline)
config = (lvddcc, :Clock => (:init => 1900u"yr", :step => 1u"d")),
stop = 20u"yr",
kind = :line,
colors = [1, 2],
names = [],
)

#### Evaluation

We have visualized how the simulated LotkaVolterra and LotkaVolterraDD systems compare to the the original dataset. Let us obtain a metric for how well the simulations fit the original dataset using the evaluate() function in Cropbox. The evaluate() function supports numerous different metrics for evaluation. Here, we will calculate the root-mean-square error (RMSE) and modeling efficiency (EF).

Here are the evaluation metrics for LotkaVolterra. The numbers in the tuples correspond to hare and lynx, respectively.

evaluate(LotkaVolterra, pelts_subset;
index = :Year => :y,
target = [:Hare => :N, :Lynx => :P],
config = (lvcc, :Clock => (:init => 1900u"yr", :step => 1u"d")),
stop = 20u"yr",
snap = 365u"d",
metric = :rmse,
)
(15.715248930777141, 11.762987336092772)
evaluate(LotkaVolterra, pelts_subset;
index = :Year => :y,
target = [:Hare => :N, :Lynx => :P],
config = (lvcc, :Clock => (:init => 1900u"yr", :step => 1u"d")),
stop = 20u"yr",
snap = 365u"d",
metric = :ef
)
(0.5507274075784151, 0.633879015444955)

Here are the evaluation metrics for LotkaVolterraDD:

evaluate(LotkaVolterraDD, pelts_subset;
index = :Year => :y,
target = [:Hare => :N, :Lynx => :P],
config = (lvddcc, :Clock => (:init => 1900u"yr", :step => 1u"d")),
stop = 20u"yr",
snap = 365u"d",
metric = :rmse
)
(20.047359944114255, 10.691595466753402)
evaluate(LotkaVolterraDD, pelts_subset;
index = :Year => :y,
target = [:Hare => :N, :Lynx => :P],
config = (lvddcc, :Clock => (:init => 1900u"yr", :step => 1u"d")),
stop = 20u"yr",
snap = 365u"d",
metric = :ef
)
(0.2688916022924638, 0.6975355070970235)