Ex. 1:

Write down Newton’s formula for finding the minimum of \(f(x)\) = \((3x^4 - 4x^3)\)/12 in the range of [-10, 10]. Then, implement it in R.

Solution:

First find the first two derivatives of f(x)

\(f'(x) = (12x^3 - 12x^2)/12 = x^2(x - 1) = x^3 - x^2\) \(f"(x) = 3x^2 - 2x\)

Per Newton’s formula \(x_k+1 = x_k - (x_k^3 - x_k^2)/(3x_k^2 - 2x_k)\)

The minimum is 1 as demonstrated below.

f <- function(x) {x - (x^3 - x^2) / (3 * x^2 - 2 * x)}
f1 <- function(X_0) {
  values <- c()
  for(i in 1:10) {
    if(i == 1){
      values[[i]] <- f(X_0)
    }
    else{
      values[[i]] <- f(values[[i - 1]])
    }
  }
  print(values)
}

f1(2)
## [[1]]
## [1] 1.5
## 
## [[2]]
## [1] 1.2
## 
## [[3]]
## [1] 1.05
## 
## [[4]]
## [1] 1.004348
## 
## [[5]]
## [1] 1.000037
## 
## [[6]]
## [1] 1
## 
## [[7]]
## [1] 1
## 
## [[8]]
## [1] 1
## 
## [[9]]
## [1] 1
## 
## [[10]]
## [1] 1
f1(5)
## [[1]]
## [1] 3.461538
## 
## [[2]]
## [1] 2.445307
## 
## [[3]]
## [1] 1.782962
## 
## [[4]]
## [1] 1.36611
## 
## [[5]]
## [1] 1.127755
## 
## [[6]]
## [1] 1.023598
## 
## [[7]]
## [1] 1.00104
## 
## [[8]]
## [1] 1.000002
## 
## [[9]]
## [1] 1
## 
## [[10]]
## [1] 1

Ex. 2:

Explore optimize() in R and try to solve the previous problem.

Solution:

The minimum value using the optimize function is close to 1 as seen below.

f <- function(x) {(3*x^4 - 4*x^3) / 12}
xmin <- optimize(f, interval = c(-10,10), tol = 0.0001)
xmin
## $minimum
## [1] 0.9999986
## 
## $objective
## [1] -0.08333333

Ex. 3:

Use any optimization algorithm to find the minimum of f(x, y) = (x - 1)^2 + 100(y - x2)2 in the domain -10 <= x, y <= 10. Discuss any issues concerning the optimization process.

Solution:

We will use the Newton method for a higher dimension to solve this problem.

\(x_t+1 = x_t - H^-1 V f(x,y)\) \(x_0\) = \[\begin{bmatrix} x\\y \end{bmatrix}\] = \[\begin{bmatrix} 0\\0 \end{bmatrix}\] Starting value $V f(x,y) = \[\begin{bmatrix} (x-1)^2 \\ 100(y - x^2)^2 \end{bmatrix}\]

We find the Hessian matrix by finding the second derivative and plugging in starting values.

\(f_xx = 1200*x^2 + 2 − 400*y\) \(f_xy = -400*x\) \(f_yy = 200\) \(f_yx = -400*x\)

\(H = V^2 f(x,y)\) = \[\begin{bmatrix}2 & -0\\ -0 & 200\end{bmatrix}\] \(H^-1\) = \[\begin{bmatrix}0.5 & 0\\ 0 & 0.005\end{bmatrix}\]

With t = 0, \(x_1 = x_0 - H^-1 V f(x,y)\)

\(x_1\) = \[\begin{bmatrix}0 & 0\end{bmatrix}\] \(-\) \[\begin{bmatrix}0.5 & 0\\ 0 & 0.005\end{bmatrix}\]

\begin{bmatrix}1 & 0\end{bmatrix}

\(x_1\) = \[\begin{bmatrix}0 & 0\end{bmatrix}\] \(-\) \[\begin{bmatrix}0.5 & 0\end{bmatrix}\] \(x_1\) = \[\begin{bmatrix}-0.5 & 0\end{bmatrix}\]

From here f(x, y) will converge into a single point as you go through more values of of t.

Ex. 4:

Explore the optimr package for R and try to solve the previous problem.

library(optimr)

f_xy <- function(param) {(param[1] - 1)^2 + 100 * (param[2] - param[1]^2)^2}

options(scipen = 999)
optimr(c(-10,10), f_xy, method = "L-BFGS-B")
## $par
## [1] 0.9997999 0.9995999
## 
## $value
## [1] 0.00000004002673
## 
## $counts
## function gradient 
##       71       71 
## 
## $convergence
## [1] 0
## 
## $message
## [1] "CONVERGENCE: REL_REDUCTION_OF_F <= FACTR*EPSMCH"