Aproximarea prin derivare#

Simplificat, încercăm să găsim un polinom de un anumit grad fixat care să treacă printre punctele ce reprezintă evaluări ale unei funcții date și care să aproximeze funcția cât mai bine. Vom aborda acum subiectul puțin mai riguros.

O astfel de metodă poartă denumirea de aproximare în sensul celor mai mici pătrate. Vom detalia în subcapitolele următoare, momentan focalizându-ne pe cunoștințele deja dobândite.

Regresia liniară#

Fie \(M=\Big\{(x_k,y_k)\in\mathbb{R}^2\,|\,k\in\overline{0,n}\Big\}\), unde \(n\in\mathbb{N}^*\), o mulțime de \(n+1\) puncte distincte.

Vrem să găsim o dreaptă care să aproximeze cât mai bine funcția de la care provin aceste evaluări. Așadar, să considerăm polinomul \(P_1\in\mathbb{P}_1\) de grad \(1\), unde \(P_1:\mathbb{R}\rightarrow\mathbb{R}\):

\[ \boxed{P_1(x)=a_0+a_1x}\,,\,\,a_0,a_1\in\mathbb{R} \]

Pentru a considera problema rezolvată, este suficientă aflarea coeficienților \(a_0\) și \(a_1\) astfel încât \(P_1\) să aproximeze "cât mai bine" funcția numerică \(M\). Pentru aceasta, să definim următoarea funcție ce are ca scop evaluarea erorii, \(\varepsilon:\mathbb{R}^2\rightarrow\mathbb{R}\):

\[ \boxed{ \varepsilon(a_0,a_1)=\sum_{i=0}^n\Big[y_i-P_1(x_i)\Big]^2=\sum_{i=0}^n\Big[y_i-(a_1x_i+a_0)\Big]^2 } \]

Scopul metodei este de a minimiza această eroare. Din această formulă se explică proveniența denumirii metodei. Evident, pentru a minimiza \(\varepsilon(a_0,a_1)\), se vor aplica metodele de analiză matematică deja cunoscute, adică se vor egala derivatele parțiale în funcție de \(a_0\), respectiv în funcție de \(a_1\), cu \(0\). Se obține astfel următorul sistem:

\[ \begin{cases} \dfrac{\partial \varepsilon}{\partial a_0}(a_0,a_1)=0\\\\ \dfrac{\partial \varepsilon}{\partial a_1}(a_0,a_1)=0 \end{cases}\Rightarrow \begin{cases} -2\sum_{i=0}^n(y_i-a_1x_i-a_0)=0\\ -2\sum_{i=0}^n(y_i-a_1x_i-a_0)x_i=0 \end{cases} \]

Separând termenii, obținem următoarele două ecuații pe care le vom considera fundamentale:

\[ \begin{cases} \sum_{i=0}^ny_i=a_0n+a_1\sum_{i=0}^nx_i\\ \sum_{i=0}^nx_iy_i=a_0\sum_{i=0}^nx_i+a_1\sum_{i=0}^nx_i^2 \end{cases} \]

Vom extrage din prima ecuație valoarea lui \(a_0\), acesta rămânând însă în funcție de \(a_1\):

\[ a_0=\frac{\sum_{i=0}^ny_i-a_1\sum_{i=0}^nx_i}{n} \]

Dacă notăm cu \(\overline{x}\) și \(\overline{y}\) media valorilor \(x_k\), respectiv \(y_k\) (adică \(\overline{x}=\frac{\sum_{i=0}^nx_i}{n}\) și \(\overline{y}=\frac{\sum_{i=0}^ny_i}{n}\)), atunci formula anterioară se poate rescrie ca:

\[ \boxed{a_0\in\overline{y}-a_1\overline{x}} \]

Înlocuind acum în ecuația \(\sum_{i=0}^nx_iy_i=a_0\sum_{i=0}^nx_i+a_1\sum_{i=0}^nx_i^2\), se obține:

\[ a_1=\frac{\sum_{i=0}^nx_iy_i-\overline{y}\sum_{i=0}^nx_i+a_1\overline{x}\sum_{i=0}^nx_i}{\sum_{i=0}^nx_i^2} \]

Continuarea calculului

Se poate demonstra faptul că următoarele două relații sunt adevărate:

\[ \begin{cases} \sum_{i=0}^n(x_i-\overline{x})(y_i-\overline{y})=\sum_{i=0}^nx_iy_i-n\overline{x}\overline{y}\\ \sum_{i=0}^n(x_i-\overline{x})^2=\sum_{i=0}x_i^2-n\overline{x}^2 \end{cases} \]

Aceste demonstrații rămând drept un exercițiu pentru cititor.

În acest context, se obține următoarea formulă finală pentru \(a_1\):

\[ \boxed{ a_1=\frac{\sum_{i=0}^n(x_i-\overline{x})(y_i-\overline{y})}{\sum_{i=0}^n(x_i-\overline{x})^2} }\]

Vizualizarea algoritmului#

Pornim de la o dreaptă cunoscută, să zicem \(f(x)=2x-3\), iar apoi generăm aleator o mulțime de \(n=20\) puncte astfel încât \(M=\Big\{(x_k,y_k)\in\mathbb{R}^2\,|\,k\in\overline{0,n},\,y_k=f(x_k)\pm50\%f(x_k)\Big\}\). Suplimentar, vom alege valori ale lui \(x_k\) astfel încât \(x_k\in[10,90]\).

Nu vom plictisi cititorul cu calculele numerice, însă vom prezenta punctele alături de funcția inițială și de polinomul construit, \(P_1(x)\approx2.4175x-13.8107\) (după cum putem observa, acest polinom este foarte apropiat de funcția inițială):

Mulțimea \(M\) de puncte alese aleator

Mulțimea \(M\) de puncte peste care s-a suprapus funcția inițială (cu portocaliu) și aproximarea sa (cu mov)

Cod ilustrativ#

Implementarea CMMP (cu drepte) se poate face prin funcția aprox_CMMP. Aceasta primește:

x - abscisele din mulțimea \(M\);
y - ordonatele din mulțimea \(M\);

Codul de mai jos returnează [a0, a1], valori utilizate în polinomul de grad \(1\) ce aproximează funcția (numerică) \(M\):

\[ P_1(x)=a_0+a_1x \]

function [a0, a1] = aprox_CMMP(x, y)
    n = length(x);
    sx = sum(x);
    sy = sum(y);
    sxy = sum(x.*y);
    sx2 = sum(x.^2);
    a0 = (sx2*sy-sxy*sx) / (n*sx2-sx.^2);
    a1 = (n*sxy-sx*sy) / (n*sx2-sx.^2);
end

Aproximarea cu un polinom oarecare#

Să generalizăm acum agloritmul deja cunoscut.

Fie \(M=\Big\{(x_k,y_k)\in\mathbb{R}^2\,|\,k\in\overline{0,n}\Big\}\), unde \(n\in\mathbb{N}^*\), o mulțime de \(n+1\) puncte distincte. Căutăm polinomul \(P_\mu\in\mathbb{P}_\mu\) de grad \(\mu\in\mathbb{N}^*\) \(\mu\lt n\) (dacă am accepta \(\mu\geq n\), am obține o interpolare polinomială), unde \(P_\mu:\mathbb{R}\rightarrow\mathbb{R}\):

\[ \boxed{ P_\mu(x)=a_0+a_1x+\dots+a_\mu x^\mu=\sum_{k=0}^\mu a_kx^k }\,,\,\,a_0,\dots,a_\mu\in\mathbb{R} \]

Ca să putem afla valorile \(a_0,\dots,a_n\), vom generaliza funcția de evaluare a erorii deja descrisă anterior, adică definim \(\varepsilon:\mathbb{R}^{n+1}\rightarrow\mathbb{R}\):

\[ \boxed{ \varepsilon(a_0,\dots,a_\mu)=\sum_{i=0}^n\Big[y_i-P_\mu(x_i)\Big]^2 }\,= \sum_{i=0}^ny_i^2-2\sum_{i=0}^ny_iP_\mu(x_i)+\sum_{i=0}^n\Big[P_\mu(x_i)\Big]^2 \]

Folosind această formulă, dezvoltăm polinomul deja cunoscut:

\[\begin{split} \varepsilon(a_0,\dots,a_\mu)&=\sum_{i=0}^ny_i^2-2\sum_{i=0}^n\left(\sum_{j=0}^\mu a_jx_i^j\right)y_i+\sum_{i=0}^n\left(\sum_{j=0}^\mu a_jx_i^j\right)^2\\ &=\sum_{i=0}^ny_i^2-2\sum_{j=0}^\mu a_j\left(\sum_{i=0}^ny_ix_i^j\right)+\sum_{j=0}^\mu\sum_{k=0}^\mu a_ja_k\left(\sum_{i=0}^nx_i^{j+k}\right) \end{split}\]

Scopul algoritmului este de a găsi un polinom \(P_\mu\) ce minimizează funcția \(\varepsilon(a_0,\dots,a_\mu)\). Așadar, derivatele parțiale ale lui \(\varepsilon\) în raport cu \(a_j\), \(\forall j\in\overline{0,\mu}\), trebuie egalate cu \(0\). Nu este greu de dovedit că acestea vor lua forma:

\[ \frac{\partial\varepsilon}{\partial a_j}(a_0,\dots,a_\mu)=-2\sum_{i=0}^ny_ix_i^j+2\sum_{k=0}^\mu a_k\sum_{i=0}^nx_i^{j+k}=0,\,\forall j\in\overline{0,\mu} \]

Astfel, luând fiecare valoarea în parte posibilă a lui \(j\in\overline{0,\mu}\), se vor forma \(\mu+1\) ecuații liniare cu \(\mu+1\) necunoscute (anume coeficienții polinomului, \(a_j\)). Ecuațiile componente vor lua forma:

\[ \boxed{\sum_{k=0}^\mu a_k\sum_{i=0}^nx_i^{j+k}=\sum_{i=0}^ny_ix_i^j}\,,\,\,\forall j\in\overline{0,\mu} \]

...trivial.

Licență#

The book "Metode Numerice", written by Valentin-Ioan Vintilă, is licensed under CC BY-NC-SA 4.0