Gnuplot la procedura di “fit”

Lo scopo del nostro esperimento è quello di cercare la relazione che lega due grandezze fisiche tra di loro:

Supponiamo di aver riportato, per comodità, le coppie di grandezze misurate (corredate dal rispettivo errore assoluto) in un foglio elettronico.

Se cerco di capire quale funzione lega la grandezza Y alla grandezza X, è sicuramente una buona idea riportare in grafico Y in funzione di X. Per eseguire il grafico con gnuplot per prima cosa occorre copiare i dati dal foglio elettronico…

… e incollarli su un editor di testo (tipo “blocco note” di Windows)

L'editor di testo deve essere in grado di generare e salvare file con estensione .txt, altri tipi di file di testo non sono compresi da gnuplot

Nel caso in cui al posto del punto decimale ci sia la virgola, occorre rimpiazzarla col punto…

E’ comodo avere nel file una intestazione che ricordi quali grandezze sono riportate in ciascuna colonna, ma gnuplot non vuole righe di testo tra i dati numerici, quindi “commentiamo” la prima riga aggiungendo all’inizio il “cancelletto” #

Adesso salviamo il file con un nome semplice, ad esempio “dati.txt”

Apriamo gnuplot e fissiamo le etichette degli assi del grafico (tanto per non dimenticarcene):

gnuplot> set xlabel "i (mA)" gnuplot> set ylabel "V (V)"

Se il file di testo non si trova nella stessa cartella di gnuplot, occorre cambiare directory (oppure spostare il file nella stessa cartella…)

gnuplot> cd "Desktop" gnuplot> !ls dati.ods dati.txt

eseguiamo il grafico:

gnuplot> plot [0:1] [0:2.5] "dati.txt" u 1:3 gnuplot> set grid

X

Y

il comando !ls visualizza la lista di file nella directory selezionata, per essere sicuri che ci sia davvero il file voluto…

la scrittura u 1:3 serve a dichiarare che la colonna 1 del file è la variabile X, mentre la colonna 3 è la variabile Y

ecco il risultato:

si tratta di un grafico sperimentale, quindi occorre inserire anche gli errori di misura:

gnuplot> plot [0:1] [0:2.5] "dati.txt" u 1:3:2:4 with xyerrorbars

la scrittura 2:4 indica che la colonna 2 del file contiene gli errori sulla X e la colonna 4 quelli sulla Y

ΔX

ΔY

with xyerrorbars errori su entrambi gli assi

Si può notare che nel grafico la barra di errore sulle X si nota appena.

Se ΔX/X ⪻ ΔY/Y, l’errore sulle X risulta trascurabile in confronto a quello delle Y. Si omette quindi di rappresentare nel grafico l’errore su X

Nel proseguo supporremo sempre di poterci mettere nella condizione (da verificare sempre!), che l’errore relativo di una delle due grandezze sia trascurabile, rispetto all’errore relativo dell’altra grandezza.

Scegliere sempre come variabile indipendente quella con errore relativo trascurabile

gnuplot> plot [0:1] [0:2.5] "dati.txt" u 1:3:4 with yerrorbars

la colonna 4 contiene l’errore su Y

ΔY

with yerrorbars errore solo su asse Y

Il grafico ottenuto suggerisce che i punti sperimentali sono disposti su una retta di equazione:

y = ax + b Problema: come facciamo a determinare la migliore retta che passa per i punti sperimentali?

(xi , yi )

FIT

trovare la cur va migliore che si adatta ai dati

y = y(x)

Alla base di molti esperimenti c’è la necessità di trovare o verificare la relazione esistente tra una coppia di grandezze fisiche, al fine di determinare la legge che regola il fenomeno. Dobbiamo prima di tutto scegliere il tipo di curva (intuito, analogia, teoria, ...) Successivamente trovare un criterio quantitativo che permetta di ricavare una stima dei parametri della curva accanto alla stima dell’errore associato Nel nostro esempio la procedura di fit dovrà farci ottenere:

a ± ∆a

b ± ∆b

In questo modo otterremo un fascio di rette che approssima i nostri dati, di cui la migliore è:

y = ax + b

Un modo semplice per eseguire il fit lineare è quello di disegnare a mano sul grafico la retta che passa più vicino ai punti sperimentali, insieme alle due rette estreme, che comunque attraversano le barre di errore

Dalla lettura del grafico siamo poi in grado di risalire all’equazione matematica delle rette e quindi stimare sia i parametri a e b che gli errori associati

Metodo dei “minimi quadrati” X

Y

σ

x1

y1

σ1

x2

y2

σ2

…

…

…

xi

yi

σi

…

…

…

xN

yN

σN

Ipotesi: • errore relativo su x trascurabile rispetto a quello su y. • grandezza yi distribuita statisticamente con deviazione standard σi : yi ± σi • la relazione che lega le due grandezze è del tipo: y=f(x) • la funzione f(x) (retta, esponenziale, ecc…) è descritta da alcuni parametri che voglio ottenere dal fit

Definisco χ2 la somma:

2

χ =

N 2 ! [yi − f (xi )] i=1

σi2

Tra tutte le funzione del tipo f(x), quella che approssima meglio i dati sperimentali è quella i cui parametri rendono minimo il χ2

Senza entrare nei dettagli della teoria statistica, proviamo a giustificare la precedente affermazione, facendo riferimento al nostro esempio:

χ2 =

f (x) = ax + b

N ! [yi − (axi + b)]2 i=1

yi f (xi ) = axi + b ta

ret

}.

Mi aspetto che la retta migliore passi più vicino possibile a tutti i punti sperimentali. Quindi i parametri a e b saranno tali che, per ogni xi (con i da 1 a N), la grandezza:

σi

}. y − f (x )

t s e b

xi

i

σi2

i

[yi − (axi + b)]2 σi2 sia la più piccola possibile. Te n e n d o c o n t o d i t u t t i i p u n t i sperimentali è ragionevole cercare i parametri a e b in modo che la somma χ2 risulti minima.

Quanto piccolo mi aspetto sia il χ2 ? Se gli errori sulle y non sono eccessivamente sovrastimati, mi aspetto che:

[yi − (axi + b)]2 !1 2 σi Se gli N addendi del χ2 fossero tutti indipendenti tra loro, mi aspetterei che χ2 ∼ N Questo non è però vero perché a e b sono ricavati a partire dalle N misure iniziali, questo comporta che χ2 ∼ N-2. In generale chiamo gradi di libertà il numero di misure effettuate meno il numero dei parametri che devo ricavare dal fit. Mi aspetto che χ2 sia circa uguale al numero di gradi di libertà. In effetti, per comodità, si introduce il χ2 ridotto:

χ2red

1 = N DF

N ! [yi − (axi + b)]2 i=1

σi2

Con NDF= numero gradi di libertà = numero misure effettuate - numero parametri della funzione Un test della bontà di un fit è ritenuto il fatto che χ2red∼1

Generalmente minimizzare il in funzione dei parametri da stimare implica molti calcoli matematici. Gnuplot può farlo per noi, adottando un procedimento numerico iterativo per ottenere i parametri del fit. Vediamo come fare:

sempre indicare il nome assegnato alla funzione seguito dal nome del file tra virgolette

gnuplot> f(x)=a*x+b gnuplot> fit f(x) "dati.txt" u 1:3:4 via a,b

dichiaro le colonne del file che devono essere utilizzate per eseguire il fit: colonna 1 -> X colonna 3 -> Y colonna 4 -> ΔY

la scrittura via a,b dichiara quali sono i parametri del fit

Alla fine del fit la finestra del terminale di gnuplot visualizza il risultato ottenuto. Lo stesso risultato viene riassunto in un file di testo generato dal programma col nome fit.log il programma restituisce come voluto a±Δa e b±Δb

attenzione! Il programma non scrive i risultati in modo corretto per quanto riguarda le cifre significative, inoltre mancano, come ovvio, le unità di misura

a = (2, 23 ± 0, 02) V/mA = (2, 23 ± 0, 02) · 103 V/A b = (−0, 009 ± 0, 012) V

N.B. La stima del parametro b ottenuta è consistente con l’ipotesi che la migliore retta passi per l’origine. Le due grandezze X e Y sono quindi direttamente proporzionali

Se viene richiesto di verificare la diretta proporzionalità tra due grandezze fisiche, oppure appare dal grafico (come in questo caso) che i punti sperimentali si allineano su una retta passante per l’origine, invece di eseguire il fit con la funzione f(x)=ax+b, conviene passare alla funzione f(x)=ax

gnuplot> f(x)=a*x gnuplot> fit f(x) "dati.txt" u 1:3:4 via a

a = (2, 22 ± 0, 01) · 103 V/A

adesso rimane un solo parametro: a

Riportiamo in grafico i dati sperimentali e la funzione ottenuta dal fit. Oltre ad essere richiesta per completezza, questa operazione offre anche un riscontro visivo della bontà del fit

gnuplot> plot [0:1] [0:2.5] "dati.txt" u 1:3:4 with yerrorbars, f(x)