Pagina 1 di 2

Minimi quadrati

MessaggioInviato: 28/10/2011, 19:49
da anonymous_ed8f11
Buongiorno, sono alle prese con un problema pratico dove devo utilizzare degli strumenti dell'analisi numerica, e vorrei chiedervi un consiglio su come fare.

Ho studiato in maniera analitica il rendimento di una turbina ed è emerso che la curva che ne descrive l'andamento in funzione della velocità di rotazione è una parabola passante per l'origine (a turbina frema il rendimento è nullo).
Poi in laboratorio ho raccolto dei dati sperimentali ed in effetti l'andamento rispecchia approssimativamente il modello teorico.

Vorrei trovare la curva che descrive il fenomeno e ho pensato di farlo in matlab utilizzando il metodo dei minimi quadrati. Ho delle coppie di punti $(x_n,y_n)$ nel piano cartesiano, e cerco di trovare il polinomio di secondo grado che più si avvicina al loro andamento.
Tutto ciò non è complicato e l'ho già fatto con il comando $"polyfit"(x,y,2)$, soltanto che la parabola trovata non passa per l'origine a causa di errori sperimentali.
Dunque vorrei imporre il passaggio per l'origine prima di trovare la soluzione ai minimi quadrati, ma non so come fare :?

Per caso sapete se in matlab esiste qualche funzione per questo, oppure qualche stratagemma per risolvere il problema?
Grazie in anticipo, Lorenzo

Re: Minimi quadrati

MessaggioInviato: 29/10/2011, 16:57
da dissonance
Io farei una semplice sottrazione. Proprio nel modo più brutale possibile:
Codice:
zstar=polyfit(x, y, 2);
z=zstar-zstar(1);

Re: Minimi quadrati

MessaggioInviato: 29/10/2011, 17:43
da anonymous_ed8f11
Haha, è vero, così passa di sicuro per l'origine :P
Il problema è che il ragionamento dell'interpolazione ai minimi quadrati va un po'a farsi benedire :)

Comunque in questo caso la parabola trovata col metodo dei minimi quadrati per fortuna non si allontana molto dall'origine, quindi se uso questo stratagemma nessuno se ne accorgerà! Ti ringrazio molto perchè non ci avevo proprio pensato.

Ad ogni modo, se a qualcuno viene in mente un metodo rigoroso sarei curioso di capire come fare, anche per una questione di interesse mio personale

Re: Minimi quadrati

MessaggioInviato: 29/10/2011, 20:32
da dissonance
Ma non stiamo barando, facendo così. E' una procedura che si usa, o almeno, così mi risulta (non ho mai studiato granché di matematica applicata, purtroppo). Il concetto è che c'è un "rumore di fondo" presente nelle misure e lo dobbiamo sottrarre via prima di potere usare i risultati.

Re: Minimi quadrati

MessaggioInviato: 29/10/2011, 23:49
da anonymous_ed8f11
Ho provato con il tuo metodo ma purtroppo non va bene :(
Il fatto è che per x=0 il polinomio vale -0.4, e dal momento che la curva descrive un rendimento, compreso sempre tra 0 e 1, risulta essere troppo da aggiungere e "sballa" tutti i dati.

A questo punto penso toccherà fare un piccolo script. Il generico polinomio il cui grafico passa per l'origine è "ax^2+bx+0", e quindi in qualche modo dovrei andare a ricondurmi a risolvere un sistema ai minimi quadrati dove però le incognite sono soltanto $a$ e $b$.
Domani mattina cerco i quaderni di analisi numerica e provo a buttare giù qualcosa, che ora non mi ricordo le funzioni a memoria..

Re: Minimi quadrati

MessaggioInviato: 30/10/2011, 09:29
da Deckard
Puoi usare un solver per l'ottimizzazione (es. AMPL o non so se anche matlab o mathematica possano farlo) che ti minimizzi lo scarto ai minimi quadrati:
$ min_{a,b} sum_(i) (y_i - (a*x_i^2 + b*x_i))^2 $

Così non hai bisogno di implementare algoritmi o cose del genere. Comunque in matlab ci saranno sicuramente librerie con all'interno metodi già implementati come ad es. Gauss–Newton.

Re: Minimi quadrati

MessaggioInviato: 30/10/2011, 10:23
da anonymous_ed8f11
Ho trovato che il comando per cercare il minimo di una funzione a più variabili è $"fminsearch(f,v)"$ con $f$ la funzione di cui cerco il minimo e $v$ un punto vicino al quale voglio cercare il minimo stesso.

Ho scritto la funzione utilizzando il prodotto riga colonna la posto della sommatoria:
Codice:
f=inline('(y-(a.*x.^2+b.*x).^2)*ones','a','b')
dove x de y sono vettori riga 1x8 (devo inperpolare 8 valori) e ones è un vettore colonna 8x1 con tutti gli elementi che valgono 1.

Ora però se provo a minimizzare la funzione ottengo errori :(
Codice:
>> fminsearch(f,[0,0])
??? Error using ==> inline.subsref at 14
Not enough inputs to inline function.

Error in ==> fminsearch at 205
fv(:,1) = funfcn(x,varargin{:});

Penso sia dovuto a come ho definito la funzione, perchè nell'esempio che ho trovato veniva usato il comando $"function"$, ma nel mio contesto se definivo la funzione $f$ col comando function mi venivano fuori altri errori :s

Re: Minimi quadrati

MessaggioInviato: 30/10/2011, 11:03
da dissonance
No, no, non usare quella cosa là!!! :-) Quel comando serve per minimizzare funzioni non lineari e usa algoritmi molto più incasinati del linearissimo metodo dei minimi quadrati. Guarda, secondo me fai prima a scriverti una piccola function che fa quello che vuoi tu. Risolvere il problema dei minimi quadrati è in ultima analisi risolvere un sistema di equazioni lineari, non occorrono grandi manovre.

Re: Minimi quadrati

MessaggioInviato: 30/10/2011, 11:49
da anonymous_ed8f11
Il fatto è che non so bene come fa matlab a risolvere i minimi quadrati. So che per risolvere un sistema $Ax=b$ la matrice $A$ viene divisa con la fattorizzazione QR, cioè $A=QR$ con Q quadrata e R rettangolare. Poi si calcola il residuo e lo si minimizza, ma a qual punto non so che algoritmo viene usato per minimizzarlo..

Re: Minimi quadrati

MessaggioInviato: 30/10/2011, 14:31
da Deckard
Ops, scusami ma non so perché ti ho indicato il metodo per risolvere i minimi quadrati in generale, ovviamente in un caso lineare come questo usare qualcosa come Gauss-Newton è fuori luogo come dice dissonance.
Supponiamo di avere $m$ coppie $(x_i,y_i)$; è sufficiente costruire la matrice $m times 2$ $X$ formata da elementi $X_{ij} = x_i^j$ (ovvero una matrice costituita da due colonne: la prima contente tutte le $x_i$ la seconda le $x_i$ elevate al quadrato) e il vettore delle incognite $beta = (b, a)^T$. Quello che vorresti fare è risolvere il sistema $X beta = y$ dove $y = (y_1, ... , y_m)^T$. Ovviamente tale sistema non avrà soluzione. Vuoi quindi trovare il vettore $beta$ che più si avvicina a tale soluzione (che minimizza gli scarti quadratici). Con qualche semplice passaggio algebrico si può dimostrare che tale vettore è $beta = (X^T X)^{-1} X^T y$.
Non c'è bisogno quindi di usare nessun algoritmo, costruisci la matrice $X$ e il vettore $y$ e sei a posto.