In michaelfrancenelson/rmd.utils:

require(here)
require(rmd.utils)
knitr::opts_chunk$set(echo = TRUE)

exp_fun = function(x, a, b) { return(a * exp(-b * x))}
ricker_fun = function(x, a, b) { return(a * x * exp(-b * x))}

Learning objectives, technical sills and concepts:

Random number generation
Distribution Functions
Building custom functions
Fitting deterministic models
Calculating residuals
Defining and quantifying error

The Ricker Function

The Ricker funciton is a modification fo the exponential function. The formula is:

$f(x) = a \times x \times e^{-b \times x}$

You can see the general form of the curve (with some hints as to how to visually parameterize it) in McGarigal's figure:

knitr::include_graphics(find_file("mcGarigal_ricker_fun.PNG"))

Notice that the x-coordinate of the peak depends on the parameter b, the initial slope depends on the parameter a, while the height of the maximum depends on both parameters.

We can create a Ricker function in R:

ricker_fun = function(x, a, b) 
{
  return(a * x * exp(-b * x))
}

Here's a plot of the shape using the simplest possible parameters: a = 1 and b = 1:

curve(
  ricker_fun(x, 1, 1), 
  from = 0, to = 5, add = FALSE, 
  main = "Ricker function: a = 1, b = 1",
  ylab = "f(x)", xlab = "x")

Note the use of the curve() function:
The from and to arguments tells curve() what range of x-values it should include in the plot.
The add = FALSE argument value creates a new plot.

The exponential function

The exponential function is the building block of the Ricker.

Here's McGarigal's infographic:

knitr::include_graphics(find_file("mcGarigal_exponential_fun.PNG"))

We can create an exponential function in R:

I'll let you build it! You can use the Ricker function code I wrote above as a template.

exp_fun = function(x .......

NOTE: McGarigal defines the Ricker funciton with a negative sign before the parameter `b` but no negative sign before `b` in the exponential. For consistency, you need to include the negative sign with parameter `b` in your exponential function definition.

Plot of the exponential function

curve(
  exp_fun(x, 0.1, 0.5), 
  from = 0, to = 10, add = FALSE, 
  main = "Exponential function",
  ylab = "f(x)", xlab = "x")

What code could you use to create a similar plot?

Simulating data with different stochastic distributions

We've used curve() to plot deterministic models using linear, logistic, Ricker, and exponential functions using parameter values we selected.

We can also create simulated data that use these functions a determinsitic model and r's random deviates functions to generate the noise or 'errors'.

Simulated data on a line.

First, let's choose 50 uniformly-distributed random x-values within the interval 2 to 10:

# Seed the RNG so we can reproduce our results
set.seed(1234567)

# Specify the x-range and number of poitns:
n_pts = 50
x_min = 2
x_max = 10

# Generate the x-values
x_sim = runif(n_pts, min = x_min, max = x_max)

Next, we can choose an intercept and slope for our determinsitic model and generate the 'predicted' y values:

param_intercept = 2.3
param_slope = 0.67
y_pred = param_intercept + x_sim * param_slope
plot(x_sim, y_pred)

Normal errors 1

Now we can add some normally-distributed noise to generate our 'observed' y-values:

error_mean = 0
error_sd = 0.25

y_observed = 
  y_pred + 
  rnorm(
    n = n_pts, 
    mean = error_mean, 
    sd = error_sd)
plot(x_sim, y_observed)

Normal errors 2

We could also use a more sophisticated stochastic model. For example, we could make the variablity larger with increasing values of x:

error_mean = 0
error_sd = 0.1

y_observed_2 = 
  y_pred + 
  rnorm(
    n = n_pts, 
    mean = error_mean, 
    sd = error_sd * x_sim)
plot(x_sim, y_observed)
plot(x_sim, y_observed_2)

What did I do differently to make the points more variable as x increased?

Exponentially-distributed errors.

To produce this plot, I generated exponentially-distributed errors using rexp(). I used a rate parameter of 1.2. I'll let you think about how you might create such a plot yourselves.

y_observed_3 = 
  y_pred +
  rexp(n_pts, rate = 1.2) - 1 / 1.2
plot(x_sim, y_observed)
plot(x_sim, y_observed_3)

Choosing a model

Consider the three plots

par(mfrow = c(3, 1))
plot(x_sim, y_observed)
plot(x_sim, y_observed_2)
plot(x_sim, y_observed_3)

If you didn't know the process that generated these data, how would you choose deterministic and stochastic models?

The deterministic functions look linear, so we'll focus ont he stochastic models.
A first step is to examine histograms of the residuals:

par(mfrow = c(3, 1))
hist(y_observed - y_pred, main = "sim data 1", xlab = "observed y=values")
hist(y_observed_2 - y_pred, main = "sim data 2", xlab = "observed y=values")
hist(y_observed_3 - y_pred, main = "sim data 3", xlab = "observed y=values")

We'll apply these concepts to a real data set below.

Marbled salamander dispersal data

Background

The data for this exercise represent the dispersal of juvenile marbled salamanders from their natal ponds to neighboring ponds.

The data were derived from a long-term study of marbled salamanders in western Massachusetts in which a cluster of 14 vernal pools were monitored continuously between 1999-2004. All juveniles were marked upon leaving their natal ponds. Subsequent recaptures at non-natal ponds were used to determine dispersal rates between ponds for first-time breeders (ftb).

In the data set provided, the dispersal rates have been standardized to account for several factors, including the propensity for dispersal from each pond and the available distances between ponds owing to the unique configuration of ponds. The data set includes three variables:

dist.class = distance class, based on 100 m intervals;
disp.rate.ftb = standardized dispersal rate for first-time breeders, which can be interpreted as a relative dispersal probability.
disp.rate.eb = standardized dispersal rate for experienced breeders, which can be interpreted as arelative dispersal probability.

Obtain the data file here:

https://michaelfrancenelson.github.io/eco_602_634_2020/data/salamander_dispersal.csv

You'll try fitting three deterministic functions to the data:

linear
exponential
Ricker

Instructions

Your assignment is to examine the data set provided (dispersal.csv) and find at least three alternative mathematical functions to describe the apparent relationship between juvenile dispersal and distance. The overall goals are for you to

Examine the data set. Plot the relationship between juvenile dispersal (disp.rate.ftb) and distance (dist.class).
Visually parameterize linear, exponential, and Ricker functions via trial-and-error to fit the data.
Note: Be patient, we'll eventually use R's power to automatically fit model parameters. For now, the goal is graphical intuition.
Calculate the residuals and propose a distribution to model the errors.

Lab Questions

question_source_files = find_file("Q", search_path = find_file("lab_05_un", directory = TRUE), return_all = TRUE)
tmp = add_moodle_quiz_questions(question_source_files, include_solution = TRUE, include_metadata = FALSE)