Аппроксимация функции полиномом
Введение
В вычислительной математике и обработке экспериментальных данных одной из ключевых задач является восстановление аналитической зависимости по дискретному набору наблюдений. Часто истинный закон, описывающий физический процесс, неизвестен или слишком сложен для прямого анализа, а исходные данные содержат неизбежные погрешности измерений. В таких случаях возникает необходимость найти относительно простую функцию, которая наилучшим образом в смысле минимума среднеквадратичного отклонения описывала бы имеющиеся точки. Эту задачу решает классический метод наименьших квадратов.
Импорт библиотек
В данном примере нам понадобится библиотека LinearAlgebra для операций с матрицами и Printfдля табличного вывода результатов.
using LinearAlgebra, Printf
Определение функции
Определим исследуемую функцию. В данном примере мы исследуем sinc-функцию (кардинальный синус).
f(x) = x == 0 ? 1.0 : sin(x)/x
Таблица значений функции
Для того чтобы подобрать полином, который наилучшим образом описывает поведение исследуемой функции, сначала необходимо задать конкретные значения этой функции в конкретных точках. Определим шаг и интервал, и построим таблицу значений функции.
a = 0
b = 16
шаг = 1.0
X = collect(a:шаг:b)
Y = [f(x) for x in X]
println("Таблица исходных значений функции:\n")
for i in eachindex(X)
@printf("X = %.4f, Y = %.6f\n", X[i], Y[i])
end
Метод наименьших квадратов
Реализуем метод наименьших квадратов для подбора коэффициентов полинома заданной степени. Для этого необходимо решить систему нормальных уравнений:
Где:
В результате получим массив коэффициентов искомого полинома, который наилучшим образом (в смысле среднеквадратичного отклонения) приближает исходные данные.
function МНК(x, y, степень)
M = [x[i]^j for i in 1:length(x), j in 0:степень]
return (M' * M) \ (M' * y)
end
Определение полиномов
В общем виде полином степени для аппроксимации функции методом наименьших квадратов выглядит следующим образом:
Определим и сравним полиномы разной степени для аппроксимации исследуемой функции. В результате получим таблицу коэффициентов для полиномов, которые в дальнейшем будем использовать для вычисления значений в новых точках и визуального сравнения с исходной функцией. Это позволяет оценить, как увеличение степени полинома влияет на качество аппроксимации.
степени = [1, 3, 5, 7, 9]
П = Dict()
println("\nКоэффициенты полиномов наилучшего среднеквадратичного приближения:")
for d in степени
К = МНК(X, Y, d)
П[d] = К
println("\nСтепень = $d")
for i in 1:length(К)
@printf("a%d = %.6f\n", i-1, К[i])
end
end
Сравним исходную функцию с полученными полиномами для детального анализа качества аппроксимации. В результате мы увидим, как каждый полином приближает исходную функцию на всём интервале
ш = 0.1
X = collect(a:ш:b)
function Значение(К, x)
sum(К[i] * x^(i-1) for i in 1:length(К))
end
println("\nТаблица значений функции и полиномов:")
for x in X
@printf("\nx = %.4f | f(x)=%.6f", x, f(x))
for d in степени
значение = Значение(П[d], x)
@printf(" | P%d=%.6f", d, значение)
end
end
Аппроксимация полиномом 1-го порядка
Построим графики для визуального сравнения функции и её аппроксимации полиномом первого порядка.
d = 1
график = plot(X, [f(x) for x in X], label="f(x)=sin(x)/x", lw=2)
plot!(X, [Значение(П[d], x) for x in X], label="Полином $d-го порядка",
xlabel="X", ylabel="Y", lw=2, legend=:right)
display(график)
.png)
Полином первого порядка это линейная функция. Очевидно, совпадения с исследуемой функцией нет.
Аппроксимация полиномом 3-го порядка
Построим графики для визуального сравнения функции и её аппроксимации полиномом третьего порядка.
d = 3
график = plot(X, [f(x) for x in X], label="f(x)=sin(x)/x", lw=2)
plot!(X, [Значение(П[d], x) for x in X], label="Полином $d-го порядка",
xlabel="X", ylabel="Y", lw=2, legend=:right)
display(график)
.png)
Полином 3-го порядка это кубическая функция, совпадение становится ближе.
Аппроксимация полиномом 5-го порядка
Построим графики для визуального сравнения функции и её аппроксимации полиномом пятого порядка.
d = 5
график = plot(X, [f(x) for x in X], label="f(x)=sin(x)/x", lw=2)
plot!(X, [Значение(П[d], x) for x in X], label="Полином $d-го порядка",
xlabel="X", ylabel="Y", lw=2, legend=:right)
display(график)
.png)
При аппроксимации полиномом 5-го порядка совпадение становится ещё ближе.
Аппроксимация полиномом 7-го порядка
Построим графики для визуального сравнения функции и её аппроксимации полиномом седьмого порядка.
d = 7
график = plot(X, [f(x) for x in X], label="f(x)=sin(x)/x", lw=2)
plot!(X, [Значение(П[d], x) for x in X], label="Полином $d-го порядка",
xlabel="X", ylabel="Y", lw=2, legend=:right)
display(график)
.png)
При аппроксимации полиномом 7-го порядка совпадение становится ещё ближе, чем при аппроксимации полиномом 5-го порядка.
Аппроксимация полиномом 9-го порядка
Построим графики для визуального сравнения функции и её аппроксимации полиномом девятого порядка.
d = 9
график = plot(X, [f(x) for x in X], label="f(x)=sin(x)/x", lw=2)
plot!(X, [Значение(П[d], x) for x in X], label="Полином $d-го порядка",
xlabel="X", ylabel="Y", lw=2, legend=:right)
display(график)
.png)
При аппроксимации полиномом 9-го порядка, исследуемая функция практически совпадает с полиномом. Таким образом можно сделать вывод о том, что чем выше порядок полинома, тем выше его совпадение с исследуемой функцией.
Заключение
В данном примере была численно реализована и исследована классическая задача полиномиальной аппроксимации методом наименьших квадратов, который является фундаментальным инструментом в самых разных областях деятельности, где требуется восстанавливать зависимости по зашумленным или дискретным данным:
-
Физика и инженерия — обработка результатов экспериментов (например, аппроксимация вольт-амперных характеристик, калибровка датчиков, сглаживание сигналов).
-
Экономика и финансы — построение трендов временных рядов, прогнозирование цен, анализ рыночных индикаторов.
-
Биология и медицина — моделирование дозовых зависимостей, анализ роста популяций, обработка данных клинических испытаний.
-
Геофизика и астрономия — выделение полезного сигнала на фоне шумов, интерполяция спутниковых данных.
-
Машинное обучение — обобщение данных (линейная и полиномиальная регрессия) как основа многих более сложных алгоритмов.
Таким образом, рассмотренный в примере подход не ограничивается чисто учебной задачей, а составляет базу для множества прикладных численных методов. Выбор Engee позволяет в дальнейшем масштабировать вычисления на большие объёмы данных без потери производительности, что особенно важно при решении реальных инженерных и научных задач.