library(tidyverse)
library(Deriv)
library(scales)
# Formula for polynomial curve: 60x^4 + 5x + 6
#
# I invented this by using WolframAlpha to fit a polynomial linear model that
# when through a set of specific points using this command (with a maximum
# exponent of 4):
#
# LinearModelFit[{{0, 6}, {0.5, 12}, {-0.8, 24}, {-0.9, 48}}, {1, x, x^4}, x]
#
# It yielded a messy model with like 78.3x^4, 10.1x + 3.44, so I put that in
# Desmos.com, moved the y intercept to 6, and tinkered with the other coefficients
# until it looked okay. Finally, I switched x and y so it's horizontal
response_curve <- function(y) 60*y^4 + 5*y + 6
response_deriv <- Deriv(response_curve)
# We can use these functions to find the x value (wage) that corresponds with a
# given y value (level of effort)
response_curve(0.5)
#> [1] 12.25
# We can also find the slope of the response curve at that point:
response_deriv(0.5)
#> [1] 35
# Armed with a slope an a single point (here slope = 35, point = (12.25, 0.5)),
# we can create a formula for the line that is tangent to the response curve at
# that point.
#
# Formula for y-intercept give a slope and a single point:
# y - y1 = m(x = x1)
#
# Because we've flipped the curve sidways and switched x and y, we need to
# rearrange this in terms of x:
#
# x = (y - y1) / m + x1
#
tangent_intercept <- function(x, m, x1, y1) (x - x1) / m + y1
# Here we invert tangent_intercept() by dividing by the negative derivative
# I DON'T KNOW WHY THIS WORKS BUT IT DOES
flipped_intercept <- function(x, effort) {
# m needs to be inverted because of switched x/y thing, so it's 1/deriv
tangent_intercept(x, m = 1 / response_deriv(effort),
x1 = effort, y1 = response_curve(effort)) /
# THIS NEGATIVE DIVISION IS A MYSTERY
-response_deriv(effort)
}
# Create a line based on the slope and intercept of the flipped derivative
tangent_line <- function(x, effort) {
(1 / response_deriv(effort) * x) + flipped_intercept(0, effort)
}
# FINALLY we can plot the flipped polynomial curve. This should, in theory, be
# as simple as rearranging our x = 60y^4 + 5y + 6 formula to be in terms of y
# and then using stat_function() to plot it, BUT polynomials are demonic beasts
# and doing this is beyond my abilities.
#
# GO HERE and see what I mean:
# https://www.wolframalpha.com/input/?i=x+%3D+60y%5E4+%2B+5y+%2B+6+in+terms+of+y
#
# So to cheat, we make a little data frame with a range of x and y values and
# then plot that with geom_line()
flipped_curve <- tibble(y = seq(0, 1, 0.01)) %>%
mutate(x = response_curve(y))
# Finally, we can plot all these things!
ggplot(data = tibble(x = 0:50), aes(x = x)) +
# Add flipped polynomial curve
geom_line(data = flipped_curve, aes(x = x, y = y),
size = 0.5, color = "grey30") +
# Add tangent lines
# If I wanted to be super fancy, I could figure out the width of the base of
# the triangle formed by each tangent line so that I could get ± the same x
# distance to the left and right of each point and have each tangent segment
# use the same horizontal distance, but that sounds too hard, so I just picked
# a bunch of numbers (here 5, 2, and 7) until it looked okay
stat_function(fun = tangent_line, args = list(effort = 0.5),
color = "darkred", size = 1,
xlim = c(response_curve(0.5) - 5, response_curve(0.5) + 5)) +
stat_function(fun = tangent_line, args = list(effort = 0.25),
color = "darkred", size = 1,
xlim = c(response_curve(0.25) - 2, response_curve(0.25) + 2)) +
stat_function(fun = tangent_line, args = list(effort = 0.8),
color = "darkred", size = 1,
xlim = c(response_curve(0.8) - 7, response_curve(0.8) + 7)) +
# Add dots at tangent points
annotate(geom = "point", x = response_curve(0.5), y = 0.5) +
annotate(geom = "point", x = response_curve(0.25), y = 0.25) +
annotate(geom = "point", x = response_curve(0.8), y = 0.8) +
# Labels and scale stuff
labs(x = "Wage per hour", y = "Work effort from employee") +
scale_x_continuous(labels = dollar, expand = c(0, 0)) +
scale_y_continuous(labels = percent_format(accuracy = 1), expand = c(0, 0)) +
coord_cartesian(xlim = c(0, 50), ylim = c(0, 1)) +
theme_minimal()
Created on 2019-02-21 by the reprex package (v0.2.1)