Základy algoritmizace a programování

Petr Sváček
Ústav technické matematiky

Cviceni č. 6

Aproximace metodou nejmenších čtverců
Pro zadanou tabulku hodnot 𝑥𝑖,𝑦𝑖 pro 𝑖=1,…𝑛 hledáme aproximaci této tabulky pomocí vhodné funkce 𝑝(𝑥) (např. polynomu 1. nebo 2. stupně).

Snažíme se minimalizovat kvadratickou odchylku (zde v druhé mocnině a pro p(x) polynom 2. stupně)
což lze zapsat maticově 𝛿(𝑝)=‖𝔸𝑢⃗ −𝑦⃗ ‖^2 kde 𝑢⃗  označuje vektor neznámých koeficientů u = [a; b; c], vektor 𝑦⃗ =( y_i) a matice 𝔸 je daná 
v MATLAB zápise

 
  A = [ones(n,1), x, x.^2 ];
 

Minima se nabývá pro vektor u, který je řešením soustavy rovnic

(𝔸^T 𝔸)𝑢⃗ =𝔸^T 𝑦⃗ 

Stejný postup lze použít pro lineární polynom, v takovém případě je matice A dána jako A = [ones(n,1), x]. Pro kubický polynom bychom matici A naopak rozšířili.

Rozmyslete si jak bychom postupovali pro data 𝑧𝑖, o kterých předpokládáme, že mají exponenciální růst, tedy 𝑧(𝑥)=˙𝑒^𝑎𝑥+𝑏. Užijeme-li logaritmus, dostaneme ln𝑧𝑖≈𝑎𝑥𝑖+𝑏.
a) Je daná tabulka hodnot x, y zatížených chybou měření. Tyto data proložte lineární funkcí tak, aby kvadratická odchylka bylo co nejmenší. Graf dat i získané funkce znázorněte.
Např. náhodná data - tabulka x, y

  x = (-4:5)'; n = length(x);
  coefA = round(40 * rand(1)) / 2 - 4;
  coefB = round(40 * rand(1)) / 2 - 9;
  y = (coefA * x + coefB + 10 * coefA * rand(n,1));


Lze aproximovat dle

   AA = [ones(n,1), x];
  
  reg = (AA' * AA) \ (AA' * y);
  
  h = 0.01; xx = x(1):h:x(end);
  
  yy = reg(2) * xx + reg(1);
  plot(x, y, 'bo', xx, yy);