Методы решения инженерных задач
Регрессионный анализ
Регрессионный анализ применяется для описания зависимости между различными величинами в виде функций. Методы регрессионного анализа применяются в химической технологии для описания экспериментальных данных (очень часто используется для описания физико-химических свойств веществ, участвующих в процессе). Также применяется для подстройки математических моделей по экспериментальным данным.
Полином
Полином является одной из наиболее часто используемых функций для описания данных. Данная функция простая, от нее легко определяется производная и интеграл и есть строгий алгоритм для определения параметров модели.
Полином n-ой степени:
Расписывая суммирование полином первой степени будет являться линейной функцией:
Второй степени — параболой:
и т. д.
Для проведения вычисления с полиномами в julia используется пакет Polynomials.jl . Параметры полинома можно задавать вручную, однако проще найти функцию по экспериментальным данным.
using Polynomials
#исходные данные
X = [1.0, 2.0, 3.0, 4.0, 5.0, 6.0]
Y = [0.3, 0.9, 1.2, 2.3, 2.8, 2.5]
Для определения параметров полинома по экспериментальным данным используется функция fit. Если использовать ее с двумя аргументами, то параметр полинома будет на единицу меньше количества экспериментальных данных. В данном случае количество неизвестных будет равно количеству параметров.
pf0 = fit(X,Y)
using Plots #подключим пакет для построения графиков
xx = 1:0.1:6 #данные в диапазоне исследуемых значений с маленьким шагом
scatter(X,Y)
plot!(xx, pf0.(xx)) #! - рисуем на этом же графике, . в названии функции - применение функции поэлементно
С увеличением степени полинома увеличивается количество экстремумов. Однако физические величины обычно не имеют множества экстремумов, поэтому нужно внимательно к выбору порядка полинома.
Степень полинома можно задать третьим аргументом функции fit.
pf1 = fit(X,Y,1) #линейная функция (полином 1 степени)
pf2 = fit(X,Y,2) #парабола (полином 2 степени)
pf3 = fit(X,Y,3) #кубическая парабола (полином 3 степени)
scatter(X,Y, label = "")
plot!(xx, pf1.(xx), label = "линейная")
plot!(xx, pf2.(xx), label = "парабола")
plot!(xx, pf3.(xx), label = "кубическая парабола")
Нелинейные по парметрам функции.
Все виды используемых для описания каких-либо экспериментальных данных можно разделить на 2 больших класса:
- Линейные по параметрам функции можно предаставить в виде , тоесть функцию модно представить как сумму произведения параметра на какую-либо функцитю которая зависит только от аргумента x
- Нелинейные по парамтерам функии не получается предавтить таком виде и параметр входит в функцию . Например в функции параметр b нельзя вытащить за функцию экспененты. Однако можно преобразовать саму функцию и находить не само значение y, а его логарифм.
С точки зрения регрессионного анализа данная классификация важна, потому что используются разные алгоритмы для опредения параметров.
Рассмотрим алгоритм для определения параметров полинома 2 степени.
Функция, описывающая данные должна быть как можно ближе к ним. Обычно находят значения параметров минимизируя квкадрат отклонения экспериментальных значений и значений по найденных по функции .
В точке эстремума производные по параметрам будут равны нулю
Подставляя выражение линейной по парамтерам функции получим:
Продифференцировав получим:
Однако в случае нелинейных по парамтерам функций получить в общем виде значение производной нельзя. Поэтому применяются другие алгоритмы для поиска парамтеров.
Перегруппировав члены выражения:
Решение этой системы уравнений лучше проводить в матричном виде:
Для этого вводят понятие матрицы входных переменных:
Реализуем ее для наших данных и полинома второй степени. значит , ,
Φ = [1 X[1] X[1]^2;
1 X[2] X[2]^2;
1 X[3] X[3]^2;
1 X[4] X[4]^2;
1 X[5] X[5]^2;
1 X[6] X[6]^2]
Информационная матрица:
:
I = transpose(Φ) * Φ
Вектор где - вектор экспериментальных данных:
Yₘ = [Y[1]; Y[2]; Y[3]; Y[4]; Y[5]; Y[6]]
Матрица B:
B = transpose(Φ) * Yₘ
Матрица параметров:
A= I^-1 * B
Сравнивая с результататми, полученными из пакета Polynomials.jl. Параметры идентичные:
fit(X,Y,2)
Нелинейные по параметрам функции
Для нелинейных по парамтерам функций применяются алгоритмы, отличные от рассмотернного выше. Применяются разные виды алгоритмов, например алгоритм Левенберга — Марквардта. Основным отличием от описанного выше алгоритма является обязательное использование начального приближения - значений параметров с корторых стартует алгоритм.
Недостатком нелинейных по параметрам функций ялвяется зависимость решения от значений начального приближения. Алгоритм может найти параметры которые будут описывать данные хорошо, а может вернуть параметры которые дают большое отклонение или выдать ошибку.
Для опредделения парамтеров будем использовать пакет LsqFit.
using LsqFit
Функция, которой будут описаны экспериментальные данные, должна также зависеть от параметров (p в аргументах функции). Т.к. параметров у функции может быть несколько, то p является массивом, содержащим все паратмеры.
Примедем пример, для описания представленных выше данных функций:
, где a и b - параметры
myfunc(x, p) = @. p[1] * exp(p[2] * x)
Для определения параметров используется функция curve_fit. Аргументы функции:
- myfunc - название функции, которой описываем данные
- x, y - данные которые описываются
- [0.1, 1000.0] - начальное приближение
В случае задания начального приближения в виде [0.1, 1000.0], для параметра a используется значение 0.1, для b - 1000. Ошибка возникает в связи с тем, что на первом этапе идет расчет значений по функции с начальными приближениями. Величина дает очень большое число за пределами вычислитеной возможности типа Float64.
curve_fit(myfunc, X, Y, [1.1, 1000.0])
С начальными приближениями [1.1, 10.0] значения параметров нашлись. Однако все равно найденные значения дают плохое описание, что видно на рисунке ниже. В результате функция выдает стуктуру (несколько значений), в которой содержится несколько переменных под своими названиями.
a = curve_fit(myfunc, X, Y, [1.1, 10])
Например посмотреть сошелся ли алгоритм можно по элементу стрктуры converged. Оно дает логическое да (true) или нет (false)
a.converged
pfail = curve_fit(myfunc, X, Y, [1.1, 10.0]).param
p = curve_fit(myfunc, X, Y, [1.1, 0.1]).param
mf2(x) = myfunc(x, p)
scatter(X,Y, label="")
plot!(xx, mf2(xx), label = "\$ y(x) = a_1 e^{a_2 x}\$")
plot!(xx, myfunc(xx, pfail), label = "неправильные параметры")
Сплайны
Сплайн - кусочная функция. Каджый дапазон между точками описывается функцией одного вида, но с разным набором параметров. Спалайны бывают разного вида. Например линейный сплайн соединяет экспериментальные точки прямыми. Параболический и кубический сплайны соединяются соответственно полиномом первой и второй степени.
Дополнительным условием для параболического и кубического сплайна является равенство производной в точках соприкосновения двух функций. Таким образом сплайн будет функцией непрерывной (соединяется строго в точках которые мы задали) и гладкой (производная от функции не имеет разрывов).
using Dierckx
Степень сплайна задается параметром функция Spline1D создает одномерный (зависящий только от одной переменной) сплайн по указанным данным
spl1 = Spline1D(X,Y, k=1)
spl2 = Spline1D(X,Y, k=2)
spl3 = Spline1D(X,Y, k=3)
scatter(X,Y, label="")
plot!(xx,spl1(xx), label="линейная интерполяция")
plot!(xx,spl2(xx), label="параболический сплайн")
plot!(xx,spl3(xx), label="кубический сплайн")