优化设计
在这个例子中,我们将演示几个简单优化问题的解决方案:使用包查找双参数函数的最小值 Optim,梯度下降的显式实现和线性规划问题中具有约束的最优解的搜索。
让我们做一个加载库的准备单元。
In [ ]:
Pkg.add(["GLPK", "LinearAlgebra", "StatsBase", "JuMP", "Optim"])
In [ ]:
# Загрузим необходимые библиотеки
using Plots
plotlyjs()
Out[0]:
几个变量函数的优化
http://julianlsolvers.github.io/Optim.jl/v0.9.3/user/minimization/
求函数的最小值 使用图书馆 Optim.
指定起点 .
In [ ]:
using Optim
f_optim(x) = (1.0 - x[1])^2 + 100.0 * (x[2] - x[1]^2)^2
results_auto = Optim.optimize(f_optim, [0.0, 0.0])
print(results_auto)
实现梯度下降
https://nextjournal.com/DeepLearningNotes/Ch06StochasticGradientDescent
本练习的目的是通过梯度下降实现线性回归过程。 首先,运行准备单元。
In [ ]:
using Plots, LinearAlgebra, StatsBase
此任务的输入数据将在以下单元格中创建。
In [ ]:
regX = rand(100)
regY = 50 .+ 100 * regX + 2 * randn(100);
scatter( regX, regY, fmt = :png, legend=:bottomright, label="data" )
Out[0]:
这个问题的确切解决方案可以表示如下。
In [ ]:
ϕ = hcat(regX.^0, regX.^1);
δ = 0.01;
θ = inv(ϕ'ϕ + δ^2*I)*ϕ'regY;
让我们定义一个二次损失函数并计算梯度。 在这个函数之后,我们可以绘制优化函数的图形,但我们会在练习的最后绘制它。
In [ ]:
function loss(X,y,θ)
cost = (y .- X*θ)'*(y .- X*θ);
grad = - 2 * X'y + 2*X'X*θ;
return (cost, grad)
end
Out[0]:
In [ ]:
function sgd(loss, X, y, θ_start, η, n, num_iters)
#=
Перечень аргументов:
loss -- функция, которую мы оптимизируем,
она принимает датасет (X, y) и
текущий набор параметров θ;
функция должна возвращать
значение потерь и градиент при наборе параметров θ
theta0 -- исходная оценка вектора параметров θ
η -- шаг обучения (learning rate)
n -- размер батча
num_iters -- макс.количество выполнения SGD
Возвращает:
theta -- вектор параметров по завершению SGD
path -- оценки оптимального θ на каждом шаге обучения
=#
data = [(X[i,:],y[i]) for i in 1:length(y)]
θ = θ_start
path = θ
for iter in 1:num_iters
s = StatsBase.sample(data, n)
Xs = hcat(first.(s)...)'
ys = last.(s)
cost, grad = loss(Xs,ys,θ)
θ = θ .- ( ( η * num_iters/(num_iters+iter) * grad) / n )
path = hcat(path,θ)
end
return (θ,path)
end
Out[0]:
具有优化进度的图形输出:
In [ ]:
# Поверхность функции потерь
xr = range(40, stop=60, length=10);
yr = range(90, stop=110, length=10);
X = hcat(ones(length(regX)), regX);
f(x,y) = loss(X,regY,[x y]')[1][1];
plot(xr, yr, f, st=:wireframe, fmt = :png, label="Функция потерь", color="blue")
# Пути оптимизации
X = hcat(ones(length(regX)), regX);
y = regY;
(t,p) = sgd(loss, X, y, [40 90]', 0.1, 1, 500);
scatter!( p[1,:], p[2,:], f.(p[1,:],p[2,:]),
label="Путь оптимизации", ms=3, markerstrokecolor = "white", camera=(40,30) )
Out[0]:
让我们显示优化结果:
In [ ]:
println("Координаты найденной точки, в которой функция принимает минимальное значение:")
print( p[:,end] )
用约束条件找到最优解
假设我们面临这样的优化问题:
可接受值的范围如图所示(x和y是连续的,自变量):
我们给出了一个优化目标(红色箭头)和一个方向(最大化),绿点是最优解,其坐标为x=3.75和y=1.25。 目标函数的最优值为23.75。
让我们使用JuMP库Julia进行优化来解决这个问题。
In [ ]:
using JuMP, GLPK
# Зададим модель
m = Model(); # m будет значить "модель"
# Выберем оптимизатор
set_optimizer(m, GLPK.Optimizer); # выберем решатель GLPK (GNU Linear Programming Kit)
# Определим
@variable(m, x>=0);
@variable(m, y>=0);
# Определим ограничения
@constraint(m, x+y<=5);
@constraint(m, 10x+6y<=45);
# Определим целевую функцию
@objective(m, Max, 5x+4y);
# Запустим решатель
optimize!(m);
# Вывод
println(objective_value(m)) # Значение целевой функции в оптимальной точке
println("x = ", value.(x), "\n","y = ",value.(y)) # Координаты оптимальных x и y