python1

Ciao a tutti, riprendo a scrivere di test statistici per capire se diverse serie di dati sono o meno uguali, usando in particolare il test ANOVA. Come detto la volta passata il test si basa su una serie di calcoli che solitamente sono comodamente svolti da funzioni già predefinite in vari strumenti software: la volta scorsa ho parlato di excel, ora mi concentro sul Python.

In realtà il Python non è uno strumento software per analisi statistiche ma è un completo linguaggio di programmazione. Come tale può fare molto più che analisi statistiche e per questo motivo approfitto di questo post per far conoscere (solo) un po’ questo linguaggio che è un formidabile strumento per tutti gli studenti che si avvicinano alla programmazione, alla fisica, alla matematica o più in generale all’analisi dei dati.

I test statistici sono quindi un scusa per parlare di Python!

Perchè Python

Come ho appena scritto il Python è un linguaggio di programmazione: è un linguaggio ad oggetti, interpretato (cioè il programma viene convertito in codice binario dall’interprete ogni volta che viene eseguito il codice, diversamente dai programmi compilati che vengono convertiti in codice binario una volta sola), multipiattaforma (cioè le stesse istruzioni senza nessuna modifica possono essere eseguite su tanti sistemi operativi diversi) ed è ampiamente supportato in rete. Forse fra tutte le sue caratteristiche quest’ultima è quella che mi sembra più interessante. Infatti è possibile trovare moltissime librerie per Python. L’ampia diffusione del linguaggio oltre a fornire molte librerie permette di trovare in rete moltissimi esempi, tutorial e suggerimenti per i possibili problemi che si incontrano.

Python inoltre è distribuito gratuitamente (osservazione ovvia ma la ribadisco). Ciò rende Python uno strumento potente come Matlab ma gratuito e totalmente supportato in rete. Il sito principale è https://www.python.org/downloads/, una valida alternativa, che include qualche libreria in più dello standard è https://www.continuum.io/downloads.

In fine uno strumento, distribuito con Python, utilissimo nelle analisi o nella programmazione è la shell IPython: è un terminale che permette di eseguire a linea di comando le istruzioni di Python. Ciò permette di lavorare in modo interattivo con il programma potendo fare analisi in tempo reale, un banalissimo esempio: se in C dovessi disegnare la retta y=12*x+5 dovrei scrivere un programmino con il suo main, inizializzare l’ambiente grafico  dare un range di variazione ad x calcolare y, passarle al comando grafico, compilare e ed eseguire. Con IPython scrivi istruzioni simili ma le esegui volta per volta verificando l’out del computer. Ottenere il grafico della curva è così estremamente più comodo e veloce. E così, se dovessi disegnare una seconda curva, in C devo modificare il main compilare ed eseguire, con IPython basta scrivere la nuova funzione e disegnare. Questo è un ulteriore strumento che rende Python estremamente simile a Matlab ed anche comodo nell’analisi dei dati in tempo reale oltre che nella programmazione.

Qualche Riferimento

I tutorial per iniziare sono moltissimi, a diversi livelli e lingue; giro questi solo perchè li ho usati, ma di certo ce ne sono di migliori, basta cercarli in rete… avrete l’imbarazzo della scelta:

Aggiungo https://ipython.org/, link utilissimo per trovare esempi di fisica, analisi dei segnali, matematica e statistica con una descrizione delle librerie necessarie.

Partenza rapida con Python

Volendo usare Python come linguaggio di scripting per semplici funzioni statistiche possiamo saltare le caratteristiche salienti del linguaggio (strutture dati, oggetti, moduli). Certo così si perde tanto di questo strumento e per questo vi invito ad approfondire con i tutorial e tutti gli altri strumenti disponibili. Supponiamo di usare IPython per provare tutti i comandi descritti sotto. Partiamo quindi dalle librerie.

Per iniziare a lavorare con le funzioni matematiche dobbiamo importare prima di tutto la libreria Numpy. Il comando è import numpy as np. Due riferimenti: https://docs.scipy.org/doc/numpy-dev/user/quickstart.html oppure http://www.engr.ucsb.edu/~shell/che210d/numpy.pdf

Importata la libreria Numpy abbiamo a disposizione molte funzioni matematiche e l’oggetto array (vettore). Numpy infatti permette di lavorare facilmente con i vettori. Per avere un’idea di cosa abbiamo a disposizione basta fare help(np) o dir(np): nel primo caso python elencherà metodi e variabili disponibili nella libreria e modi di uso sintetici, il secondo fornirà solo un elenco di metodi e variabili.

Per creare un array il modo più semplice è definirlo: a=np.array([1,43,28,18,13,8,6,5]). Questo comando ha creato un vettore, chiamato a, con 8 elementi che sono i numeri interi indicati in parentesi quadra. Se digito in IPython a ottengo l’array indicato. Se voglio lavorare con un singolo elemento del vettore lo indirizzo con la sua posizione: a[0]=1 a[3]=28. Posso anche indicare un intervallo di elementi: per approfondire vedere slicing dei vettori.

Altri metodi per creare dei vettori, più rapidi del precedente, possono fare ricorso alla funzione arange(start,stop,incr,dtype=..) che fornisce l’intervallo di numeri compreso fra start, stop, al passo incr. del formato dtype: a=np.arange(-2,3,0.8,dtype=float) fornisce l’array([-2,-1.2,-0.4,0.4,1.2,2,2.8]). I vettori di numpy hanno diverse proprietà, una molto utile è la possibilità di filtrare gli elementi di un vettore che rispondono come true ad un criterio, esempio: a[a>1] fornisce l’array array([1.2,2,2.8]) perché solo gli ultimi tre elementi di a rispondono come true al criterio a>1. Su questi vettori sono inoltre definite diverse operazioni statistiche np.median(a), np.average(a).

Un modo per avere vettori con molti elementi distribuiti non solo a passo costante è la funzione random (np.random()). Con random possiamo ottenere array di numeri casuali distribuiti secondo diverse densità di probabilità: a=np.random.normal(media,sigma,campioni), a=np.random.lognormal(media,sigma,campioni), a=np.poison(lamda).

Se volessimo fare i grafici di queste distribuzioni dovremmo accedere ad altre due librerie (basterebbe la prima ma usarle entrambe è più comodo): import matplotlib.pyplot as plt e import seaborn as sns (per la libreria matplotlib modulo pyplot: http://matplotlib.org/users/pyplot_tutorial.html per seaborn: https://stanford.edu/~mwaskom/software/seaborn/tutorial.html). Importate le librerie basta scrivere i comandi sotto per ottenere una prima vista delle distribuzioni:

import numpy as np
import matplotlib.pyplot as plt
import seaborn as sns
x=np.random.normal(3,0.1,2000)
y=np.random.normal(3.3,0.3,2000)
sns.distplot(x,hist=False)
sns.distplot(y,hist=False)
sns.plt.show()

distribuzioni

Sia chiaro che i grafici si possono arricchire di particolari incrementando il numero di impostazioni. Se necessario si può definire un modulo privato a cui passo solo i vettori da disegnare e che esegua automaticamente tutte le configurazioni scelte prima.

Fra i metodi messi a disposizione da numpy non ci sono test statistici. Per ottenere funzioni più complesse occorre usare la libreria Scipy (sempre import scipy as sc). Scipy ha un modulo chiamato stats che mette a disposizione diverse funzioni statistiche (oltre i tutorial un semplice help(sc.stats) è sempre utile). Ormai il gioco è chiaro, basta richiamare la funzione sc.stats.f_oneway(v1,v2,v3,…) all’interno del modulo per sc.stats per eseguire il test anova fra i vettori V1,V2 ,V3 etc. L’out di sc.stats.f_oneway(v1,v2,v3,…) sono due numeri, il primo indica il valore del test F, il secondo il pvalue. Come scritto nel primo post un pvalue basso (inferiore convenzionalmente a 0.05, ma è appunto convenzionale) indica che l’ipotesi nulla (le medie delle distribuzioni messe a confronto sono uguali) va rifiutata (almeno una distribuzione è diversa). Questo risultato è appunto una stima soggetta ad un errore che è indicato dal pvalue (cioè ho il 5% di probabilità che in realtà le distribuzioni siano uguali).

Importare dati in Python

Affichè le funzioni precedenti abbiano un un risvolto pratico occorre chiaramente saper importare dei dati reali (esperimenti, test di produzione, andamento di titoli, dati di vendita e tanto altro).  Supponiamo di avere un file csv nel formato sotto, dove A,B,C indicano il “nome” di tre serie di dati (sono gli stessi dati del mio precedente post dove i valori indicavano le misure dei diametri dei bulloni prodotti dai tool A B e C).

A 54.6 45.7 56.7 37.7 48.3
B 53.4 57.5 54.3 52.3 64.5
C 56.7 44.7 56.5 50.6 49.5

Voglio leggere il file usando le sole librerie fino ad ora usate (numpy); numpy mette a disposizione diversi metodi per leggere file (riferimento https://docs.scipy.org/doc/numpy/reference/routines.io.html#text-files): essendo il file con un formato misto di stringhe e numeri usiamo genfromtxt per cui il comando per leggere il dato diventa:

import numpy as np
data=np.genfromtxt('path',delimiter=';')

A questo punto ho importato nella variabile data il file in oggetto, ciò che ottengo è (mi scuserete ho scritto dt invece che data):

array_data

datat[0,:] mi da la prima riga, data[:,0] la prima colonna. Se invece voglio accedere ad un singolo elemento di data uso data[1,:][2], in questo caso l’out di python sarà 57.5 (come l’array di python parte da 0). La presenza dell’elemento nan (Not a Number) è dovuta alla presenza delle stringhe iniziali nel file csv che non sono ben gestite dalla funzione genfromtxt. Un vettore così fatto non permetterebbe di utilizzare nessuno dei metodi disponibili, in quanto un elemento del vettore non è definito. Si può ovviare a questo problema eliminando dal vettore l’elemento nan. Un metodo generale si affida alla funzione di numpy isfinite(array) (come si intuisce dal nome della stessa funzione, essa fornisce true se l’elemento del vettore non è nan altrimenti fornisce False). Quindi visto che in numpy posso selezionare gli elementi di un vettore che sono True secondo un certo criterio va da se che il comando dt[0,:][np.isfinite(dt[0,:])] seleziona solo gli elementi non nan del vettore:

array_data_2

A questo punto applico semplicemente il test anova sui vettori ripuliti per ottenere, il dato riportato sotto (le distribuzioni sono statisticamente uguali), come si intuisce dai boxplot riportati ottenuti graficando i vettori v0,v1,v2 con la funzione di seaborn boxplot:

Anova

boxplot

La strada seguita non è la migliore per importare dati csv ed elaborarli, l’ho seguita per mostrare qualche funzione in più di python. Come dicevo inizialmente anche per questa esigenza (importare comodamente dati csv ed elaborarli) esiste una libreria in python: Pandas (riferimenti qui e qui anche se probabilmente, dovendo lavorare intensamente in questo campo, conviene investire su qualche buon testo).

E poi arriva Pandas

Semplicemente riporto la sequenza:

import pandas as pd
dt=pd.read_csv('path','delimiter')
dt.head()

Tutto il file adesso è caricato in dt. Per vedere che formato ha, il file appena letto, posso fare dt.head() ottenendo l’intestazione del file. Dal header del file posso selezionare una colonna indicandone il nome ad esempio dt[‘Tool’], posso sapere quante volte si ripete un singolo valore della colonna usando dt[‘Tool’].value_counts(). Infine posso facilmente filtrare la tabella ad esempio se scrivo dt[dt[‘Tool’]==’A’] seleziono le righe della tabella il cui valore nella colonna Tool vale A. Detto questo posso facilmente costruire i tre vettori che contengono le misure dei tre tool in modo da poter confrontare i tre vettori. Per ottenere il primo vettore scrivo v0=dt[‘Dati’][dt[‘Tool’]==’A’], per il secondo v1=dt[‘Dati’][dt[‘Tool’]==’B’] e cosi via. Di seguito un’immagine di quanto detto fino ad ora:

dt

dt2

dt3

Cosi ottenuti i tre vettori v0, v1, v2 li posso passare ad stats.f_oneway(v0,v1,v2) per ottenere un F value di 2.16 ed un pvalue=0.16.

Avere librerie adeguate per eseguire certi lavori è fondamentale infatti il lavoro con pandas è stato molto più semplice, ma comunque avremmo superato i problemi anche senza la libreria adeguata (la comodità di un linguaggio di programmazione). In fine una curiosità: pandas non è stata sviluppata in ambito accademico, ma presso una banca d’affari per semplificare la lettura e gestione dei dati finanziari, questo solo per sottolineare ancora una volta la flessibilità e versatilità di python .

Spero che in questa carrellata su python vi sia venuta un pò di curiosità, se invece non ci sono riuscito, date ancora una chance a Python leggendo altro materiale in rete, cercate ancora… qualcuno meglio di me vi mostrerà la sua utilità!

CC BY-NC-SA 4.0
This work is licensed under a Creative Commons Attribution-NonCommercial-ShareAlike 4.0 International License.