Implementazione del problema della approssimazione ai minimi quadratiCamillo BoscoCorso di Analisi Numerica A.A. 2004-2005

Definizione del problema (richiamo)CASO DISCRETO: siano dati m punti (x1,y1),,(xm,ym). Si vuole determinare un polinomio p*(x) appartenente a Pn, con m>n ed m,n appartenenti a N, tale che: risulti minima rispetto ai coefficienti del polinomio. I valori w1, . . . , wm sono delle costanti positive dette pesi.A scopo didattico tali costanti assumono tutte valore 1 (peso uniforme).

Fitting Lineare: una istanza del problemaNel caso in cui n=1 vogliamo determinare il polinomio di primo grado p(x)=ax+b (geometricamente una retta) che costituisce la migliore approssimazione ai minimi quadrati. Vogliamo cio, noti m punti (xi,yi) i=1m trovare i valori di a e b che minimizzano:Come minimizzare tale funzione ? Calcolando le derivate prime di F rispetto ad a e rispetto a b e ponendole uguali a zero.Si ottiene cos un sistema di due equazioni in due incognite che risolto consente di esprimere a e b (le incognite !) in funzione delle coordinate degli n punti.

Le soluzioni del sistema

Esempio 1 : fitting_lineare.m

Ricordiamo la teoria !Risolvere il problema discreto ai minimi quadrati equivale a determinare la soluzione di un sistema lineare sovradeterminato (m>n) ottenuto calcolando i valori p*(xi).Qualora tale sistema Ax=b non abbia soluzione si determina il vettore soluzione:Tale vettore, che minimizza la somma dei quadrati delle componenti del vettore resto R = Ax-b, soluzione del sistema: ATAx = ATbIl problema viene quindi ricondotto alla risoluzione di un sistema lineare classicoRisolvere il problema discreto ai minimi quadrati equivale a determinare la soluzione di un sistema lineare sovradeterminato (m>n) ottenuto calcolando i valori p*(xi).Qualora tale sistema Ax=b non abbia soluzione si determina il vettore soluzione

Ripassiamo il metodo di Jacobi !E un metodo iterativo per la risoluzione di un sistema lineare quadrato (m=n).Al passo k+1-esimo le componenti del vettore soluzione sono cos definite: i=1nMATLAB risolve un sistema lineare utilizzando loperatore \ nella forma x=A\b.Tale operatore corrisponde alla funzione mldivide. Matlab utilizza un metodo diretto: il metodo di eliminazione di Gauss

Esempio 2 : min_quad.m

Definizione del problema (caso continuo)CASO CONTINUO: sia w(x) una funzione positiva, continua ed integrabile in [a,b]. Si vuole determinare p*(x) appartenente a Pn in modo tale da minimizzare: La funzione w(x) detta funzione peso. A scopo didattico si sceglie w(x)=1 cio peso uniforme ed unitario.

Una istanza del problemaNel caso in cui w(x)=1 e consideriamo lo spazio V=C[a,b], ovvero lo spazio delle funzioni continue in [a,b], vogliamo trovare i coefficienti a0,,an del polinomio p*(x) che approssima f appartenente a V.Dobbiamo minimizzare:Come minimizzare tale funzione ? Calcolando le derivate prime di F rispetto ai valori a0,,an e ponendole uguali a zero. Otteniamo il sistema:

Esempio 3 : min_quad_continuo1.mIn tal caso [a,b] = [0,1] il sistema diventa:PROBLEMA:La matrice dei coefficienti di tale sistema :i,j= 0nSi tratta della matrice di Hilbert. E una matrice malcondizionata !

Soluzione !Per evitare di ottenere una matrice malcondizionata e difficilmente invertibile si sceglie una base differente da quella canonica in Pn. Linsieme S={1,x,x2,,xn}, base canonica di Pn, costituito da polinomi non ortogonali rispetto a nessun prodotto interno.Vogliamo trasformare la base canonica in una famiglia di polinomi ortogonali.Il procedimento di ortogonalizzazione di Gram-Schmidt ci consente di generare due successioni di polinomi:

Ortogonalizzazione di Gram-Schmidt

Le due successioniLa successione costituita da polinomi ortogonali. Pertanto a norma di definizione: La successione costituita da polinomi ortonormali. Quindi: Scegliendo tali successioni la matrice del sistema di equazioni diagonale, quindi non pi malcondizionata !

Esempio 4 : min_quad_continuo2.mIl polinomio di migliore approssimazione ai minimi quadrati nel caso continuo calcolato come segue:1) Si genera linsiemeutilizzando il procedimento di ortogonalizzazione di Gram-Schmidt2) Si calcolano i coefficienti generalizzati di Fourier:j= 0n3) Si costruisce il polinomio di approssimazione:

Teoria dei polinomi ortogonali-Esistono delle famiglie di polinomi ortogonali utilizzate in diversi ambiti della analisi numerica.Ciascuna famiglia caratterizzata da una propriet principale: ciascun polinomio pn della famiglia ortogonale a tutti i polinomi di grado minore o uguale a n.Nel caso della approssimazione si osserva che, scegliendo w(x)=1 in [-1,1], si ottiene la famiglia dei polinomi di Legendre:-Scegliendo invece w(x)=1/sqrt(1-x2) in [-1,1], si ottiene la famiglia dei polinomi di Chebichev:

Esempio 5 : polinomiLegendre.mGenera i primi n polinomi di Legendre secondo la seguente ricorrenza: (Rif. Naldi-Pareschi-Russo pag. 262)L0 (x)= 1L1 (x)= x

(n+1) Ln+1 (x)=(2(n+1)-1)xLn(x)-nLn-1(x)

Esempio 6 : polinomiChebichev.mGenera i primi n polinomi di Chebichev secondo la seguente ricorrenza:T0 (x)= 1T1 (x)= x

Tn+1 (x)=2xTn(x)-Tn-1(x)