Optimizing

Examples

Example

Find any local maximum, minimum, or saddle points of the function f(x,y)=3x-(x^3)-2(y^2)+(y^4).

Solution

f[x_,y_]=3x-(x^3)-2(y^2)+(y^4)

Solve the equations where first partial derivatives are zero:

sol=Solve[{D[f[x,y],x]==0,D[f[x,y],y]==0},{x,y}]

Compute the Hessian:

D[f[x,y],x,x]
D[f[x,y],x,y]
D[f[x,y],y,y]
Hess[x_,y_]=D[f[x,y],x,x]*D[f[x,y],y,y]-((D[f[x,y],x,y])^2)

Create the table with the points, second partial derivative with respect to x, and Hessian:

{x,y,D[f[x,y],x,x],Hess[x,y]}/.sol//TableForm

Conclusion:
(-1,0), (1,-1), and (1,1)  are saddle points (Hessian is negative).
(-1,1) and (1,1) are local mins (since (Subscript[f, xx])>0, the function is concave up near the points).
(1,0)  is a local max (since (Subscript[f, xx])<0, the function is concave down near the points).

Here is a sketch of the function, to verify our results:

surface=Plot3D[f[x,y],{x,-2,2},{y,-2,2},AxesLabel->{"x","y","z"}]

Here I'll add small spheres at the critical points to help you see what's going on:

points=Graphics3D[{Red,
Table[Sphere[{x,y,f[x,y]},0.2]/.sol[[i]],{i,1,Length[sol]}]
}];
Show[surface,points]

Example

Find any local maximum, minimum, or saddle points of the function f(x,y)=2(x^3)+(xy^2)+5(x^2)+(y^2).

Solution

f[x_,y_]=2(x^3)+x*(y^2)+5(x^2)+(y^2)

Solve the equations where first partial derivatives are zero:

sol=Solve[{D[f[x,y],x]==0,D[f[x,y],y]==0},{x,y}]

Compute the Hessian:

D[f[x,y],x,x]
D[f[x,y],x,y]
D[f[x,y],y,y]
Hess[x_,y_]=D[f[x,y],x,x]*D[f[x,y],y,y]-((D[f[x,y],x,y])^2)

Create the table with the points, second partial derivative with respect to x, and Hessian:

{x,y,D[f[x,y],x,x],Hess[x,y]}/.sol//TableForm

Conclusion:
(-1,-2), (-1,2)  are saddle points (Hessian is negative).
(0,0) is a local min (since (Subscript[f, xx])>0, the function is concave up near the point).
(-5/3,0)  is a local max (since (Subscript[f, xx])<0, the function is concave down near the point).

Here is a sketch of the function, to verify our results:

surface=Plot3D[f[x,y],{x,-2,1},{y,-3,3},AxesLabel->{"x","y","z"}];
points=Graphics3D[{Red,
Table[Sphere[{x,y,f[x,y]},0.2]/.sol[[i]],{i,1,Length[sol]}]
}];
Show[surface,points]

Example: optimizing over a closed bounded domain

f[x_,y_]=2x^2+x+y^2-2

We want to find the absolute maximum and absolute minimum of f(x,y) on the domain (x^2)+(y^2)≤4.

In the interior:

Solve[{D[f[x,y],x]==0,D[f[x,y],y]==0},{x,y}]

pt1={-1/4,0}

On the boundary:

g[x_]=f[x,y]/.(y^2)->4-(x^2)

Solve[g'[x]==0,x]

pt2={x,(Sqrt[4-(x^2)])}/.x->-1/2
pt3={x,-(Sqrt[4-(x^2)])}/.x->-1/2

The function g(x)  also has a maximum at the endpoint where x=2 (g(x) is a parabola opening up, so the endpoints could be maximums):

g[-2]
g[2]

pt4={x,-(Sqrt[4-(x^2)])}/.x->2

Evaluate the function at the points we have found:

f[pt1[[1]],pt1[[2]]]
f[pt2[[1]],pt2[[2]]]
f[pt3[[1]],pt3[[2]]]
f[pt4[[1]],pt4[[2]]]

So f(x,y) has an absolute maximum of 8 which occurs at (2,0) on the boundary, and an absolute minimum of -17/8 which occurs at (-1/4,0) in the interior. Here's a plot to verify these results. Notice the built-in command RegionFunction that can plot surfaces subject to contraints!

Plot3D[f[x,y],{x,-3,3},{y,-3,3},RegionFunction->Function[{x,y,z},x^2+y^2≤4]]

Linear Regression

Derive the Regression Formula

data={{1,12.4442},{2,5.62473},{3,10.2803}}

points=ListPlot[data,PlotStyle->{Red,PointSize[0.02]}]

To fit a straight line to these data, we would have picture something like the following:

f[x_]=-1*x+10
fit=Plot[f[x],{x,0,30},PlotStyle->{Blue,Thick}]

Show[points,fit]

To find a way to do this without guessing, we need something to mimimize. Let's choose to minimize the square of the vertical distance between the line and all the points:

lines=Graphics[{Darker[Green],Thick,
Line[{data[[1]],{data[[1,1]],f[data[[1,1]]]}}],
Line[{data[[2]],{data[[2,1]],f[data[[2,1]]]}}],
Line[{data[[3]],{data[[3,1]],f[data[[3,1]]]}}]}
];
Show[points,fit,lines]

Our goal is to figure of the values of m and b gives a line y=m x +b that minimize the follwoing quantity, working in general with N data points ((Subscript[x, i]),(Subscript[y, i])).

Minimize: F=(∑_(i=1)^N)((((Subscript[y, i])-m (Subscript[x, i])-b))^2)

Solve ((∂F)/(∂m))=0  ((∂F)/(∂b))=0 for m and b.
 ((∂F)/(∂m))=(∑_(i=1)^N)2((Subscript[y, i])-m (Subscript[x, i])-b)(-(Subscript[x, i]))=0                  ((∂F)/(∂b))=(∑_(i=1)^N)2((Subscript[y, i])-m (Subscript[x, i])-b)(-1)=0
Rewrite these equations as
m (∑_(i=1)^N)((Subscript[x, i])^2)+b (∑_(i=1)^N)(Subscript[x, i])=(∑_(i=1)^N)(Subscript[y, i])(Subscript[x, i])
m (∑_(i=1)^N)(Subscript[x, i])+b (∑_(i=1)^N)1=(∑_(i=1)^N)(Subscript[y, i])
Solution found using Cramer's Rule or back-substitution, using  (∑_(i=1)^N)1=N,

m=((N(∑_(i=1)^N)(Subscript[y, i])(Subscript[x, i])-((∑_(i=1)^N)(Subscript[x, i]))((∑_(i=1)^N)(Subscript[y, i])))/(N(∑_(i=1)^N)((Subscript[x, i])^2)-((((∑_(i=1)^N)(Subscript[x, i])))^2)))         
 b=((((∑_(i=1)^N)((Subscript[x, i])^2))((∑_(i=1)^N)(Subscript[y, i]))-((∑_(i=1)^N)(Subscript[x, i]))((∑_(i=1)^N)(Subscript[y, i])(Subscript[x, i])))/(N(∑_(i=1)^N)((Subscript[x, i])^2)-((((∑_(i=1)^N)(Subscript[x, i])))^2)))

It's a bit tricky to use the second derivative test to prove these values of m and b produce a minimum, but the geometry of the situation should be enough to convince you they provide a minimum and not a maximum.

So, for the case above with only three points, we get:

n:=Length[data]
Sx:=(∑_(i=1)^n)data[[i,1]]
Sx2:=(∑_(i=1)^n)((data[[i,1]])^2)
Sy:=(∑_(i=1)^n)data[[i,2]]
Sxy:=(∑_(i=1)^n)data[[i,1]]*data[[i,2]]

m=((n*Sxy-Sx*Sy)/(n*Sx2-(((Sx))^2)))
b=((Sx2*Sy-Sx*Sxy)/(n*Sx2-(((Sx))^2)))

fit=Plot[m*x+b,{x,0,30},PlotStyle->{Blue,Thick}]

Show[points,fit]

This is what the built-in Mathematica command Fit does:

Fit[data, {1,x},x]

Which is the same m  and b that we found.

Example

data={{1,12.4442},{2,5.62473},{3,10.2803},{4,12.1662},{5,7.9472},{6,-3.01402},{7,4.01769},{8,10.0809},{9,10.8153},{10,3.85094},{11,-3.9047},{12,-11.658},{13,3.28664},{14,-7.78576},{15,4.18959},{16,3.19379},{17,1.81433},{18,-10.3094},{19,-8.66676},{20,-3.66591},{21,-12.9382},{22,-11.1233},{23,-17.9311},{24,-11.4085},{25,-10.3824},{26,-22.748},{27,-15.2114},{28,-11.5746},{29,-27.3296},{30,-29.734}}

points=ListPlot[data,PlotStyle->{Red,PointSize[0.02]}]

Here are the formulas we worked out by trying to minimize the sum of the distances from all the points to the line:

n:=Length[data]
Sx:=(∑_(i=1)^n)data[[i,1]]
Sx2:=(∑_(i=1)^n)((data[[i,1]])^2)
Sy:=(∑_(i=1)^n)data[[i,2]]
Sxy:=(∑_(i=1)^n)data[[i,1]]*data[[i,2]]

m=((n*Sxy-Sx*Sy)/(n*Sx2-(((Sx))^2)))
b=((Sx2*Sy-Sx*Sxy)/(n*Sx2-(((Sx))^2)))

fit=Plot[m*x+b,{x,0,30},PlotStyle->{Blue,Thick}]

Show[fit,points]

This process is in the built-in command Fit:

Fit[data, {1,x},x]