Optimization
MACS 33001
University of Chicago
\(f(x) = -x^2\), \(x \in [-3, 3]\)
\(f(x) = x^3\), \(x \in [-3, 3]\)
Example: MLE for a normal distribution
set.seed(1234)
y <- rnorm(n = 10000, mean = 0.25, sd = 10)
head(y)## [1] -11.82 3.02 11.09 -23.21 4.54 5.31# define the log-likelihood function for a normal distribution
log.like <- function(mu, sigma.2, y){
part1 <- - (1 / (2 * sigma.2)) * sum((y - mu)^2)
part2 <- - (length(y)/2) * log(sigma.2)
out <- part1 + part2
return(out)
}# define search space
(search_space <- expand.grid(mu = seq(-2, 2, by = .01),
sigma2 = seq(8^2, 12^2, by = .1)) %>%
as_tibble)## # A tibble: 321,201 x 2
## mu sigma2
## <dbl> <dbl>
## 1 -2 64
## 2 -1.99 64
## 3 -1.98 64
## 4 -1.97 64
## 5 -1.96 64
## 6 -1.95 64
## 7 -1.94 64
## 8 -1.93 64
## 9 -1.92 64
## 10 -1.91 64
## # ... with 321,191 more rows# evaluate parameters
system.time((search_space <- search_space %>%
mutate(logLik = map2_dbl(mu, sigma2, log.like, y = y))))## user system elapsed
## 11.87 3.59 15.65# which parameter values maximize the likelihood function
filter(search_space, logLik == max(logLik))## # A tibble: 1 x 3
## mu sigma2 logLik
## <dbl> <dbl> <dbl>
## 1 -1.55 142. Inf# draw a heatmap of the log-likelihood values
ggplot(search_space, aes(mu, sigma2, fill = logLik)) +
geom_raster() +
geom_contour(aes(z = logLik)) +
scale_fill_continuous(type = "viridis") +
labs(x = expression(mu),
y = expression(sigma^2))
mean(y)## [1] 0.311var(y)## [1] 97.5Implementation in R
\[y = x^3 + 2x - 5\]
# define the function
func2 <- function(x) {
x^3 - 2 * x - 5
}
# draw a plot of the function
range <- data_frame(x = c(-5, 5))
f0 <- range %>%
ggplot(aes(x)) +
stat_function(fun = func2, size = .5, linetype = 2) +
geom_hline(yintercept = 0) +
geom_vline(xintercept = 0) +
labs(x = expression(x),
y = expression(f(x)))
f0
uniroot(func2, interval = c(2, 3))## $root
## [1] 2.09
##
## $f.root
## [1] -0.000115
##
## $iter
## [1] 5
##
## $init.it
## [1] NA
##
## $estim.prec
## [1] 6.1e-05Implementation in R
# define the function
obj_func <- function(x) {
x^4 - 3 * x^3 + 2
}
# draw a plot of the function
range <- data_frame(x = c(-1, 3.25))
f0 <- range %>%
ggplot(aes(x)) +
stat_function(fun = obj_func, size = .5, linetype = 2) +
geom_hline(yintercept = 0) +
geom_vline(xintercept = 0) +
labs(x = expression(x),
y = expression(f(x)))
f0
Implementation in R
\[\nabla f(x) = 4x^3 - 9x^2 = 0\]
# define the gradient
gradient <- function(x) {
(4 * x^3) - (9 * x^2)
}
f0 <- range %>%
ggplot(aes(x)) +
stat_function(fun = gradient, size = .5, linetype = 2) +
geom_hline(yintercept = 0) +
geom_vline(xintercept = 0) +
labs(x = expression(x),
y = expression(f(x)))
f0
Example: one-dimensional utility
politician.utility <- function(policy.content){
politician.support <- -(policy.content - 1)^2 + 8
return(politician.support)
}range <- data_frame(x = c(-2, 4))
range %>%
ggplot(aes(x)) +
stat_function(fun = politician.utility, size = .5, linetype = 2) +
geom_hline(yintercept = 0) +
geom_vline(xintercept = 0) +
labs(title = "Legislator's Utility Function",
x = "Policy Ideology",
y = "Legislator Utility From Policy")
Example: two-dimensional utility
politician.utility.2d <- function(params){
# split parameters
economic.content <- params[1]
social.content <- params[2]
# calculate utility
politician.support <- (-(economic.content - 1)^2 + 8 ) + (-(social.content + 2)^2 + 8)
return(politician.support)
}# generate data values
expand.grid(economic.substance = seq(from=-8,to=8, by=.2),
social.substance = seq(from=-8, to=8, by=.2)) %>%
as_tibble %>%
mutate(utility = map2_dbl(economic.substance, social.substance,
~ politician.utility.2d(c(.x, .y)))) %>%
ggplot(aes(economic.substance, social.substance, fill = utility)) +
geom_raster() +
geom_contour(aes(z = utility)) +
scale_fill_continuous(type = "viridis") +
labs(title = "Utility Based on Policy Content on Two Dimensions",
x = "Economic Policy Content",
y = "Social Policy Content")