Приближение кривой методом наименьших квадратов
В данном примере продемонстрировано, как выполнить нелинейную задачу по приближению (аппроксимации) кривой методом наименьших квадратов (МНК) с помощью рабочего процесса проблемно-ориентированной оптимизации.
Мы воспользуемся функциями библиотеки JuMP.jl для формулировки оптимизационной задачи, библиотекой нелинейной оптимизации Ipopt.jl, библиотекой генерации случайных чисел Random.jl, библиотекой для работы с распределениями вероятностей Distributions.jl и библиотекой Plots.jl для визуализации результатов.
Установка библиотек
Если в вашем окружении не установлена последняя версия пакета JuMP, раскомментируйте и запустите ячейку ниже:
Pkg.add(["Ipopt", "Distributions", "JuMP"])
#Pkg.add("JuMP");
Для запуска новой версии библиотеки после завершения установки нажмите на кнопку «Личный кабинет»:
Затем нажмите на кнопку «Стоп»:
Перезапустите сессию нажатием кнопки «Старт Engee»:
Описание задачи
Задача оптимизации состоит в том, чтобы выполнить нелинейную аппроксимацию кривой методом наименьших квадратов.
Модельное уравнение зависимости для этой задачи:
где , , и — это неизвестные параметры, — отклик, а t — время.
Цель задачи состоит в том, чтобы найти значения и , которые минимизируют целевую функцию:
Подключение библиотек
Подключите библиотеку JuMP:
using JuMP;
Подключите библиотеку нелинейного решателя Ipopt:
using Ipopt;
Подключите библиотеку генерации случайных чисел Random:
using Random;
Подключите библиотеку для работы с распределениями вероятностей Distributions:
using Distributions;
Подключите библиотеку построения графиков Plots:
using Plots;
Генерация данных
Для данной задачи мы сгенерируем искусственные зашумленные данные.
Задайте произвольное зерно генератора случайных чисел:
Random.seed!(0);
Создайте переменные A_true и r_true, содержащие истинные значения параметров модели. Эти переменные используются для генерации искусственных данных ydata.
В качестве истинных значений используйте , , и :
A_true = [1, 2];
r_true = [-1, -3];
Переменные A_true и r_true являются эталоном, который алгоритм оптимизации должен попытаться воссоздать. Сравнение результатов оптимизации с истинными значениями позволяет нам оценить точность и эффективность процесса аппроксимации кривой.
Сгенерируйте 200 случайных значений tdata в диапазоне от 0 до 3 в качестве данных времени. Также создайте переменную noisedata, содержащую шум, которую мы добавим к итоговым данным ydata.
tdata = sort(3 * rand(200));
noisedata = 0.05 * randn(length(tdata));
ydata = A_true[1] * exp.(r_true[1] * tdata) + A_true[2] * exp.(r_true[2] * tdata) + noisedata;
Выведите полученные точки на график:
plot(tdata, ydata, seriestype=:scatter, color=:red, legend=:right, label="Исходные данные")
xlabel!("t")
ylabel!("Отклик")
Как видно из графика, данные содержат шум. Поэтому решение задачи, вероятно, не будет точно совпадать с истинными параметрами и .
Создание задачи оптимизации
Создайте оптимизационную задачу с помощью функции Model() и укажите название решателя в скобках:
model = Model(Ipopt.Optimizer)
Создайте переменные и , содержащие по два значения, – , , и :
@variable(model, A[1:2]);
@variable(model, r[1:2]);
Создайте нелинейную целевую функцию, минимизирующую сумму квадратов:
@NLobjective(model, Min, sum((A[1]*exp(r[1]*t) + A[2]*exp(r[2]*t) - y)^2 for (t,y) in zip(tdata, ydata)));
Установка соответствующих начальных значений – важный шаг для эффективного и точного решения задач оптимизации, существенно влияющий на скорость, качество и надежность процесса оптимизации.
Укажите начальные значения , , и .
set_start_value.(A, [1/2, 3/2]);
set_start_value.(r, [-1/2, -3/2]);
Выведите начальные значения и :
println("Начальные значения A: ", start_value.(A))
println("Начальные значения r: ", start_value.(r))
Решение задачи
Решите оптимизационную задачу:
optimize!(model)
Сохраните значения и в переменных:
A_sol = value.(A);
r_sol = value.(r);
Выведите результаты оптимизации:
println("Оптимизированные значения A: ", A_sol)
println("Оптимизированные значения r: ", r_sol)
Визуализация и анализ результатов
Рассчитайте результат аппроксимации кривой:
fitted_response = A_sol[1] * exp.(r_sol[1] * tdata) + A_sol[2] * exp.(r_sol[2] * tdata);
Выведите результаты на график:
plot!(tdata, fitted_response, color=:blue, label="Приближённая кривая")
title!("Приближённый отклик")
Вы можете сравнить полученные результаты с изначально заданными значения , , и . А функция round позволит округлить избыточно точные значения.
Сравните результаты с изначальными значениями:
percent_diffs = ((A_sol .- A_true) ./ A_true) .* 100
avg_percent_diff = mean(percent_diffs)
println("Процентная разница: ", round.(percent_diffs, digits=2))
println("Средняя процентная разница: ", round.(avg_percent_diff, digits=2))
Сравните результаты с изначальными значениями:
percent_diffs = ((r_sol .- r_true) ./ r_true) .* 100
avg_percent_diff = mean(percent_diffs)
println("Процентная разница: ", round.(percent_diffs, digits=2))
println("Средняя процентная разница: ", round.(avg_percent_diff, digits=2))
Таким образом, результаты в среднем отличаются от изначальных значений на , а результаты на .
Заключение
В данном примере мы решили нелинейную задачу по приближению (аппроксимации) кривой с помощью метода наименьших квадратов (МНК), используя проблемно-ориентированный подход. Мы также визуализировали и произвели сравнение результатов с истинными параметрами модели.