W praktyce dane otrzymywane w wielu eksperymentach wcale ne muszą mieć charakteru liniowego.
Załóżmy, że poszukiwaną funkcją jest wielomian:
\(\displaystyle{ p_{m}(x) = \sum_{k=0}^{m} a_{k}x^{k} \ \ (1) }\)
stopnia \(\displaystyle{ m \leq n-1. }\)
Stosując zasadę najmniejszych kwadaratów - mamy znaleźć współczynniki \(\displaystyle{ a_{0}, a_{1}, ...,a_{m}, }\) dla których forma kwadratowa:
\(\displaystyle{ E(a_{0},a_{1}, ..., a_{m}) = \sum_{i=1}^{n}[ p_{m}(x_{i}) -y_{i}]^2 = \sum_{i=1}^{n} \left [ \sum_{k=0}^{m} (a_{k}x^{k}_{i})-y_{i} \right ]^2 \ \ (2) }\)
osiąga minimum.
\(\displaystyle{ E }\) osiąga minimum (warunek wystarczający pomijamy), jeśli
\(\displaystyle{ E'_{|a_{j}}(a_{0},a_{1},...,a_{m})=0, \ \ j=0,1,..., m.}\)
\(\displaystyle{ E'_{|a_{0}} = \sum_{i=1}^{n}2 \left [ \sum_{k=0}^{m} (a_{k}x^{k}_{i})-y_{i} \right ] = 0 }\)
\(\displaystyle{ E'_{|a_{1}} = \sum_{i=1}^{n}2 \left [ \sum_{k=0}^{m} (a_{k}x^{k}_{i})-y_{i} \right ](x_{i}) = 0 }\)
\(\displaystyle{ ......................................}\)
\(\displaystyle{ E'_{|a_{m}} = \sum_{i=1}^{n}2 \left [ \sum_{k=0}^{m} (a_{k}x^{k}_{i})-y_{i} \right ](x^{m}_{i}) = 0 }\)
Przekształcając równoważnie ten układ równań, otrzymujemy \(\displaystyle{ (m+1) }\) rownań normalnych z \(\displaystyle{ (m+1) }\) niewiadomymi \(\displaystyle{ a_{0},
a_{1}, ..., a_{m}:}\)
\(\displaystyle{ a_{o}\cdot n + a_{1} \sum_{i=1}^{n}x_{i} + ...+ a_{m} \sum_{i=1}^{n}x^{m}_{i} = \sum_{i=1}^{n}y_{i} }\)
\(\displaystyle{ a_{o}\sum_{i=1}^{n} x_{i} + a_{1} \sum_{i=1}^{n}x^2_{i} + ...+ a_{m} \sum_{i=1}^{n}x^{m+1}_{i} = \sum_{i=1}^{n}y_{i}x_{i} }\)
\(\displaystyle{ a_{o}\sum_{i=1}^{n} x^2_{i} + a_{1} \sum_{i=1}^{n}x^3_{i} + ...+ a_{m} \sum_{i=1}^{n}x^{m+2}_{i} = \sum_{i=1}^{n}y_{i}x^2_{i} }\)
\(\displaystyle{ ...............................................}\)
\(\displaystyle{ [ a_{o}\sum_{i=1}^{n} x^{m}_{i} + a_{1} \sum_{i=1}^{n}x^{m+1}_{i} + ...+ a_{m} \sum_{i=1}^{n}x^{2m}_{i} = \sum_{i=1}^{n}y_{i}x^{m}_{i} }\)
Układ równań normalnych o macierzy Hilberta ma dokładnie jedno rozwiązanie, ale gdy wartości \(\displaystyle{ x_{i}, \ \ i=1,2,..,n }\) są
równomiernie rozłożone to macierz ta, jest macierzą źle uwarunkowaną. Wtedy wielomian aproksymujący jest obarczony błędami
zaokrągleń.
Patrz na przykład Janina i MIchał Jankowscy. Przegląd metod i algorytmów numerycznych. Część 1. WNT Warszawa 1981.
Program w MATLAB
Kod: Zaznacz cały
function mnk(x,y,m)
% Wielomian aproksymujący w sensie najmniejszych kwadratów
n=length(x);
n=length(y);
for k=1:2*m+1
c(k)=sum(x.^(k-1));
end
% Współrzędne wektora b układu równań normalnych
for k=1:m+1
b(k)=sum(y.*x.^(k-1));
end
% Elementy macierzy A układu równań normalnych
for i=1:m+1
for j=1:m+1
A(i,j)=c(j+i-1);
end
end
fprintf('\n')
disp('Macierz rozszerzona ukladu ,[A b]=')
fprintf('\n')
disp([A b'])
z=A\b';
disp('Wspolczynnikami a0,a1,...,an wielomianu aproksymującego są')
fprintf('\n')
disp(z')
% Wartości wielomianu w węzłach xi, i=1,2,...,n
disp('_____________________________________________________')
disp(' xi yi p(xi) |yi-p(xi)| ')
disp('______________________________________________________')
for i=1:n
s=z(1);
for k=2:m+1
s=s+z(k)*x(i)^(k-1);
end
p(i)=s;
err(i)=abs(y(i)-p(i));
fprintf('%6.2f %6.2f %12.6f %12.6f\n',x(i),y(i),p(i),err(i))
end
err=sum(err.*err);
fprintf('\n E(a0,a1,...,an)= %12.6g\n',sum(err))
[x' y' p']
Kod: Zaznacz cały
>> x=[0:0.5:2.5];
>> y=[0 0.20 0.27 0.30 0.32 0.33];
>> mnk(x,y,3)
Macierz rozszerzona ukladu ,[A b]=
6.0000 7.5000 13.7500 28.1250 1.4200
7.5000 13.7500 28.1250 61.1875 2.2850
13.7500 28.1250 61.1875 138.2813 4.3375
28.1250 61.1875 138.2813 320.5469 9.0237
Wspolczynnikami a0,a1,...,an wielomianu aproksymującego są
0.0033 0.4970 -0.2738 0.0511
_____________________________________________________
xi yi p(xi) |yi-p(xi)|
______________________________________________________
0.00 0.00 0.003333 0.003333
0.50 0.20 0.189762 0.010238
1.00 0.27 0.277619 0.007619
1.50 0.30 0.305238 0.005238
2.00 0.32 0.310952 0.009048
2.50 0.33 0.333095 0.003095
E(a0,a1,...,an)= 0.000292857
ans =
0 0 0.0033
0.5000 0.2000 0.1898
1.0000 0.2700 0.2776
1.5000 0.3000 0.3052
2.0000 0.3200 0.3110
2.5000 0.3300 0.3331
Kod: Zaznacz cały
>> x=[-3.0:0.4:3.3];
>> y=sin(x);
>> mnk(x,y,4)
Macierz rozszerzona ukladu ,[A b]=
1.0e+04 *
0.0016 0.0000 0.0054 -0.0000 0.0331 -0.0000
0.0000 0.0054 -0.0000 0.0331 0 0.0016
0.0054 -0.0000 0.0331 0 0.2388 0
-0.0000 0.0331 0 0.2388 0 0.0062
0.0331 0 0.2388 0 1.8648 -0.0000
Wspolczynnikami a0,a1,...,an wielomianu aproksymującego są
-0.0000 0.8552 0.0000 -0.0928 -0.0000
_____________________________________________________
xi yi p(xi) |yi-p(xi)|
______________________________________________________
-3.00 -0.14 -0.060136 0.080984
-2.60 -0.52 -0.592535 0.077033
-2.20 -0.81 -0.893326 0.084830
-1.80 -0.97 -0.998142 0.024295
-1.40 -0.99 -0.942615 0.042835
-1.00 -0.84 -0.762376 0.079095
-0.60 -0.56 -0.493058 0.071585
-0.20 -0.20 -0.170291 0.028378
0.20 0.20 0.170291 0.028378
0.60 0.56 0.493058 0.071585
1.00 0.84 0.762376 0.079095
1.40 0.99 0.942615 0.042835
1.80 0.97 0.998142 0.024295
2.20 0.81 0.893326 0.084830
2.60 0.52 0.592535 0.077033
3.00 0.14 0.060136 0.080984
E(a0,a1,...,an)= 0.0685988
ans =
-3.0000 -0.1411 -0.0601
-2.6000 -0.5155 -0.5925
-2.2000 -0.8085 -0.8933
-1.8000 -0.9738 -0.9981
-1.4000 -0.9854 -0.9426
-1.0000 -0.8415 -0.7624
-0.6000 -0.5646 -0.4931
-0.2000 -0.1987 -0.1703
0.2000 0.1987 0.1703
0.6000 0.5646 0.4931
1.0000 0.8415 0.7624
1.4000 0.9854 0.9426
1.8000 0.9738 0.9981
2.2000 0.8085 0.8933