# ---
# jupyter:
#   jupytext:
#     text_representation:
#       extension: .r
#       format_name: light
#       format_version: '1.5'
#       jupytext_version: 1.14.4
#   kernelspec:
#     display_name: R
#     language: R
#     name: ir
# ---

# # Extended Technology Application

# ## *Maximum Sustainable Harvest*

# In certain situations, biologists are able to determine what is called a *reproduction curve*. This is a function
#
# > $y = f(P)$
#
# such that if $P$ is the population after $P$ years, then $f(P)$ is the population a year later, at time $t + 1$. Such a curve is shown below.

# &emsp; The line $y = P$ is significant because if it ever coincides with the curve $f(P)$, then we know that the population stays the same from year to year. Here the graph of $f$ lies mostly above the line, indicating that the population is increasing.  
# &emsp; Too many deer in a forest can deplete the food supply and eventually cause the population to decrease for lack of food. Often in such cases and with some controversy, hunters are allowed to "harvest" some of the deer. Then with a greater food supply, the remaining deer population might actually prosper and increase.  
# &emsp; We know that a population $P$ will grow to a population $f(P)$ in a year. If this were a population of fur-bearing animals and the population were increasing, then hunters could "harvest" the amount
#
# > $f(P) - P$
#
# each year without shrinking the initial population $P$. If the population were remaining the same or decreasing, then such a harvest would deplete the population.  
# &emsp; Suppose that we want to know the value of $P_0$ that would allow the harvest to be the largest. If we could determine that $P_0$, we could let the population grow until it reached that level and then begin harvesting year after year the amount $f(P_0) - P_0$.

# &emsp; Let the harvest function $H$ be given by
#
# > $H(P) = f(P) - P$.
#
# Then
# > $H'(P) = f'(P) - 1$.
#
# Now, if we assume that $H'(P)$ exists for all values of $P$ and that there is only one critical value, it follows that the *maximum sustainable harvest* occurs at that value $P_0$ such that
#
# > $H'(P_0) = f'(P_0) - 1 = 0$
#
# and
# > $H''(P_0) = f''(P_0) < 0$.
#
# Or, equivalently, we have the following.
# > `THEOREM`  
# > ---
# > The **maximum sustainable harvest** occurs at $P_0$ such that  
# >
# > &emsp; &emsp; $f'(P_0) = 1$ and $f''(P_0) < 0$,  
# > 
# > is given by  
# > 
# > &emsp; &emsp; $H(P_0) = f(P_0) - P_0$.

# ## EXERCISES
#
# For each reproduction curve in Exercises 1&ndash;3, do the following.
#
# 1. Graph the reproduction function, the line $y = P$, and the harvest function using the same viewing window.
# 2. Find the population at which the maximum sustainable harvest occurs. Use both a graphical solution and a calculus solution.
# 3. Find the maximum sustainable harvest.

# ### 1. $f(P) = P(10 - P)$, where $P$ is measured in thousands.

# #### 1. Graphing

f <- function (P) P * (10 - P)
y <- function (P) P
H <- function (P) f(P) - P
exercise1.plot <- function (xlim.max = 15, ylim.max = 30) {
    plot(f, xlim = c(0, xlim.max), ylim = c(0, ylim.max), col = 2, lwd = 2, xlab = expression("P"), ylab = expression("f(P)")) # red
    plot(y, xlim = c(0, xlim.max), ylim = c(0, ylim.max), col = 3, lwd = 2, add = TRUE) # green
    plot(H, xlim = c(0, xlim.max), ylim = c(0, ylim.max), col = 4, lwd = 2, add = TRUE) # blue
    grid(lty = 2)
    legend("topright",
           legend = c(expression(f), expression(y), expression(H)),
           col = c(2, 3, 4),
           lty = c(1, 1, 1)
          )
}
exercise1.plot()

# #### 2. Find population at which maximum sustainable harvest occurs
#
# ##### Graphical approach:

exercise1.plot()
max <- optimize(H, c(0, 100), maximum = TRUE)
points(max$maximum, max$objective, col = 4, lwd = 2)
print(max)

# ##### Calculus approach:
#
# > $f(P) = P(10 - P)$  
# > $f'(P) = 10 - 2P$   
# > $H(P) = f(P) - P$  
# > $H'(P) = f'(P) - 1$  
# > $H'(P) = 10 - 2P - 1$  
# > $H'(P) = -2P + 9$  
# > $0 = -2P + 9$  
# > $\color{blue}{P = 4.5}$

# +
df <- function (P) 10 - (2 * P)
exercise1.plot()
a <- max$maximum[1]
max.x <- 15

points(a, f(a), col = 2, lwd = 2)

# Draw tangent
slope <- function (x) df(a) * (x - a) + f(a)
lines(c(0, max.x), c(slope(0), slope(max.x)), lty = 5)
text(a - a / 10, f(a) + a / 3, "m = 1")

# Draw Maximum Sustainable Harvest
lines(c(a, a), c(0, y(a)), lty = 5)
lines(c(a, a), c(y(a), f(a)), lwd = 2, col = 6)
text(a + a * 0.8, (f(a) - y(a)) - a * 0.75, "Maximum sustainable harvest", col = 6)
print(f(a))
# -

# ### 2. $f(P) = -0.025P^2 + 4P$, where $P$ is measured in thousands. This is the reproduction curve in the Hudson Bay area for the snowshoe hare, a fur-bearing animal.

# #### 1. Graphing 

f <- function (P) (-0.025 * P ^ 2) + (4 * P)
y <- function (P) P
H <- function (P) f(P) - P
exercise2.plot <- function (xlim.max = 175, ylim.max = 175) {
    plot(f, xlim = c(0, xlim.max), ylim = c(0, ylim.max), col = 2, lwd = 2, xlab = expression("P"), ylab = expression("f(P)")) # red
    plot(y, xlim = c(0, xlim.max), ylim = c(0, ylim.max), col = 3, lwd = 2, add = TRUE) # green
    plot(H, xlim = c(0, xlim.max), ylim = c(0, ylim.max), col = 4, lwd = 2, add = TRUE) # blue
    grid(lty = 2)
    legend("topright",
           legend = c(expression(f), expression(y), expression(H)),
           col = c(2, 3, 4),
           lty = c(1, 1, 1)
          )
}
exercise2.plot()

# #### 2. Find population at which maximum sustainable harvest occurs
#
# ##### Graphical approach:

exercise2.plot()
max <- optimize(H, c(0, 100), maximum = TRUE)
points(max$maximum, max$objective, col = 4, lwd = 2)
print(max)

# ##### Calculus approach:
#
# > $f(P) = -0.025P^2 + 4P$  
# > $f'(P) = -0.05P + 4$   
# > $H(P) = f(P) - P$  
# > $H'(P) = f'(P) - 1$  
# > $H'(P) = -0.05P + 4 - 1$  
# > $H'(P) = -0.05P + 3$  
# > $0 = -0.05P + 3$  
# > $\color{blue}{P = 60}$

# #### 3. Maximum sustainable harvest

# +
df <- function (P) -0.05 * P + 4 
exercise2.plot()
a <- max$maximum[1]
max.x <- 175

points(a, f(a), col = 2, lwd = 2)

# Draw tangent
slope <- function (x) df(a) * (x - a) + f(a)
lines(c(0, max.x), c(slope(0), slope(max.x)), lty = 5)
text(a - a / 10, f(a) + a / 5, "m = 1")

# Draw Maximum Sustainable Harvest
lines(c(a, a), c(0, y(a)), lty = 5)
lines(c(a, a), c(y(a), f(a)), lwd = 2, col = 6)
text(a + a * 0.65, (f(a) - y(a)) + a * 0.5, "Maximum sustainable harvest", col = 6)
print(f(a))
# -

# ### 3. $f(P) = -0.01P^2 + 2P$, where $P$ is measured in thousands. This is the reproduction curve in the Hudson Bay area for the lynx, a fur-bearing animal.

# #### 1. Graphing

f <- function (P) (-0.01 * P ^ 2) + (2 * P)
y <- function (P) P
H <- function (P) f(P) - P
exercise3.plot <- function (xlim.max = 200, ylim.max = 150) {
    plot(f, xlim = c(0, xlim.max), ylim = c(0, ylim.max), col = 2, lwd = 2, xlab = expression("P"), ylab = expression("f(P)")) # red
    plot(y, xlim = c(0, xlim.max), ylim = c(0, ylim.max), col = 3, lwd = 2, add = TRUE) # green
    plot(H, xlim = c(0, xlim.max), ylim = c(0, ylim.max), col = 4, lwd = 2, add = TRUE) # blue
    grid(lty = 2)
    legend("topright",
           legend = c(expression(f), expression(y), expression(H)),
           col = c(2, 3, 4),
           lty = c(1, 1, 1)
          )
}
exercise3.plot()

# #### 2. Find population at which maximum sustainable harvest occurs
#
# ##### Graphical approach:

exercise3.plot()
max <- optimize(H, c(0, 100), maximum = TRUE)
points(max$maximum, max$objective, col = 4, lwd = 2)
print(max)

# ##### Calculus approach:
#
# > $f(P) = -0.01P^2 + 2P$  
# > $f'(P) = -0.02P + 2$   
# > $H(P) = f(P) - P$  
# > $H'(P) = f'(P) - 1$  
# > $H'(P) = -0.02P + 2 - 1$  
# > $H'(P) = -0.02P + 1$  
# > $0 = -0.02P + 1$  
# > $\color{blue}{P = 50}$

# #### 3. Maximum Sustainable Harvest

# +
df <- function (P) -0.02 * P + 2
exercise3.plot()
a <- max$maximum[1]
max.x <- 200

points(a, f(a), col = 2, lwd = 2)

# Draw tangent
slope <- function (x) df(a) * (x - a) + f(a)
lines(c(0, max.x), c(slope(0), slope(max.x)), lty = 5)
text(a - a / 10, f(a) + a / 5, "m = 1")

# Draw Maximum Sustainable Harvest
lines(c(a, a), c(0, y(a)), lty = 5)
lines(c(a, a), c(y(a), f(a)), lwd = 2, col = 6)
text(a + a * 0.95, (f(a) - y(a)) + a * 0.85, "Maximum sustainable harvest", col = 6)
print(f(a))
# -

# For each reproduction curve in Exercises 4 and 5, do the following.
#
# 1. Graph the reproduction function, the line $y = P$, and the harvest function using the same viewing window.
# 2. Graphically determine the population at which the maximum sustainable harvest occurs.
# 3. Find the maximum sustainable harvest.

# ### 4. $f(P) = 40\sqrt{P}$, where $P$ is measured in thousands. Assume that this is the reproduction curve for the brown trout population in a large lake.

#

f <- function (P) (40 * sqrt(P))
y <- function (P) P
H <- function (P) f(P) - P
exercise4.plot <- function (xlim.max = 2400, ylim.max = 2400) {
    plot(f, xlim = c(0, xlim.max), ylim = c(0, ylim.max), col = 2, lwd = 2, xlab = expression("P"), ylab = expression("f(P)")) # red
    plot(y, xlim = c(0, xlim.max), ylim = c(0, ylim.max), col = 3, lwd = 2, add = TRUE) # green
    plot(H, xlim = c(0, xlim.max), ylim = c(0, ylim.max), col = 4, lwd = 2, add = TRUE) # blue
    grid(lty = 2)
    legend("topright",
           legend = c(expression(f), expression(y), expression(H)),
           col = c(2, 3, 4),
           lty = c(1, 1, 1)
          )
}
exercise4.plot()

# #### 2. Find population at which maximum sustainable harvest occurs

exercise4.plot()
max <- optimize(H, c(0, 1500), maximum = TRUE)
points(max$maximum, max$objective, col = 4, lwd = 2)
print(max)

# +
df <- function (P) (20 / sqrt(P))
exercise4.plot()
a <- max$maximum[1]
max.x <- 2400

points(a, f(a), col = 2, lwd = 2)

# Draw tangent
slope <- function (x) df(a) * (x - a) + f(a)
lines(c(0, max.x), c(slope(0), slope(max.x)), lty = 5)
text(a - a / 10, f(a) + a / 5, "m = 1")

# Draw Maximum Sustainable Harvest
lines(c(a, a), c(0, y(a)), lty = 5)
lines(c(a, a), c(y(a), f(a)), lwd = 2, col = 6)
text(a + a * 1.35, (f(a) - y(a)) + a * 0.75, "Maximum sustainable harvest", col = 6)
print(f(a))
# -

# ### 5. $f(P) = 0.237P\sqrt{2000 - P^2}$, where $P$ is measured in thousands.

f <- function (P) (0.237 * P) * (sqrt((2000 - (P ^ 2))))
y <- function (P) P
H <- function (P) f(P) - P
exercise5.plot <- function (xlim.max = 60, ylim.max = 240) {
    plot(f, xlim = c(0, xlim.max), ylim = c(0, ylim.max), col = 2, lwd = 2, xlab = expression("P"), ylab = expression("f(P)")) # red
    plot(y, xlim = c(0, xlim.max), ylim = c(0, ylim.max), col = 3, lwd = 2, add = TRUE) # green
    plot(H, xlim = c(0, xlim.max), ylim = c(0, ylim.max), col = 4, lwd = 2, add = TRUE) # blue
    grid(lty = 2)
    legend("topright",
           legend = c(expression(f), expression(y), expression(H)),
           col = c(2, 3, 4),
           lty = c(1, 1, 1)
          )
}
exercise5.plot()

# #### 2. Find population at which maximum sustainable harvest occurs

exercise5.plot()
max <- optimize(H, c(0, 44), maximum = TRUE)
points(max$maximum, max$objective, col = 4, lwd = 2)
print(max)

# > $f(P) = 0.237P\sqrt{2000 - P^2}$  
# > $f(P) = 0.237P(2000-P^2)^{\frac{1}{2}}$  
# > $u = (2000-P^2)^{\frac{1}{2}}$  
# > $u' = {\frac{-P}{\sqrt{2000-P^2}}}$  
# > $f'(P) = 0.237P \cdot u' + u \cdot 0.237$  
# > $f'(P) = 0.237P \cdot {\frac{-P}{\sqrt{2000 - P^2}}} + \sqrt{2000 - P^2} \cdot 0.237$  
# > $f'(P) = \frac{-0.237P^2}{\sqrt{2000 - P^2}} + 0.237\sqrt{2000 - P^2}$  
# > $f'(P) = \frac{-0.237P^2 + 0.237(2000 - P^2)}{\sqrt{2000 - P^2}}$  
# > $f'(P) = \frac{-0.474P^2 + 474}{\sqrt{2000 - P^2}}$

# +
df <- function (P) ((-0.474 * P ^ 2 + 474) / sqrt(2000 - P ^ 2))
exercise5.plot()
a <- max$maximum[1]
max.x <- 60

points(a, f(a), col = 2, lwd = 2)

# Draw tangent
slope <- function (x) df(a) * (x - a) + f(a)
lines(c(0, max.x), c(slope(0), slope(max.x)), lty = 5)
text(a - a / 10, f(a) + a / 5, "m = 1")

# Draw Maximum Sustainable Harvest
lines(c(a, a), c(0, y(a)), lty = 5)
lines(c(a, a), c(y(a), f(a)), lwd = 2, col = 6)
text(a + a * 0.45, (f(a) - y(a)) - a, "Maximum sustainable harvest", col = 6)
print(f(a))
# -

# ### 6. The table below lists data regarding the reproduction of a certain animal.

# Population (in thousands)
current.population <- c(10.0, 20.0, 30.0, 40.0, 50.0)
next.population    <- c(09.7, 23.1, 37.4, 46.2, 42.6)

# #### a. Use REGRESSION to fit a cubic polynomial to these data.

model <- lm(next.population ~ current.population + I(current.population ^ 2) + I(current.population ^ 3))
print(model)

# #### b. Graph the reproduction function, the line $y = P$, and the harvest function using the same viewing window.

f <- function (P) model$coefficients[1] + model$coefficients[2] * P + model$coefficients[3] * P ^ 2 + model$coefficients[4] * P ^ 3
y <- function (P) P
H <- function (P) f(P) - P
exercise6.plot <- function (xlim.max = 60, ylim.max = 60) {
    plot(f, xlim = c(0, xlim.max), ylim = c(0, ylim.max), col = 2, lwd = 2, xlab = expression("P"), ylab = expression("f(P)")) # red
    plot(y, xlim = c(0, xlim.max), ylim = c(0, ylim.max), col = 3, lwd = 2, add = TRUE) # green
    plot(H, xlim = c(0, xlim.max), ylim = c(0, ylim.max), col = 4, lwd = 2, add = TRUE) # blue
    grid(lty = 2)
    legend("topright",
           legend = c(expression(f), expression(y), expression(H)),
           col = c(2, 3, 4),
           lty = c(1, 1, 1)
          )
}
exercise6.plot()

# #### c. Graphically determine the population at which the maximum sustainable harvest occurs.

# +
df <- function (P) -0.003324 * P ^ 2 + 0.143072 * P - 0.03381
exercise6.plot()

max <- optimize(H, c(0, 45), maximum = TRUE)
points(max$maximum, max$objective, col = 4, lwd = 2)

a <- max$maximum[1]
max.x <- 60

points(a, f(a), col = 2, lwd = 2)

# Draw tangent
slope <- function (x) df(a) * (x - a) + f(a)
lines(c(0, max.x), c(slope(0), slope(max.x)), lty = 5)
text(a - a / 10, f(a) + a / 5, "m = 1")

# Draw Maximum Sustainable Harvest
lines(c(a, a), c(0, y(a)), lty = 5)
lines(c(a, a), c(y(a), f(a)), lwd = 2, col = 6)
text(a + a * 0.4, (f(a) - y(a)) + a * 0.9, "Maximum sustainable harvest", col = 6)
print(f(a))

StarcoderdataJupyterScriptsDedupFiltered