From a075c778cba02032a54c8fc984111e2341913dc7 Mon Sep 17 00:00:00 2001 From: Leonardo Robol Date: Thu, 5 Nov 2009 18:29:18 +0100 Subject: [PATCH] =?UTF-8?q?Aggiunta=20lezione=20di=20Gemignani=20del=204?= =?UTF-8?q?=20novembre=20che=20in=20realt=C3=A0=20non=20avevo=20mai=20aggi?= =?UTF-8?q?unto=20perch=C3=A8=20mi=20ero=20dimenticato=20di=20aggiungere?= =?UTF-8?q?=20il=20file=20capitolo3.tex..=20ops!?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit --- CalcoloScientifico.tex | 1 + capitolo3.tex | 163 +++++++++++++++++++++++++++++++++++++++++++++++++ 2 files changed, 164 insertions(+) create mode 100644 capitolo3.tex diff --git a/CalcoloScientifico.tex b/CalcoloScientifico.tex index 7dc7b82..4cb3426 100644 --- a/CalcoloScientifico.tex +++ b/CalcoloScientifico.tex @@ -115,6 +115,7 @@ % \newcommand{\deter}[1]{\mathrm{det}(#1)} \newcommand{\rk}[1]{\mathrm{rk}(#1)} +\newcommand{\diag}[1]{\mathrm{diag}(#1)} %% %% Il comando per avere l'insieme delle parti diff --git a/capitolo3.tex b/capitolo3.tex new file mode 100644 index 0000000..909860d --- /dev/null +++ b/capitolo3.tex @@ -0,0 +1,163 @@ +\chapter{Decomposizione in valori singolari} + +In questa sezione ci occuperemo di un tipo di fattorizzazione basato sulle proprietà spettrali, e quindi differente +dalle fattorizzazioni viste fino ad ora (come la QR o la LU). Cominceremo esponendo un caso concreto in cui +ci si trova nella necessità di utilizzarla. + +\section{Problemi lineari ai minimi quadrati} +Consideriamo il seguente problema; sia $A$ una matrice $m \times n$ con $m \geq n$ e $b$ un vettore +di $\R^m$. Vogliamo trovare $x$ tale che $Ax = b$. In generale questo problema è risolubile solo se +$b \in \imm{A}$ vista come applicazione lineare, ma non è conveniente affrontare il problema risolvendo +effettivamente il sistema lineare. Possiamo quindi riformulare la nostra richiesta in questo modo: cerchiamo +$x$ tale che +\[ + x = \min_{x \in \R^n}( || Ax - b ||_2 ) +\] +Se $x$ è effettivamente soluzione del sistema lineare avremo $||Ax - b||_2 = 0$, ed in caso contrario +avremo un $x$ tale che $Ax$ sia una buona\footnote{dove buona significa che non si potrebbe trovare +un altro $x'$ tale ceh $Ax'$ sia più vicino (nel senso della norma 2) a $b$ di quanto lo sia già $Ax$} +approssimazione di $b$. + +Analizziamo alcune tecniche per la risoluzione di questo problema + +\subsection{Fattorizzazione QR} +Consideriamo la fattorizzazione QR di $A = QR$. Osserviamo che in questo caso $Q \in \matr{\R}{n}{n}$ e +$R \in \matr{\R}{m}{n}$. Osservando che minimizzare la norma 2 è la stessa cosa che minimizzarne il +quadrato si ottiene +\[ + ||Ax - b||_2^2 = ||QRx - b||_2^2 = ||Q(Rx - \trasp Qb)||_2^2 = ||Rx - \trasp Qb||_2^2 +\] +dove l'ultimo passaggio è giustificato dal fatto che $Q$ è unitaria. +Consideriamo che la matrice $R$ e il vettore $\trasp Qb$ avranno una struttura di questo tipo: +\[ + R = \left[ \begin{array}{c} + \hat R \\ + \trasp 0 + \end{array} \right] +\qquad \trasp Qb = \left[ \begin{array}{c} + w_1 \\ + w_2 + \end{array} \right] +\] +dove $\hat R$ è una matrice quadrata $n \times n$, $w_1$ un vettore con $n$ componenti e $w_2$ un vettore +di $\R^{m-n}$. Possiamo quindi riscrivere la relazione precedente in questo modo +\[ + ||Ax - b||_2^2 = ||\hat Rx - w_1||_2^2 + ||w_2||^2 +\] +Abbiamo ricondotto il problema a minimizzare la norma di $\hat Rx - w_1$. Se $\hat R$ è invertibile allora +la soluzione al problema è evidentemente unica ed è la soluzione del sistema lineare $\hat Rx = w_1$. +Osserviamo che $\hat R$ è invertibile se e solo se lo è $\rk{A} = n$; per ora assumeremo che +questo sia vero, riservandoci di analizzare il caso più generale in seguito. + +In conclusione la fattorizzazione QR ci ha permesso di ricondurre il problema alla risoluzione di un sistema +lineare quadrato di ordine $n$. Analizziamo ora un diverso approccio. + +\subsection{Sistema delle equazioni normali} +Questo sistema deriva dalla seguente osservazione. Supponiamo di avere $x$ soluzione del sistema $Ax = b$. Allora +$x$ deve essere anche soluzione del sistema quadrato $\trasp AAx = \trasp Ab$, che si ottiene semplicemente moltiplicando +a sinistra la precedente relazione per $\trasp A$. +Osserviamo ancora una volta che scrivendo la fattorizzazione QR di $A$ e supponendo che $A$ abbia rango massimo si +ottiene lo stesso risultato di prima: +\[ + A = QR \qquad \trasp AAx = b \iff \trasp RRx = \trasp R\trasp Qb \iff \trasp R\hat R = \trasp Rw_1 \iff + \hat Rx = w_1 +\] +e quindi si ha l'unicità della soluzione e la risolubilità del sistema lineare. In generale però, il calcolo +di $\trasp AA$ potrebbe esserer oneroso e, soprattutto, nessuno ci garantisce in una situazione generale +che il rango di $A$ sia massimo. + +\subsection{La decomposizione in valori singolari} +Possiamo pensare di voler migliorare l'idea delle equazioni normali. In altre parole, vogliamo trovare il +modo di evitare il calcolo di $\trasp AA$. Consideriamo dunque il seguente +\begin{te}[Decomposizione in valori singolari] +Sia $A \in \matr{\R}{m}{n}$. Esistono $U \in \matr{\R}{m}{n}$ e $V \in \mat{\R}{n}$ ortogonali e +$\Sigma \in \matr{\R}{m}{n}$ tali che +\[ + \Sigma = \left[ \begin{array}{cccc} + \sigma_1 & & & \\ + & \sigma_2 & & \\ + & & \ddots & \\ + & & & \sigma_n \\ \hline + \multicolumn{4}{c}{\trasp 0} \\ + \end{array} \right]; \qquad \textrm{dove} \: \sigma_i \geq \sigma_j \ \forall i < j \qquad \textrm{e} \qquad +A = U\Sigma\trasp V +\] +\end{te} + +\begin{os} + Osservando che $\Sigma$ è quasi diagonale si capisce che in un certo senso questa è una decomposizione spettrale. + In particolare si può scrivere $\trasp AA$ come + \[ + \trasp AA = \trasp V\Sigma\trasp U U \Sigma V = \trasp V \left[ \begin{array}{cccc} + \sigma_1^2 & & & \\ + & \sigma_2^2 & & \\ + & & \ddots & \\ + & & & \sigma_n^2 \\ \hline + \multicolumn{4}{c}{\trasp 0} \\ + \end{array} \right] V + \] +e quindi se prendiamo $\lambda_i$ gli autovalori di $\trasp AA$ (che sono reali perché la matrice $\trasp AA$ +è simmetrica) ordinati in modo decrescente abbiamo che $\sigma_i = \sqrt{\lambda_i}$ per ogni $i = 1 \ldots n$. +\end{os} +Inoltre questa fattorizzazione mette in luce altre proprietà della matrice $A$. In particolare: +\begin{os} + Se $A$ ha rango $k$, allora $\sigma_1 \ldots \sigma_k$ sono tutti diversi da $0$ mentre $\sigma_{k+1} \ldots +\sigma_n$ sono $0$. Questo si nota dal fatto che $\trasp AA$ ha lo stesso rango di $A$ ed il suo rango + è pari al numero\footnote{inteso come numero moltiplicato per la molteplicità algebrica} di autovalori diversi da $0$. + L'ordinamento che abbiamo richiesto sui $\sigma_i$ ci permette di concludere. +\end{os} + +Possiamo ora osservare come cambia il problema lineare ai minimi quadrati quando introduciamo questa +nuova fattorizzazione. +\[ + ||Ax - b||_2^2 = ||U\Sigma\trasp Vx - b||_2^2 = ||\Sigma\trasp Vx - \trasp Ub||_2^2 +\] +e ponendo $z = \trasp V x$ e $\trasp Ub = \left[ \begin{array}{c} w_1 \\ w_2 \end{array} \right] $ si ottiene +\[ + \Bigg|\Bigg|\Sigma z - \left[ \begin{array}{c} w_1 \\ w_2 \end{array} \right] \Bigg|\Bigg|_2^2 = + \Bigg|\Bigg|\left[ \begin{array}{ccc} \sigma_1 & & \\ & \ddots & \\ & & \sigma_n \end{array} \right] z - w_1\Bigg|\Bigg|_2^2 + ||w_2||_2^2 +\] +Se il rango di $A$ è $n$ allora la matrice $\diag{\sigma_1 \ldots \sigma_n}$ è invertibile e quindi per +minimizzare il primo fattore è sufficiente risolvere un sistema diagonale di ordine $n$. + +Se però $A$ è di ordine $k < n$ allora si può ripartizionare $\trasp Ub$ in una parte con $k$ componenti +e una parte con le restanti $n-k$. Riscrivendo l'espressione si ottiene che si può minimizzare il primo fattore +risolvendo un sistema di $k$ equazioni e scegliendo arbitrariamente le restanti $n-k$. +\begin{os} Come regola generale si preferisce scegliere queste ultime componenti tutte nulle in modo + da minimizzare la norma della soluzione. \`E infatti generalmente preferibile ottenere una soluzione di norma piccola. +\end{os} +Dalle relazioni sopra possiamo anche ottenere una espressione esplicita della soluzione $z$. In particolare +\[ + z_j = \frac{w_{1,j}}{\sigma_j} = \frac{\trasp U_j b }{\sigma_j} \quad \text{dove} \ U_j \ \text{è la $j$-esima +colonna di} \ U +\] +Ora possiamo ricordare che $x = Vz$ e quindi se $A$ è di rango massimo +\[ + x = Vz = \sum_{k=1}^{n} \frac{V_kU_kb}{\sigma_k} = \underbrace{\left( \sum_{k=1}^{n} \frac{V_kU_k}{\sigma_k} \right)}_{A^+} b +\] +e chiamiamo $A^+$ \emph{pseudo-inversa} di $A$. Più precisamente +\[ + A^+ = V \Sigma^{-1} \trasp U +\] +Si verifica subito infatti che $AA^+ = I$ e $A^+A = I$\footnote{ma attenzione! Non è vero che $AA^+ = A^+A$ perchè +queste sono addirittura due matrici di ordine diverso. Il loro prodotto dà sempre l'identità, ma non dello stesso +spazio.}. +Questa pseudo-inversa di chiama anche \emph{inversa di Moore-Penrose}. + +Se $A$ non è di rango massimo possiamo ugualmente definire la pseudo-inversa ``fermando'' la somma +a $k$ dove $k$ è il rango di $A$. Abbiamo dunque +\[ + A^+ = \sum_{j=1}^{k} \frac{V_jU_j}{\sigma_j} +\] +e ancora una volta vale l'eguaglianza $x = A^+b$ dove $x$ è il vettore di rango minimo che risolve (nel modo migliore) + il problema lineare ai minimi quadrati. + +\subsection{Interpolazione polinomiale} +Abbiamo visto nel corso di analisi numerica che dati $(x_i,y_i)$ con $i = 0 \ldots n$ è sempre possibile +trovare un polinomio di grado al più $n$ tale che $p(x_i) = y_i$ per tutti gli $i$. + +Nelle applicazioni concrete si ha però a che fare con grandi moli di dati sperimentali e raramente +il grado del polinomio che si vorrebbe usare per creare un modello matematico ha lo stesso ordine +di grandezza. Ci chiediamo: ``\`E allora possibile trovare il polinomio di grado $m$ << $n$ che meglio +approssimi la distribuzinoe dei nostri dati sperimentali?''; la risposta è sì e si tratta di risolvere +un problema lineare ai minimi quadrati con una matrice $m \times n$. -- 2.1.4