Akademia Górniczo-Hutnicza
Katedra Robotyki i Mechatroniki
Identyfikacja i analiza sygnałów
Laboratorium 8
Metoda najmniejszych kwadratów
Metoda najmniejszych kwadratów w identyfikacji układów
dynamicznych
Metoda najmniejszych kwadratów została wprowadzona przez Legendre a w 1805. Gauss,
który twierdził, że używał jej od 1794 r., wsparł ją w 1809 założeniem o rozkładzie błędów
normalnym (zwanym też rozkładem Gaussa-Laplace'a). Od Gaussa pochodzi nazwa: Methode
der kleinsten Fehlerquadrate (po polsku: metoda najmniejszych kwadratów błędów).
Początkowo była stosowana do obliczeń geodezyjnych, określających wielkość najbardziej
prawdopodobną z wielu nie całkiem zgodnych pomiarów. Stała się podstawą teorii błędów
pomiarów, używanej początkowo w astronomii i geodezji, obecnie we wszystkich pomiarach
fizycznych. Legła też u podstaw statystyki. Szerokie jej stosowanie wpłynęło na uproszczenie
nazwy, która jest nieco myląca, ale międzynarodowa (pl: Metoda najmniejszych kwadratów,
niem.: Methode der kleinsten Quadrate, ang.: Method of least squares, ros.: Mietod
najmieńszych kwadratow).
Opiera się na postulacie Legendre a. W postaci najprostszej postulat ten brzmi tak: wartością
najbardziej prawdopodobną, otrzymaną z szeregu wyników tak samo dokładnych pomiarów,
jest taka, od której obliczone odchylenia tych wyników, po podniesieniu do drugiej potęgi i
zsumowaniu dają wielkość najmniejszą z możliwych. Czyli przyjęcie do obliczenia odchyleń
wielkości dowolnej innej, niż najbardziej prawdopodobna, da sumę ich drugich potęg
(kwadratów) większą. Z postulatu Legendre'a wynika, że najbardziej prawdopodobną
wielkością z szeregu jednakowo dokładnych pomiarów jednej wielkości jest ich średnia
zwykła. W przypadku pomiarów niejednakowo dokładnych postulat ten brzmi podobnie,
stosuje się jednak do odchyleń równoważonych "wagami", tj wartość ma tym większą wagę
im bardziej dokładny jest pomiar. W tym przypadku najbardziej prawdopodobną okazuje się
wielkość zwana średnią ważoną. Gdy w zadaniu jest wiele niewiadomych, a nie są dostępne
bezpośredniemu pomiarowi, muszą być obliczane jako funkcje wielu innych mierzonych
wielkości. Wówczas do obliczeń stosuje się jeszcze bardziej rozwinięty aparat tej metody.
Wsparta założeniem o rozkładzie błędów normalnym i nazywana w Polsce rachunkiem
wyrównawczym, daje ona też liczbowo określone miary błędności wyników, jako ich tzw.
błędy średnie (które są przybliżeniami teoretycznych wielkości statystycznych  odchyleń
standardowych). Leżący u podstaw tej metody postulat Legendre'a nie wynika z żadnej ścisłej
matematycznej teorii. Jakkolwiek poczyniono wiele prób, by udowodnić jego słuszność i
uzasadnić stosowanie, wszystko spełzło na niczym. Pozostaje nadal tym, czym był od
początku  założeniem matematyka o genialnej intuicji. Geodeta i astronom, Tadeusz
Banachiewicz, napisał, że stosuje się metodę minimum sumy kwadratów :nie dlatego, abyśmy
uważali ją za matematycznie pewną, ale dlatego, że nikt dotychczas nie wskazał lepszej
metody. Tak naprawdę lepsze metody istnieją (np. regresja medianowa, albo metody oparte o
głębokość regresji), ale nie są tak proste obliczeniowo, więc się nie przyjęły
Definicja metody najmniejszych kwadratów
Metoda najmniejszych kwadratów (bardziej odpowiednia, ale nieużywana nazwa: metoda
minimum sumy kwadratów błędów)  metoda statystycznych estymacji i wyznaczania linii
trendu na podstawie zbioru danych w postaci par liczb.
Metoda najmniejszych kwadratów jest wykorzystywana do poszukiwania wielomianu h (x),
n
który będzie najdokładniejszą aproksymacją funkcji y(x), opisującej proces stochastyczny tak,
aby:
N
2
( yi - hn( xi )) = min
(1)
"
i= 1
gdzie: n - stopień wielomianu; N - liczba wyników doświadczenia i n<Aby znalez
ć minima tej funkcji ze względu na niewiadome współczynniki wielomianu a ,
i
różniczkuje się ją kolejno po a i przyrównuje do 0. Doprowadza to w efekcie do układu
i
równań, przy pomocy którego wyznaczamy niewiadome współczynniki wielomianu a
i.
S0a0 + S1a1 + �"+ Snan = To
S1a0 + S2a1 + �"+ Sn+ 1an = T1
(2)
�ą�ą�ą�ą�ą�ą�ą�ą�ą�ą�ą�ą
Sna0 + Sn+ 1a1 + �"+ S2nan = Tn
gdzie:
N N
Si = xij oraz Ti = xij y
,
" " j
j= 0 j= 0
Metoda najmniejszych kwadratów najczęściej jest stosowana przy regresji liniowej. Dla
regresji liniowej żądamy minimalizacji funkcji �2, która mierzy odchylenie zadanej zależności
funkcyjnej od punktów doświadczalnych. W przypadku funkcji liniowej f(x) = ax + b, funkcja
�2 sprowadza się do:
2
N
( yi - axi - b)
2
� (a,b) =
(3)
"
2
�
i= 1
i
gdzie � to odchylenie standardowe (niepewność pomiaru) danego punktu pomiarowego (w
i
zmiennej y). Aby znalez
ć minima tej funkcji ze względu na parametry a i b, różniczkuje się
po a i po b i przyrównuje do 0.
2
N
" � xi( yi - axi - b)
= 0 = - 2
(4)
"
2
" a �
i= 1
i
2
N
" � ( yi - axi - b)
= 0 = - 2
(5)
"
2
" b �
i= 1
i
Gdy odchylenie standardowe (niepewność pomiaru) wszystkich punktów pomiarowych jest
jednakowa, regresję nazywa się regresją nieważoną (klasyczną lub pierwszego rodzaju),
wówczas odchylenie standardowe może być wyłączone przed znak sumowania i upraszcza się
we wzorach na współczynniki a .
i
Ograniczenia metody najmniejszych kwadratów
Metoda najmniejszych kwadratów zawsze daje wynik o najmniejszej sumie kwadratów błę-
dów. Nie ma jednak gwarancji, że wynik ten ma jakikolwiek praktyczny sens. W szczególno-
ści, jeśli w danych występuje dużo elementów odstających, rezultaty mogą nie mieć nic
wspólnego z rzeczywistą linią trendu, czy zależnością między zjawiskami opisywanymi przez
zmienne losowe. Metoda najmniejszych kwadratów dostosowuje się bowiem do punktów naj-
bardziej oddalonych od średniej, które mogą wprowadzić największy błąd. Jeśli mamy w da-
nych pojedynczą zakłócającą obserwację (outlier) bardzo oddaloną od reszty, przyciągnie ona
do siebie linię trendu. Takie zjawisko jest niestety częste w realnych danych, nie należy więc
stosować metody najmniejszych kwadratów bez sprawdzenia (choćby na wykresie rozrzutu)
braku elementów odstających i ich usunięcia.
PRZYKA
AD 1
Metodą najmniejszych kwadratów znalez
ć równanie prostej, która najlepiej przybliża poniż-
sze dane:
x 1 2 3 4
i
y 7 6 4 1
i
Poszukiwane jest takie równanie prostej y = ax + b, aby funkcja:
f(a, b) = ((a � 1 + b) - 7)2 + ((a � 2 + b) - 6)2 + ((a � 3 + b) - 4)2 + ((a � 4 + b) - 1)2
osiągnęła minimum. W pierwszym kroku wyliczone zostaną pochodne minimalizowanej
funkcji:
" f
= 2 �" (a + b - 7) + 2 �" (2a + b - 6) �" 2 + 2 �" (3a + b - 4) �" 3 + 2 �" (4a + b - 1) �" 4 = 2 �" (30a + 10b - 35) = 0
" a
" f
= 2 �" (a + b - 7) + 2�" (2a + b - 6) + 2�" (3a + b - 4) + 2 �" (4a + b - 1) = 2�" (10a + 4b - 18) = 0
" b
Z powyższego układu równań można wyznaczyć parametry funkcji:
a = - 2,
b = 9,5.
W celu graficznego sprawdzenia poprawności rozwiązania wykonano w Matlabie następujący
rysunek:
l i n e a r r e g r e s s i o n
8
7
6
5
4
3
2
1
0
0 0 . 5 1 1 . 5 2 2 . 5 3 3 . 5 4 4 . 5 5
x
PRZYKA
AD 2
Metodą najmniejszych kwadratów znalez
ć parametry fizyczne układu oscylatora
harmonicznego jak na rysunku:
P
m
k
ą
Ł
Równanie ruchu układu można zapisać w postaci: m� + ą � + kx = P .
Ł
m�i +
Pojedynczy element pomiaru spełnia zatem równanie: ą �i + kxi = Pi
.
Tworzymy funkcję błędu e:
n
2
Ł
(m�i + ą �i + kxi - Pi) = e
.
"
i= 1
Poszukujemy następnie minimum funkcji błędu względem parametrów m, ą, k.
n
ńł " e
Ł Ł
= 0 �! 2 �" (m�i + ą �i + kxi - Pi ) �"�i = 0
"
�ł
" m
i= 1
�ł
n
" e
�ł
Ł
= 0 �! 2 �" (m�i + ą �i + kxi - Pi ) �"�i = 0
�ł
"
" ą
i= 1
�ł
n
" e
�ł
Ł
= 0 �! 2 �" (m�i + ą �i + kxi - Pi ) �" xi = 0
"
�ł
" ą
ół i= 1
y
Dokonujemy wymnożenia i pogrupowania składników równań.
n n n n
ńł
2
Ł
m �i + ą Ł Ł Ł
�i�i + k �ixi - Pi�i = 0
" " " "
�ł
i= 1 i= 1 i= 1 i= 1
�ł
n n n n
�ł
2
Ł
m �i�i + ą �i + k �ixi - Pi�i = 0
�ł
" " " "
i= 1 i= 1 i= 1 i= 1
�ł
n n n n
�ł 2
Ł
m �ixi + ą xi�i + k xi - Pixi = 0
" " " "
�ł
ół i= 1 i= 1 i= 1 i= 1
n n n n
�ł łł �ł łł
2
Ł Ł Ł Ł
�i �i�i �ixi śł Pi�i śł
" " " "
�ł �ł
i= 1 i= 1 i= 1 i= 1
m
�ł łł
�ł śł �ł śł
n n n n
2 �ł śł
�ł �ł
Ł
�i�i �i �ixi śł �" ą = Pi�i śł
" " " "
�ł śł
�ł śł �ł śł
i= 1 i= 1 i= 1 i= 1
�ł śł
k
�ł n n n śł �ł n śł
�ł �ł
2
Ł
�ł �ixi �ixi xi śł �ł Pixi śł
" " " "
�ł śł �ł śł
�ł i= 1 i= 1 i= 1 �ł �ł i= 1 �ł
Po dokonaniu powyższych czynności otrzymujemy proste równanie macierzowe X a = P,
które możemy numerycznie rozwiązać przy pomocy oprogramowania Matlab.
Zadania do samodzielnego wykonania:
1. Zbudować model układu mechanicznego z przykładu 2 w Simulinku dla następujących
parametrów:
m = 2,
c = 0,3,
k=1000.
2. Dokonać symulacji układu (najlepiej sygnałem szumu białego lub sygnałem poliharmo-
nicznym (złożenie wielu sinusoid o różnych częstotliwościach) i zarejestrować przebiegi
Ł
� , � i x.
3. Opierając się na równaniach wyprowadzonych w przykładzie 2 oraz zarejestrowanych
Ł
� �
przebiegach zmiennych , , x, zidentyfikować parametry układu m, ą, k.
4. Porównać zadane i zidentyfikowane parametry układu. (Jeśli parametry są takie same 
obliczenia są poprawne, jeśli są różne należy znalez
ć błąd w programie).
Ł
5. Dodać do przebiegów � , � i x szum na poziomie 0.1 wartości amplitudy tych przebie-
gów i dokonać ponownej identyfikacji układu i porównania parametrów.