Minimi quadrati e polinomi ortogonali

30 gennaio 2009, release B

1 Il problema ai minimi quadrati

Consideriamo lo spazio normato delle funzioni misurabili (cf. [5, p.284]) quadrato integrabili ( ^{2}( ) _{2}) dove è un intervallo della retta reale, non necessariamente limitato (cf. [5, p.386], ricordando [5, p.308]), e

2 ( ) ( )_{2} _{}^{}( ) \overline{( )}

dove \overline{}

Si dimostra (non facile!) che questo spazio euclideo è un esempio di spazio di Hilbert (cf. [5, p.388]), cioè uno spazio euclideo che è completo, separabile e infinito dimensionale (cf. [5, p.155]). Più in generale se : ( ) � è una funzione positiva allora lo spazio ( _{}^{2} ( )

2) definito come

_{}^{2} ( ) misurabili t.c. _{}^{}( ) ^{2}

è uno spazio di Hilbert dotato del prodotto scalare

( )_{2} _{}^{}( ) \overline{( )} ( )

(cf. [2, p.23]). Supponiamo di seguito che in particolare sia : ( ) � una funzione nonnegativa, con ( ) non necessariamente limitato, tale che

1. _{}^{}^{}( ) per tutti gli �;

2. _{}^{}( ) ( ) 0 per qualche funzione continua e non negativa implica 0 in ( ). per cui, visto che _{}^{2} ( ) è uno uno spazio vettoriale, contiene sia ^{}_{�} come pure le loro combinazioni lineari e quindi lo spazio dei polinomi _{}di grado (con � arbitrario). Infatti, se _{}( ) _{0}^{}_{}^{}allora per la disuguaglianza triangolare e il fatto che per ogni si ha

_{2}^{2} _{}^{}^{2} ( )

necessariamente

2 ; ; ; ; _{0}^{}_{}^{} ; ; ; ; _{2} _{0}^{}_{}^{}_{2}

e quindi _{}.

Fissata _{}^{2} ([ ]) ed �, il problema ai minimi quadrati (nel continuo) consiste nel determinare il polinomio _{}di grado tale che sia minima la quantità (cf. [1, p.204-207])

2 _{}^{}( ) _{}( ) ^{2} ( )

Essendo _{}^{2} ([ ]) uno spazio euclideo, e ( _{})_{0} una base di _{}, abbiamo visto che la soluzione del problema

2 min_{span _{0}}  _{2}

è

_{0}^{}_{}^{}_{}

dove i coefficienti _{}^{}verificano le cosidette equazioni normali

_{0}^{}( _{}_{})_{2} _{}^{}( _{})_{2} 0

La soluzione è caratterizzata dalla proprietà di ortogonalità cioè che ^{}è ortogonale a tutti gli _{}, con 1 , ovvero

( _{})_{2} ( ^{}_{})_{2} 0

Nel caso _{}( ) ^{}per 0 , le equazioni normali si riscrivono quindi come

_{0}^{}_{}^{}_{}^{}^{}( )  _{}^{}^{}( ) ( )  0

o in forma matriciale, posto per 0

: _{}^{}^{}( ) _{}_{}^{}^{}( ) ( ) come . Nota facoltativa. Quanto visto si può provare in un modo alternativo. Sia : ( ) �, _{}^{2} ( ) e definiamo

( _{0} _{1}  _{}) : _{}^{}( ) _{0}^{}_{}^{}^{2} ( )  _{}^{}[ ( ) _{0}^{}_{}^{}]^{2} ( )

Sviluppando il quadrato [ ( ) _{0}^{}_{}^{}]^{2} si vede subito che ( _{0} _{1} _{}) è infinitamente differenziabile rispetto i coefficienti ( _{0} _{1} _{}) e inoltre che una condizione necessaria affinchè ( _{0} _{1} _{}) sia un minimo per è che sia

\frac{}{_{}}0 0

Osserviamo che

0 \frac{}{_{}}_{}^{}\frac{[ ( ) _{0}^{}_{}^{}]^{2} ( )}{_{}}
_{}^{}2 [ ( ) _{0}^{}_{}^{}] ( ^{}) ( )

_{}^{}2 [ ( ) ^{}_{0}^{}_{}^{}] ( )

se e solo se

_{}^{}( ) ^{}( ) _{0}^{}_{}_{}^{}^{}( )

La quantità a primo membro coincide con ( ^{})_{2}. Per ipotesi, essendo _{}^{}^{} ( ) la quantità a secondo membro e’ ben definita. Sia _{}_{}^{}^{} ( ) , ( _{}), _{}( ^{( 1)})_{2} con 0 2 . Quindi la soluzione esiste ed unica se e solo se lo è quella del sistema lineare cioè det( ) 0.

Vediamo quale caso particolare ( ) 1, ( ) (0 1). Poichè

_{0}^{1} ^{}\frac{1}{1}

si nota che la matrice coincide con la matrice di Hilbert _{}\frac{1}{1} . La matrice di Hilbert è estremamente sensibile a piccoli cambiamenti nei coefficienti o dei valori del termine noto. Per vederlo (cf. [1, p.530-535]), sia una matrice non singolare, una norma indotta di matrice [8], e

( ) ^{1}

Si dimostra che se

(  )(  )

allora

\frac{}{}\frac{( )}{1 ( )\frac{}{}}\frac{}{}\frac{}{}

Per capirne il significato, si supponga per semplicità 0, Allora

\frac{}{}\frac{( )}{1 ( )\frac{}{}}\frac{}{}\frac{}{}( )\frac{}{}

per cui più grande è ( ) più grande è il rischio che piccole perturbazioni sul termini noto provochino grandi effetti sul calcolo della soluzione . Si noti che questa stima è solo una maggiorazione. Comunque in molti casi risulta un buon indicatore di quanto un sistema lineare sia difficile da risolvere al calcolatore, indipendentemente dal metodo utilizzato.

Tornando alla matrice di Hilbert, vediamo in Matlab valori del numero di condizionamento in norma 2.

>> help cond
COND   Condition number with respect to inversion.
COND(X) returns the 2-norm condition number (the ratio of the
largest singular value of X to the smallest).  Large condition
numbers indicate a nearly singular matrix.
COND(X,P) returns the condition number of X in P-norm:
NORM(X,P) * NORM(INV(X),P).
where P = 1, 2, inf, or 'fro.'
See also RCOND, CONDEST, CONDEIG, NORM, NORMEST.
>> for n=5:5:40, fprintf('\n \t [n]: %2.0f [cond (norma 2)]: %2.2e
[det]: %2.2e', n,cond(hilb(n)),det(hilb(n))); end
[n]:  5 [cond (norma 2)]: 4.77e+005 [det]: 3.75e-012
[n]: 10 [cond (norma 2)]: 1.60e+013 [det]: 2.16e-053
[n]: 15 [cond (norma 2)]: 8.49e+017 [det]: -2.19e-120
[n]: 20 [cond (norma 2)]: 1.91e+018 [det]: -1.12e-195
[n]: 25 [cond (norma 2)]: 1.46e+019 [det]: 8.14e-275
[n]: 30 [cond (norma 2)]: 6.19e+018 [det]: 0.00e+000
[n]: 35 [cond (norma 2)]: 1.03e+019 [det]: 0.00e+000
[n]: 40 [cond (norma 2)]: 4.55e+019 [det]: 0.00e+000
>> for n=5:5:40, fprintf('\n \t [n]: %2.0f [min autoval.]: %2.2e
[max autoval.]: %2.2e', n,min(abs(eig(hilb(n)))),max(abs(eig(hilb(n))))); end
[n]:  5 [min autoval.]: 3.29e-006 [max autoval.]: 1.57e+000

          [n]: 10 [min autoval.]: 1.09e-013 [max autoval.]: 1.75e+000
          [n]: 15 [min autoval.]: 7.58e-018 [max autoval.]: 1.85e+000
          [n]: 20 [min autoval.]: 3.05e-018 [max autoval.]: 1.91e+000
          [n]: 25 [min autoval.]: 1.70e-018 [max autoval.]: 1.95e+000
          [n]: 30 [min autoval.]: 9.49e-019 [max autoval.]: 1.99e+000
          [n]: 35 [min autoval.]: 1.68e-017 [max autoval.]: 2.01e+000
          [n]: 40 [min autoval.]: 2.22e-017 [max autoval.]: 2.04e+000
>>

Ricapitolando

1. il problema ai minimi quadrati continui ha un’unica soluzione se e solo se un certo sistema lineare ha una e una sola soluzione;

2. nel caso 1, ( ) ( 1 1) la matrice coincide con la matrice di Hilbert che ha un numero di condizionamento estremamente grande, il che comporta che il sistema lineare sarà estremamente sensibile a piccole perturbazioni (cioè malcondizionato);

3. la matrice di Hilbert sembra essere singolare, cioè avente determinante nullo, per 0; in realtè si dimostra che è definita positiva, e quindi tutti gli autovalori sono positivi da cui il determinante è non nullo essendo il prodotto degli autovalori.

Notiamo subito che il sistema è difficile da risolvere poichè abbiamo rappresentato il polinomio _{}^{}nella base 1 ^{}. Cosa succede se lo scriviamo come combinazione lineare di elementi di un’altra base _{}? Se ad esempio scegliessimo la base _{}cosicchè

( _{}_{}) _{}_{}_{}0 0 1

con _{}delta di Kronecker, finiremmo per risolvere un problema dove è una matrice diagonale e quindi senza alcuna difficoltà troveremmo

\frac{_{}}{_{}}

con _{}0.

Risulta quindi di importanza fondamentale cercare basi di polinomi _{}per cui

( _{}_{}) _{}_{}_{}0 0 1

2 Polinomi ortogonali

Sia : ( ) � una funzione nonnegativa, con ( ) non necessariamente limitato, tale che

1. _{}^{}^{}( ) per tutti gli �; 2. _{}^{}( ) ( ) 0 per qualche funzione continua e non negativa implica

0 in ( ).

La funzione è detta peso. Le funzioni peso più comuni sono (cf. [1, p.206], [10], [11],[9])

1. ( ) 1 con [ 1 1] (peso di Legendre) [14];

2. ( ) \frac{1}{\overline{1 ^{2}}} con ( 1 1) (peso di Chebyshev);

3. ( ) (1 ^{2})^{(1/2)} con ( 1 1), ( 1/2) (peso di Gegenbauer) [6];

4. ( ) (1 )^{}(1 )^{}con ( 1 1), 1, 1 (peso di Jacobi) [7];

5. ( ) exp( ) con (0 ) (peso di Laguerre) [13];

6. ( ) exp( ^{2}) con ( ) (peso di Hermite) [12];

Vediamo ora che proprietà ha un insieme di polinomi _{}_{0} in cui il grado deg di ogni polinomio sia 1 deg( _{}) e ( _{}_{})_{2} _{}_{}(con _{}il delta di Kronecker, _{}0, 0 ). Una tal famiglia triangolare di polinomi (cioè tale che deg( _{}) ) si dice ortogonale rispetto alla funzione peso nell’intervallo di riferimento.

Si può dimostrare usando la procedura di Gram-Schmidt che una tal famiglia triangolare di polinomi esiste e con la stessa procedura costruirla direttamente; inoltre è immediato osservare che ogni polinomio di grado si può scrivere univocamente come combinazione lineare di _{0} _{}. Di conseguenza se _{}_{0}^{}_{}_{}, allora per la bilinearità del prodotto scalare ( )_{2}

( _{1} _{})_{2} ( _{1} _{0}^{}_{}_{})_{2} _{0}^{}_{}( _{1} _{})_{2} 0

Inoltre (cf. [3, p.978])

Teorema 2.1 Sia _{}_{0} una famiglia triangolare di polinomi ortogonali in ( ) rispetto ad una funzione peso . Allora il polinomio _{}ha esattamente radici reali e distinte nell’intervallo aperto ( ).

Prima dimostrazione. Siano _{1} _{}tutti gli zeri di _{}interni ad ( ) tali che _{} cambia segno in ogni _{}con 0 . Questa richiesta va discussa.

• Se _{}( _{}) 0 ma non cambia segno, allora necessariamente _{}^{(1)}( _{}) 0 e quindi lo zero _{}ha molteplicità maggiore di 1.

• In virtù del fatto che nella richiesta gli zeri distinti interni all’intervallo sono e con molteplicità uguale a 1, se fosse allora le rimanenti radici sarebbero in � ( ). Mostreremo che così non è .

Supponiamo sia per assurdo e che sia

( ) : ( _{1}) ( _{})

Sia ( ) _{}( ) _{}( )

con _{}_{}che per ipotesi, non avendo zeri in ( ), ha segno costante in . Quindi

( ) _{}( ) _{}( ) _{}( ) _{}( ) _{}^{2} ( ) _{}( )

è una funzione continua che ha segno costante in ( ), non è nulla, per cui dalle proprietà della funzione peso si ha

_{}^{}_{}( ) _{}( ) ( ) 0

Ma poichè deg( _{}) , essendo ( _{})_{2} 0 per ogni polinomio di grado strettamente inferiore a , allora

_{}^{}_{}( ) _{}( ) ( ) 0

il che è assurdo.

Seconda dimostrazione (simile alla prima ma eseguita durante il corso). Siano 1 _{}tutti gli zeri di _{}interni ad ( ) con molteplicità rispettivamente _{1} _{}. Di conseguenza, per qualche numero _{}abbiamo

( ) _{}_{1}^{}( _{})^{_{}}_{}( )

avendo supposto _{1}^{}( _{})^{_{}}1 se non ci sono zeri interni ad ( ). Il polinomio _{}per costruzione non ha zeri in ( ) e quindi non si annulla mai ed essendo una funzione continua ha segno costante.

Consideriamo il polinomio

( ) _{1}^{}( _{})^{mod_{2}( _{})}

Se uno zero di _{}ha molteplicità dispari ma maggiore di 1 o uno almeno ha molteplicità pari o esiste uno zero complesso non in ( ), è facile osservare che il grado di è minore di . Osserviamo ora che qualsiasi sia un numero naturale, _{}mod_{2}( _{}) è un numero pari.

Di conseguenza avendo _{1}^{}( _{})^{_{}mod_{2}( _{})} _{}( ) segno costante e non coincidente col polinomio nullo,

0  ( _{})
_{}^{}_{}( ) ( ) ( )
_{}^{}_{}_{1}^{}( _{})^{_{}}_{}( ) _{1}^{}( _{})^{mod_{2}( _{})} ( )

_{}^{}_{}_{1}^{}( _{})^{_{}mod_{2}( _{})}  _{}( ) ( )  0

(1)

da cui la contraddizione. Potrebbe venire il dubbio su perchè qualche non possa essere o . Nella dimostrazione avrebbe quale unico effetto che _{}si annulla in o , rimanendo di segno costante in ( ). Siccome l’integrale di due funzioni non cambia se sono uguali a meno di un insieme trascurabile di punti, questo effetto è del tutto inessenziale. La conclusione è che il polinomio ortogonale _{}ha radici distinte e semplici, interne ad ( ).

Vediamo ora il Teorema di Clenshaw, che stabilisce a partire dai polinomi ortogonali _{0}, _{1} di grado 0 e 1 rispettivamente come calcolare ricorsivamente la famiglia di polinomi ortogonali.

Teorema 2.2 Sia _{}_{0} una famiglia triangolare di polinomi ortogonali in [ ] rispetto ad una funzione peso . Allora per 1

1 ( ) _{}( _{}) _{}( ) _{}_{1} ( )

dove, detto _{}il coefficiente di grado massimo di _{}, si ha

\frac{_{1}}{_{}}(2)
\frac{( _{}_{})_{2}}{( _{}_{})_{2}}  (3)
\frac{_{}( _{1} _{})_{2}}{( _{}_{})_{2}}  (4)

Dimostrazione: Facoltativa. Per qualche _{}, _{}, _{}, _{2} _{2} abbiamo

1 ( ) _{}( _{}) _{}( ) _{}_{1} ( ) _{0}^{2} _{}_{}( )

Essendo

( _{1} _{})_{2} 0 1

necessariamente, per la bilinearità del prodotto scalare,

( _{1} _{})_{2}  ( _{}( _{}) _{}_{}_{1} _{0}^{2} _{}_{}_{})_{2}

( _{}( _{}) _{}_{})_{2} ( _{}_{1} _{})_{2} _{0}^{2} _{}( _{}_{})_{2} (5)

Se

• se 2, abbiamo, essendo il grado di _{}( _{}) _{}uguale a 1 1 , ( _{1} _{})_{2} 0, sia

_{0}^{2} _{}( _{}_{})_{2} _{}

sia

( _{}_{1} _{})_{2} 0

come pure

( _{}( _{}) _{}_{})_{2}  _{}^{}_{}( _{}) _{}( )  _{}( ) ( )
_{}^{}_{}( ) _{}( _{}) _{}( ) ( )
( _{}_{}( _{}) _{})_{2}
0  (6)

per cui da (5) _{}0 per 0 2. Di conseguenza

1 ( ) _{}( _{}) _{}( ) _{}_{1} ( )  (7)

• da un confronto dei coefficienti di grado massimo nella formula ricorsiva di Clenshaw si deduce subito il valore di _{}.

• se , abbiamo ( _{}_{1} _{})_{2} 0, (( _{1} _{})_{2} 0 e quindi

0 ( _{}( _{}) _{}_{})_{2})

per cui, da _{}0 visto che 1 deg( _{}) e ( _{}_{})_{2} 0 , ricaviamo

( _{}_{})_{2}) ( _{}_{}_{})_{2}) _{}( _{}_{})_{2})

\frac{( _{}_{})_{2}}{( _{}_{})_{2}}

• se 1, da (7), ( _{}_{1} )_{2} 0, (( _{1} _{1} )_{2} 0, abbiamo

0  ( _{}( _{}) _{}_{1} )_{2} _{}( _{1} _{1} )_{2}
( _{}_{1} )_{2} _{}_{}( _{}_{1} )_{2} _{}( _{1} _{1} )_{2}
( _{}_{1} )_{2} _{}( _{1} _{1} )_{2}

(8) e quindi

\frac{_{}( _{}_{1} )_{2}}{( _{1} _{1} )_{2}}

References

[1] K. Atkinson, An Introduction to Numerical Analysis, Wiley, (1989).

[2] K. Atkinson e W. Han, Theoretical Numerical Analysis, A Functional Analysis Framework, Springer, (2001).

[3] V. Comincioli, Analisi Numerica, metodi modelli applicazioni, McGraw-Hill, (1990).

[4] G. Dahlquist e A. Bjorck , Numerical methods, Dover, (2003).

[5] A.N. Kolmogorov e S.V. Fomin, Introductory Real Analysis, Dover publications, 1970.

[6] MathWorld at Wolfram, (Gegenbauer polynomial), http://mathworld.wolfram.com/GegenbauerPolynomial.html .

[7] MathWorld at Wolfram, (Jacobi polynomial), http://mathworld.wolfram.com/JacobiPolynomial.html .

[8] Wikipedia, (Norma matriciale), http://it.wikipedia.org/wiki/Norma_matriciale .

[9] Wikipedia, (Orthogonal polynomials), http://en.wikipedia.org/wiki/Orthogonal_polynomials .

[10] Wikipedia, (Polinomi di Chebyshev), http://it.wikipedia.org/wiki/Polinomi_di_Chebyshev .

[11] Wikipedia, (Polinomi Ortogonali), http://it.wikipedia.org/wiki/Polinomi_ortogonali .

[12] Wikipedia, (Polinomi di Hermite), http://it.wikipedia.org/wiki/Polinomi_di_Hermite .

[13] Wikipedia, (Polinomi di Laguerre), http://it.wikipedia.org/wiki/Polinomi_di_Laguerre .

[14] Wikipedia, (Polinomi di Legendre), http://it.wikipedia.org/wiki/Polinomi_di_Legendre .