Аппроксимация функций. Метод наименьших квадратов
В инженерной деятельности часто возникает необходимость описать в виде функциональной зависимости связь между величинами, заданными таблично или в виде набора точек с координатами (xi, yi), i = 0, 1, 2,... , n, где n - общее количество точек. Как правило, эти табличные данные получены экспериментально и имеют погрешности (рис. 2.5)
Рис.4.2
При аппроксимации желательно получить относительно простую функциональную зависимость (например, многочлен), которая позволила бы "сгладить" экспериментальные погрешности, вычислять значения функции в точках, не содержащихся в исходной таблице.
Эта функциональная зависимость должна с достаточной точностью соответствовать исходной табличной зависимости. В качестве критерия точности чаще всего используют критерий наименьших квадратов, т.е. определяют такую функциональную зависимость f(x), при которой
S =, (4.12)
обращается в минимум.
Погрешность приближения оценивается величиной среднеквадратического уклонения
= . (4.13)
В качестве функциональной зависимости рассмотрим многочлен
Pm(x)=a0 + a1x + a2x2+...+amxm. (4.14)
Формула (4.12) примет вид
S =
Условия минимума S можно записать, приравнивая нулю частные производные S по всем переменным a0, a1, a2, … , am. Получим систему уравнений
= -= 0, или
= 0, k = 0, 1, … , m. (4.15)
Систему уравнений (4.15) перепишем в следующем виде:
a0+ a1+ … +am= , k = 0, 1, … , m (4.16)
Введем обозначения:
ck = , bk = .
Система (4.16) может быть записана так:
a0ck + a1ck+1 + … + ck+mam = bk, k = 0, 1, … , m. (4.17)
Перепишем систему (4.17) в развернутом виде:
c0a0 + c1a1 + c2a2… + cmam = b0
c1a0 + c2a1 + c3a2… + cm+1am = b1
(4.18)
cma0 + cm+1a1 + cm+2a2… + c2mam = bm
Матричная запись системы (4.18) имеет следующий вид:
Ca = b. (4.19)
Для определения коэффициентов ak, k = 0, 1, … , m, и, следовательно, искомого многочлена (4.14) необходимо вычислить суммы ck, bk и решить систему уравнений (4.18). Матрица C системы (4.19) называется матрицей Грама и является симметричной и положительно определенной. Эти полезные свойства используются при решении.
Погрешность приближения в соответствии с формулой (4.13) составит
= . (4.20)
Рассмотрим частные случаи m =1 и m = 2.
1. Линейная аппроксимация (m = 1).
P1(x) = a0 + a1x.
c0 = = n + 1; c1 = = ; c2 = ; (4.21)
b0 = = ; b1 = = . (4.22)
c0 c1 n+1
C = = ,
c1 c2
b = (b0, b1)T = (,)T.
Решение системы уравнений Ca = b найдем по правилу Крамера:
a0 = , a1 = ,
где C - определитель матрицы C, аCi - определитель матрицы Ci, полученной из матрицы C заменой i-го столбца столбцом свободных членов b, i = 1, 2.
Таким образом,
a0 = , a1 = . (4.23)
Алгоритм 4.1 (Алгоритм метода наименьших квадратов. Линейная аппроксимация).
Шаг 1. Ввести исходные данные: xi, yi, i=0, 1, 2, ... , n.
Шаг 2. Вычислить коэффициенты c0, c1, b0, b1 по формулам (4.21), (4.22).
Шаг 3. Вычислить a0, a1 по формулам (4.23).
Шаг 4. Вычислить величину погрешности
1 = . (4.24)
Шаг 5. Вывести на экран результаты: аппроксимирующую линейную функцию P1(x) = a0 + a1x и величину погрешности 1.
2. Квадратичная аппроксимация (m = 2).
P2(x) = a0 + a1x + a2x2.
c0 == n+1; c1 ==; c2 =; c3 =; c4 =. (4.25)
b0 ==; b1 ==; b2 = . (4.26)
c0 c1 c2
C = c1 c2 c3 .
c2 c3 c4
b = (b0, b1, b2)T .
Решение системы уравнений Ca = b найдем по правилу Крамера:
ai = , i = 0, 1,
где C - определитель матрицы C, аCi - определитель матрицы Ci, полученной из матрицы C заменой i-го столбца столбцом свободных членов b.
C = c0c2c4 + 2c1c2c3 - c - сc4 - cc0. (4.27)
b0 c1 c2
C1 = b1 c2 c3 = b0c2c4 + b2c1c3 + b1c2c3 - b2c- b1c1c4 - b0c. (4.28)
b2 c3 c4
c0 b0 c2
C2 = c1 b1 c3 = b1c0c4 + b0c2c3 + b2c1c2 - b1c- b0c1c4 - b2c0c3. (4.29)
c2 b2 c4
c0 c1 b0
C3 = c1 c2 b1 = b2c0c2 + b1c1c2 + b0c1c3 - b0c- b2c - b1c0c3. (4.30)
c2 c3 b2
a0 = , a1 = , a2 = . (4.31)
Алгоритм 4.2 (Алгоритм метода наименьших квадратов. Квадратичная аппроксимация).
Шаг 1. Ввести исходные данные: xi, yi, i=0, 1, 2, ... , n.
Шаг 2. Вычислить коэффициенты c0, c1, c2, c3, c4, b0, b1, b2, по формулам (4.25), (4.26).
Шаг 3. Вычислить C, C1, C2, C3 по формулам (4.27) - (4.30).
Шаг 4. Вычислить a0, a1, a2 по формулам (4.31).
Шаг 5. Вычислить величину погрешности
2 = . (4.32)
Шаг 5. Вывести на экран результаты : аппроксимирующую квадратичную функцию P2(x) = a0 + a1x + a2x2 и величину погрешности 2.
Пример 4.6.
Построим по методу наименьших квадратов многочлены первой и второй степени и оценим степень приближения. Значения yi в точках xi , i =0, 1, 2, 3, 4 приведены в таблице 2.3.
Таблица 4.1
i | ||||||
xi | ||||||
yi | -1 | |||||
Вычислим коэффициенты c0, c1, c2, c3, c4, b0, b1, b2, по формулам (4.25), (4.26):
c0 = 5; c1 = 15; c2 = 55; c3 = 225; c4 = 979;
b0 = 12; b1 = 53; b2 = 235.
1. Линейная аппроксимация (m =1).
Система уравнений для определения коэффициентов a0 и a1 многочлена первой степени P2(x) = a0 + a1x + a2x2 имеет вид
5a0 + 15a1 = 12
15a0 + 55a1 = 53
По формулам (4.23) найдем коэффициенты a0 и a1:
a0 = -2.7, a1 = 1.7.
P1(x) = a0 + a1x = -2.7 + 1.7x.
2. Квадратичная аппроксимация (m =2).
Система уравнений для определения коэффициентов a0, a1 и a2 многочлена второй степени P2(x) = a0 + a1x + a2x2 имеет вид
5a0 + 15a1 + 55a2 = 12
15a0 + 55a1 + 225a2 = 53
55a0 + 225a1 + 979a2 = 235
По формулам (4.31) найдем коэффициенты a0, a1 и a2:
a0 -2.20, a1 1.27, a2 0.07.
P2(x) = a0 + a1x + a2x2 = -2.20 + 1.27x + 0.07x2.
Сравним значения, рассчитанные для функциональной зависимости, с исходными данными. Результаты приведены в табл.2.4.
Таблица 4.2
i | ||||||
xi | ||||||
yi | -1 | |||||
P1(xi) | -1 | 0.7 | 2.4 | 4.1 | 5.8 | |
P2(xi) | -1 | 0.62 | 2.24 | 6.9 | ||
Погрешность приближения в соответствии с формулами (4.24) и (4.32) составит
1 = = 0.245.
2 = = 0.084.