Funzione `polyfit`

Octave permette di effettura interpolazione con polinomi usando:

p = polyfit(x0, y0, n)

Dove:

x0 è un vettore con i valori $x_i$ dei punti da interpolare
y0 è un vettore con i valori $y_i$ dei punti da interpolare

$n$ è il grado del polinomio

Se $n = m-1$ l'interpolazione è esatta
Se $n < m-1$ viene risolto il problema ai minimi quadrati

La funzione restituisce p, che contiene i coefficienti del polinomio

Interpolazione Risolvendo un Sistema Lineare

Si può anche ottenere un polinomio interpolante/approssimante con:

p = V \ y

Dove:

V la matrice di Vandermonde, calcolata per i valori $x_i$
y un vettore colonna con i valori $y_i$

Il risultato dipende dalle dimensioni di V:

Se V è quadrata si ottengono i coefficienti del polinomio interpolante
Se V ha più righe che colonne, viene risolto il problema ai min. quad.

Matrice di Vandermonde in Octave

La matrice di Vandermonde può essere ottenuta facilmente con:

V = vander(x0, n)

Dove:

x0 è il vettore di valori $x_i$ per cui calcolare la matrice
n è il numero di colonne da utilizzare
Se n viene omesso, allora la matrice è quadrata

Elementi di Informatica e
Applicazioni Numeriche T

Esercizio 1

Si vuole determinare il profilo di velocità di un fluido in condotta:

La condotta ha un diamtro di $15$ cm
Sono state effettuate $4$ misurazioni della velocità

Si determini il polinomio interpolante utilizzando:

Prima la funzione polyfit di Octave...
... Poi risolvendo il sistema lineare $V \beta = y$

Si visualizzi il polinomio interpolante utilizzando la funzione:

function plot_approx(p, x_l, x_u, x0, y0)

La funzione è disponibile nello start-kit sul sito del corso

Elementi di Informatica e
Applicazioni Numeriche T

Esercizio 2

Il file ch11_es2.csv (nello start-kit) contiene una serie di misurazioni:

La prima colonna nel file contiene i valori $x_i$
La seconda colonna contiene i valori $y_i$

Si leggano i dati dal file (utilizzando dlmread)

Si ottenga poi un polinomio approssimante

Si utilizzi dapprima la funzione polyfit...

Si modifichi l'utilizzo di polyfit rispetto all'esercizio 1...
...in modo che la funzione risolva il problema ai minimi quadrati

Esercizio 2

Poi si ripeta l'operazione risolvendo il sistema $V^TVx = V^Ty$

Si ricordi in Octave è facile utilizzare la matrice pseudo-inversa!
- Basta usare una divisione sinistra
Si ottenga la matrice di Vandermonde...
...nel formato opportuno a risolvere il problema ai minimi quadrati

Si utilizzi la funzione plot_approx per visualizzare i risultati

Opzionalmente:

Si calcoli il valore dell'errore $e = \sum_{i = 1}^n r^2$

Elementi di Informatica e
Applicazioni Numeriche T

Esercizio 3

Il file ch11_es3.csv (nello start-kit) contiene una serie di misurazioni:

La prima colonna nel file contiene i valori $x_i$
La seconda colonna contiene i valori $y_i$

Si ottenga una funzione (non polinomiale) approssimante

Si utilizzi il metodo dei minimi quadrati

Occorrerà costruire la matrice $\Phi$ ...
...con il valore dei termini della funzione approssimante

$\Phi_{i,j} = \phi_j(x_i)$

Esercizio 3

A tal file, di definisca la funzione:

function phi = build_phi_matrix(x)

Che deve costruire la matrice $\Phi$

La si utilizzi poi per risolvere il sistema:

$\Phi^T\Phi \beta = \Phi^T y$

Ed ottenere il valore dei parametri $\beta$
Si ricordi per utilizzare la matrice pseudo-inversa...
...È sufficiente usare la divisione sinistra

Esercizio 3

Si noti che è possibile valutare la funzione approssimante:

Calcolando $\Phi$ per i valori di $x$ a cui si è interessati
Calcolando $\Phi \beta$ , una volta che $\beta$ è noto

Questa tecnica è utilizzate nella funzione:

function plot_approx(p, x_l, x_u, x0, y0)

Che può essere utilizzata per valutare la qualità dell'approssimazione

Elementi di Informatica e
Applicazioni Numeriche T

Esercizio 4

Il file ch11_es4.csv (nello start-kit) contiene una serie di misurazioni:

La prima colonna nel file contiene i valori $x_i$
La seconda colonna contiene i valori $y_i$

Si ottenga una funzione (non polinomiale) approssimante

Si utilizzi il metodo dei minimi quadrati

Si proceda utilizzando la matrice $\Phi$ ...
...come nel caso dell'esercizio precendente

Si visualizzino i risultati utilizzando plot_approx

Elementi di Informatica e
Applicazioni Numeriche T

Esercizio 5

Si calcoli il polinomio interpolante per i dati dell'esercizio 1

Si utilizzi il polinomio in forma di Lagrange

In particolare, si implementi la funzione:

function y = f(x, x0, y0)

La funzione deve calcolare il polinomio in forma di Lagrange...
Per i punti nel vettore x
I parametri x0 e y0 sono necessari al calcolo

Linear Least Squares Minimizziamo l'errore Tipicamente: minimizziamo lo scarto quadratico medio Funzione convessa - definizione - minimo locale = minimo globale Conseguenza: possiamo trovare il minimo annullando la derivata Derivazione (non banalissima) - Otteniamo un sistema lineare - dimostrazione di unicità della soluzione basata sulla convessità - casi degeneri: X^TX singolare Esempio - Fitting di una parabola Estendibilità ad altre funzioni non lineari - Formulazione con funzioni generiche - Far vedee che non cambia nulla Il metodo gestisce facilmente anche misurazioni multiple per lo stesso punto Interpolazione e Convergenza Terzo problema: caratterizzazione di un fenomeno fisico - Errori di misura => ad una stessa x possono corrispondere y diverse! - Dubbio: che sia meglio introdurre questa cosa quando parliamo di statistica? Magari è meglio farlo adesso: si può sempre riprendere dopo In generale Vale la pena di parlare anche dell'interpolazione lineare?

Scelta del valore di x iniziale in NR Funzioni anonime Plotting? Linspace figure salvare figure Piano (Vecchia lezione): Tipi di convergenza: [SPOSTARE DOPO IL CASO NON-LINEARE] - Lineare, sublineare, superlineare Differenze con il caso lineare: - Bacini di attrazione (sono un casino: possono essere discontinui...) Un caso particolare: radici di polinomi - Per adesso focalizziamoci sul trovare una sola radice - Primo approccio: riduzione al caso non lineare + Problema: condizioni di convergenza - Secondo approccio: + Usiamo la matrice compagna + Se troviamo gli autovalori, troviamo le radici + Come è definito un autovalore? + Facile, ma c'è lambda tra le scatole + Eliminiamolo + Rinunciamo a trovare l'autovettore: ci accontentiamo di un suo multiplo (che ha lo stesso comportamento) + Alé! Possiamo applicare l'iterazione di punto fisso + Alla fine troviamo un vettore unitario con la stessa direzione di un autovettore + Come trovare l'autovalore? Rapporto tra i coefficienti di due vettori successivi (qualunque coefficiente, basta che sia non nullo) + Osservazione: riduzione ad un caso non lineare + Converge? Condizioni + Produce sempre l'autovalore con valore assoluto più grande - Trovare tutte le radici: regola di Ruffini (divisione sintetica) Fix-point iteration Fix-point iteration - Calcolo della radice quadrata? Attrattori? - Punti fissi attrattivi - Bacino di attrazione Teorema del punto fisso di Banach? Power iteration/power method - Condizioni di convergenza? - Dim. di convergenza all'autovalore maggiore? - Non è che a loro interessi tanto... Bracketing methods Teorema del punto intermedio Ricerca lineare Metodo di bisezione Regula falsi NOTA: non sono riducibili a fix-point iteration Non-bracketing methods Metodo delle secanti Metodo di Newton Radici di un polinomio? - Metodo di Horner + Newton - Companion matrix + power method? Errori? Calcolo della derivata Cancellazioni? Evitate nella regula falsi Remeber to explain: Multi variate functions applied to vectors subfunctions array concatenation memorizzazione delle matrici in memoria continue (magari nella lezione sull'ottimizzazione) trattamento degli argomenti non specificati return Gestione errori (error/warning) Debugging Prossime volte: determinant computation (recursive!) factorial computation (recursive!) matrix argument to vector operation - a good occasion to explain multiple functions in a single file sieve of erathostenes other algorithms to find prime numbers? permutations - Permute vectors - Permute rows/columns of a matrix Taylor series expansion? convolution? bubble sort / naive sort comparison functions set operations - intersection - union - substract? - powerset? continued fraction expansion Coordinate conversion (cartesian, polar) Some matrix factorizations/decompositions? - Choleski - LU decomposition - QR factorization - Singular Value Decomposition Gradient descent Newton method Fixed point iteration Bisection method - bracket a zero - convex optimization Minimization by repeated evaluation Golden section search Numeric gradient/Jacobian/Hessian Nelder-Mead? RNG algorithms (something simple) Using diagonal matrices for scaling Integration methods: - Numerical integration based on Gaussian quadrature. - Numerical integration using an adaptive vectorized Simpson’s rule. - Numerical integration using an adaptive Lobatto rule. - Numerical integration using an adaptive Gauss-Konrod rule. - Numerical integration using adaptive Clenshaw-Curtis rules. - Numerical integration of data using the trapezoidal method. Multi-dimensional integration? - Recursive appraoch - Orthogonal collocation ODE - algo? - reduction to first order Optimization? - LP, QP (null-space active-set method)... the algos are too complex - QP: basic method for convex problems? Zero-gradient - NLP via sequential quadratic programming - Least Square method? Descriptive statistics - Mean - Median - Variance - Stdev (second moment, statistic bias -- maybe this one later) - Mode? - Quartiles, quantiles - Explain histograms, too? - Histogram count computation - correlation/covariance? "True" statistics & probability theory - stat plots - unbiased estimators - some probability distributions - central limit theorem - Distributions - pdf, cdf, quantiles - One or two stat tests? - Confidence intervals? Generate permutations? Centering / Normalization Polynomials - Representation? - Horner's method? - Roots of a polynomial: Frobenius companion matrix, eigenavalues - poles? Interpolation - Nope - Use variable transformation for power law fitting Verifica stato di un elemento di un array Sort by argument Calcolo della radice quadrata Serie di Taylor (per calcolare una funzione) I/O semplice + disp + printf? Per stampare con precisione alta Goal: algoritmi + complessi ed esempi di propagazione degli errori valutazione di un polinomio (trivial) valutazione di un polinomio (Horner) Somme: - Somma "normale" - Somma gerarchica - Somma di Kahan Numeri interi - identificazione di numeri primi - fattorizzazione in numeri primi (trial division) Algoritmi di ordinamento: - Solo naive sort e bubble sort RNG - un algoritmo semplice - esempi di utilizzo di RNG in Octave Valutazione di una funzione a gradini - Lookup in una tabella chiave-valore - Valutazione in punti multipli Numeri interi ed algoritmi - Interi in Octaves - Promozione/demozione? - sieve of erathostenes (magari nella lezione 3) - fattorizzazione in numeri primi (trial division) - GCD GCD: q = gcd a = n_a q // q è il GCD di a e b b = n_b q // q è il GCD di a e b a = p b + r // divisione di a per b Quindi: n_a q = a = p b + r = p n_b q + r Quindi n_a = p n_b + r / q Deve venire un intero => r/q deve essere intero => r è divisibile per q * Si dovrebbe poter fare lo stesso ragionamento anche per la sottrazione

Elementi di Informatica e Applicazioni Numeriche T

Interpolazione in Octave

Elementi di Informatica e
Applicazioni Numeriche T