Il linguaggio Python di base (incluse le librerie standard) fornisce strumenti sufficienti per completare semplici progetti computazionali. Tuttavia, esistono librerie Python dedicate che forniscono funzionalità estese che
forniscono strumenti numerici per automatizzare operazioni ricorrenti
sono semplici da usare
sono più efficienti in termini di tempo di CPU e di esigenze di memoria degli strumenti in Python base.
Citiamo tre librerie in particolare:
Il modulo numpy, introdotto in 07_numpy, che fornisce strumenti numerici. Numpy contiene fra le altre cose molte funzioni per manipolazioni statistiche che vengono presentate in 09_numpy_statistics.
Il modulo matplotlib, introdotto in 08_plots, che permette di creare grafici.
Il modulo scipy (SCIentific PYthon), introdotto in 10_scipy, che fornisce un gran numero di algoritmi numerici.
Molti degli algoritmi numerici resi disponibili da numpy e scipy sono forniti da librerie compilate, di solida tradizione, che spesso sono scritte in Fortran o C. Vengono quindi eseguite molto più velocemente di codice scritto in puro Python (che è interpretato). Generalmente, un codice compilato è due ordini di grandezza più veloce di un codice in puro Python.
Come al solito, si può utilizzare la funzione help su ciascuno dei metodi numerici per accedere alla documentazione.
La librery NumPy (NUMerical PYthon) fornisce un nuovo tipo di struttura dei dati chiamata arrays che permette di eseguire in modo efficiente operazioni su vettori e matrici. Fornisce inoltre diverse operazioni legate all'algebra lineare (come risolvere sistemi di equazioni lineari, calcolo di autovettori e autovalori). Viene tradizionalmente importata come np.
NumPy introduce un nuovo tipo di dati detto “array”. Un array sembra molto simile a una lista ma può contenere solo elementi di un singolo tipo (mentre una lista può contenere oggetti di tipo diverso). Questo significa che gli array possono essere scritti in memoria in modo più efficiente. Gli array sono la struttura di dati migliore per i calcoli numerici in cui spesso si ha a che fare con vettori e matrici di tipo omogeneo.
Vettori, matrici bidimensionali e matrici con più di due indici sono tutti chiamati “arrays” in NumPy.
La struttura dati data che useremo più spesso è il vettore. Qualche esempio di come crearne uno:
import numpy as np
x = np.array([0, 0.5, 1, 1.5])
print(x)
type(x)
type([0, 0.5, 1, 1.5])
x
a = np.array([1+1j,2-1j])
b = np.array([1-1j,2-1j])
a*b
a+b
x = np.arange(0, 2, 0.5)
print(x)
x = np.zeros(4)
print(x)
y = np.ones(10)
print(y)
Numpy ha i metodi np.append, np.insert, np.delete che permettono di modificare il numero di elementi di un array. Tutti questi metodi per&ograsve;non agiscono in-place ma creano e ritornano un nuovo array, quindi sono computazionalmente costosi. Per questa ragione è preferibile creare un array con metodo come zeroes e ones e poi modificare i valori.
In Numpy esiste il comando linspace(x0, x1, n), che genera una lista di n elementi equidistanti fra x0 e x1 (inclusi). Il comando logspace(x0, x1, n) genera una lista di n elementi equidistanti, in scala logaritmica, fra 10\*\*x0 e 10\*\*x1; è utile per fare plot logaritmici come vedremo in 08_plots. Alcuni esempi:
np.linspace(0,10,1001)
np.linspace(-1,1,21)
Ricordate che i primi due argomenti di logspace sono i logaritmi in base 10 del punto iniziale e di quello finale. Il terzo parametro è il numero di divisioni, estremi inclusi.
np.logspace(1,4,4)
#np.logspace(1.7,4.4,40)
tipo di dati di un array¶Quando viene creato un array è possibile specificare il tipo di dati che contiene. I tipi possibili sono: 'float', 'complex', 'int', 'bool', 'str' and 'object'. Un controllo più fine si può avere usando espressioni come 'float32', 'float64', 'int8', 'int16' oppure 'int32'; l'intero specifica il numero di bytes.
array_f = np.array([1,2,3.14], dtype=float)
array_i = np.array([1.,2.2,3], dtype='int16')
array_f
array_i
np.array([1,2,123456790],dtype='int32')
np.arange(1.3,2.2,0.1,dtype=complex)
Il tipo di un array può essere modificato con il metodo astype
array_f.astype(complex)
array_f.astype('int16')
x = np.zeros(4)
x[0] = 3.4
x[2] = 4
print(x)
print(x[0])
print(x[0:-1])
lenlen(x)
x = np.arange(0, 2, 0.5)
print(x)
print(x + 10)
print(x*33)
print(x ** 2)
La maggior parte delle operazioni (+, -, *, /) agiscono su tutto un array in un solo colpo, elemento per elemento.
y1 = np.array([1.,2.,3.])
y2 = np.array([-1.,-10,+100])
y1+y2
Nel modulo numpy sono definite tutte le funzioni matematiche più comuni in modo che agiscano su tutti gli elementi di un array. Sono definite anche le costanti come np.e e np.pi. In generale è superfluo importare separatamente il modulo math.
print(np.sin(x))
print(np.e)
La funzione mod corrisponde alla operazione modulo % con molte opzioni in più (Hint: help(np.mod)).
x = np.arange(10)
print(x)
np.mod(x,6)
Il modo migliore di utilizzare una funzione in numpy è definire l'array dei punti in cui la funzione deve essere calcolata e poi applicare la funzione all'array.
Esempio: per calcolare i valori del polinomio y = 3* x**3 - 2*x**2 + 1 fra -5 e 5:
def myfunc(arr):
return 3*arr**3 -2*arr*arr +1
x = np.linspace(-5,5,101)
print(x)
res = myfunc(x)
res
check:
3*(-4.9)**3 -2*(-4.9)**2 +1.
x = np.arange(5,15)
x
x1 = [x[2],x[3],x[5]]
print(x1)
type(x1)
Per ottenere un array dal risultato precedente:
x1 = np.array(x1)
print(x1)
type(x1)
ind = [2,3,5]
x[ind] # Restituisce un array
arr = np.array([2, 3, 4, 5, 6, 7, 8, 9])
arr[arr % 2 == 1]
oppure, con più flessibilità sul modo di imporre le condizioni:
np.array([x for x in arr if x% 2 == 1])
arr = np.array([0, 1, 2, 3, 4, 5, 6, 7, 8, 9])
arr[arr % 2 == 1] = -1
arr
np.where(condition,X,Y) restituisce l'elemento di X se la condizione è vera e l'elemento di Y se la condizione è falsa.#help(np.where)
arr = np.arange(2,10)
out = np.where(arr % 2 == 1, -1, arr)
print(arr)
out
np.where(condition) restituisce una ntupla il cui primo elemento è l'array degli indici degli elementi per cui la condizione è vera.
out = np.where(arr % 2 == 1)
out
Per estrarre i corrispondenti elementi:
arr[out]
Per imporre più condizioni si può utilizzare la notazione:
arr = np.array([0, 1, 2, 3, 4, 5, 6, 7, 8, 9])
np.where( (arr>=1) & (arr<+7) )
Le parentesi tonde attorno alle condizioni sono indispensabili perchè l'operatore & ha precedenza più alta degli operatori di confronto. Le parentesi forzano Python a valutare le disuguaglianze prima di combinarle con and.
Si possono anche combinare le varie condizioni in una condizione unica con logical_and:
index = np.where(np.logical_and(a>=5, a<=10))
a[index]
uniqs, counts = np.unique(out, return_counts=True)
print("Unique items : ", uniqs)
uniqs2 = np.unique(out)
uniqs2
#help(np.unique)
a = np.array([1,2,3,2,3,4,3,4,5,6])
b = np.array([7,2,10,2,7,4,9,4,9,8])
np.intersect1d(a,b)
a tutti gli elementi che compaiono anche in b.a = np.array([1,2,3,4,5])
b = np.array([5,6,7,8,9])
np.setdiff1d(a,b)
a = np.array([1,2,3,2,3,4,3,4,5,6])
b = np.array([7,2,10,2,7,4,9,4,9,8])
np.where(a == b)
np.asarray(a == b)
np.asarray(a == b).nonzero() # restituisce gli indici degli elementi non nulli cioè di quelli non False
Python ha il metodo sort e la funzione sorted per ordinare liste. In numpy si può usare la funzione numpy.sort che restituisce una copia ordinata dell'array passato in input.
x = np.array([1,3,2])
y = np.sort(x)
print(y)
print(x)
Per ordinare un array in ordine decrescente una possibilità è:
y1 = np.sort(x)[::-1]
print(y1)
Per convertire un array in una lista o una si possono usare le funzioni python standard `list(s0)` e `tuple(s0)` che accettano una sequenza `s0` come input e ritornano rispettivamente una lista e una ntupla:
a = np.array([1, 4, 10])
a
list(a)
tuple(a)
È spesso utile misurare il tempo necessario per eseguire dei comandi. Jupyter fornisce la funzione %timeit che possiamo utilizzare per mostrare la maggiore velocità delle operazioni di Numpy su array rispetto alle operazioni di math su liste.
def forloopmethod(N):
a1 = [0]*N
for i in range(N):
a1[i] = float(i)**2 # passare a numeri reali rende l'operazione più veloce
return sum(a1)
def numpymethod(N):
a3 = np.sum(np.arange(0, N)**2, dtype='d') # dtype = 'd' vuol dire double precision
return a3
N = 1000000
%timeit forloopmethod(N)
%timeit numpymethod(N)
Ecco due modi per creare un array a due dimensioni:
x = np.array([
[1, 2, 3], [4, 5, 6]
])
x
x = np.zeros((5, 4))
x
La "forma" di una matrice può essere trovata con il comando shape (in questo caso ci sono 2 righe e 3 colonne):
x=np.array([[1, 2, 3], [4, 5, 6]])
print(x)
x.shape
La "forma" di una matrice può essere modificata con il comando/metodo reshape:
a = np.arange(6).reshape((3, 2)) # metodo
np.reshape(a, (2, 3)) # funzione
I singoli elementi possono essere recuperati con la sintassi seguente:
x=np.array([[1, 2, 3], [4, 5, 6],[7,8,9]])
x[0, 0]
x[0, 1]
x[0, 2]
x[1, 0]
Per recuperare una riga o una colonna:
x[:, 0]
x[0,:]
Per recuperare una lista di righe o colonne:
x[[0,1],:]
x[np.arange(1,3),:]
La stessa sintassi si utilizza per modificare un elemento della matrice:
x[1,0] = 157
x
A = np.zeros((3,3))
print(A)
vec = np.array([-2,5,1])
indices = np.array([0,1,2])
A[indices,indices] = vec
A
Un esempio che utilizza due array diversi per gli indici sulle righe e gli indici sulle colonne:
A = np.zeros((4,4))
ind_row = np.array([0,2])
ind_col = np.array([1,3])
vec1 = np.array([5,1])
A[ind_row,ind_col] = vec1
A
x = np.array([1, 4, 0], float)
y = np.array([2, 2, 1], float)
print("Matrices and vectors.")
print("x:")
print(x)
print("y:")
print(y)
print("Inner product of x and y:")
print(np.inner(x, y))
print("Outer product of x and y:")
print(np.outer(x, y))
print("Cross product of x and y:")
print(np.cross(x, y))
Due arrays possono essere moltiplicato nel senso usuale dell'algebra lineare usando `numpy.dot`. Ecco un esempio:
import numpy as np
import numpy.random
A = numpy.random.rand(5, 5) # genera matrice random 5 by 5
x = numpy.random.rand(5) # genera un vettore di 5 elementi
print(A)
print(x)
b=np.dot(A, x) # multiplica la matrice A per il vettore x
print(b)
Per risolvere un sistema di equazioni lineari A*x = b, dato in forma matriciale (cioè A è una matrice e x e b sono vettori, con A e b noti), possiamo usare la libreria di algebra lineare `linalg` di `numpy`:
import numpy.linalg as LA
x = LA.solve(A, b)
x
check:
print(np.dot(A,x))
print(b)
X = np.array([[1, 2], [4, 5]])
LA.det(X)
np.trace(X)
#help(np.linalg)
X = np.array([[1, 2], [4, 5]])
Xinv = LA.inv(X)
print(Xinv)
np.dot(X,Xinv)
Xt = np.transpose(X)
print(Xt)
Ricordiamo che data una matrice $M$, il vettore non nullo $v_i$ è un autovettore di $M$ e lo scalare $\lambda_i$ è il corrispondente autovalore se: $$ M\,v_i = \lambda_i \, v_i.$$ Ecco un piccolo esempio che calcola gli autovettori e autovalori di una matrice con il comando `eig`:
import numpy.linalg as LA
import numpy as np
R = np.array([[1,2,3], [1,0,1], [0,1,2]])
LA.det(R)
evalues, evectors = LA.eig(R)
print(evalues)
Controlliamo se il prodotto degli autovalori è uguale al determinante:
np.prod(evalues)
print(evectors)
Attenzione: gli autovettori sono disposti sulle colonne.
Verifica:
A1=evectors[:,0]
print(A1)
print(np.dot(R,A1))
print(evalues[0]*A1)
Ciascuno di questi comandi fornisce la propria documentazione. Per esempio, `help(LA.eig)` fornisce dettagli sulla funzione che calcola autovettori e autovalori (Avendo importato numpy.linalg come LA).
X = np.array(
[[6, 3, 7],
[2, 6, 4],
[7, 2, 5]]
)
Ordinare le colonne:
np.sort(X, axis=0)
Ordinare le righe:
np.sort(X, axis=1)
Per avere programmi veloci è necessario evitare di ciclare sugli elementi di vettori e matrici usando invece algoritmi vettorializzati. Il primo passo per convertire un algoritmo scalare in un algoritmo vettorializzato è costruire delle funzioni che accettino input vettoriali.
AVVERTENZA: è quasi sempre possibile trovare funzioni di Numpy, quindi già vettorializzate, che facciano quello che ci serve. vectorize va utilizzata come ultima risorsa quando tutte le ricerche nella documentazione e in rete si sono rivelate vane.
def Theta(x):
"""
Implemenazione scalare della funzione a gradino di Heaviside.
"""
if x >= 0:
return 1
else:
return 0
Theta(0.3)
Se passiamo a Theta un array o qualunque altro oggetto iterabile otteniamo un errore:
Theta((1,2))
Theta(np.array([-3,-2,-1,0,1,2,3]))
Per ottenere una versione vettorializzata di Theta si può utilizzare la funzione di Numpy vectorize. In molti casi vectorize riesce a creare una funzione vettoriale che opera come la funzione scalare ma su tutti gli elementi di qualunque array:
Theta_vec = np.vectorize(Theta) # Theta_vec è una funzione diversa da Theta
Theta_vec(np.array([-3,-2,-1,0,1,2,3]))
a = np.array([1,2,3])
b = a
a[2] = 22
print(b)
La copia con [:] non è sufficiente a creare una deep copy
a = np.array([1,2,3])
b = a[:]
a[2] = 22
print(b)
Bisogna usare la funzione deepcopy nlla libreria copy.
import copy
a = np.array([1,2,3])
b = copy.deepcopy(a)
a[2] = 22
print(b)
Tipicamente le funzioni di numpy possono operare su array multidimensionali. È quindi necessario poter controllare come la funzione agisce sull'array. Un caso tipico è far operare la funzione lungo un singolo asse (axis) dell'array. Un semplice esempio usando la funzione sum.
np.sum([[0, 1], [3, 5]])
np.sum([[0, 1], [3, 5]], axis=0)
np.sum([[0, 1], [3, 5]], axis=1)
Quando si debba trattare una grande massa di dati è essenziale poterli importare ed esportare su file in modo veloce.
Numpy fornisce la funzione loadtxt:
import numpy as np
dataPt, time, height, error = np.loadtxt( "../Data/CadutaLiberaDati.txt", skiprows=5 , unpack=True)
Il parametro skiprows=5 dice a loadtxt di ignorare le prime cinque righe che costituiscono l'header del file. Il parametro unpack=True dice di estrarre i dati caricandoli nell'ntupla di array a primo membro. Il comando assume che i dati siano separati da uno o più spazi.
print(height)
Se i dati sono contenuti in un Comma-Separated Values file (CSV), una alternativa molto comune, è sufficiente specificare il separatore con la keyword delimiter.
dataPt, time, height, error = np.loadtxt("../Data/CadutaLiberaDati.csv", skiprows=5 , unpack=True, delimiter=',')
Si possono importare da un file con dati testuali e numerici usando dtype=str. In questo caso i dati sono rappresentati come stringhe.
gender, weight, age = np.loadtxt("../Data/mixed_data.txt", skiprows=2 , unpack=True, dtype=str)
print(gender)
print(weight)
print(age)
È poi possibile convertire parte degli array in un formato numerico usando il metodo astype:
age_n = age.astype(float)
print(age_n)
Esistono metodi ancora più raffinati come genfromtxt. Usate help(np.genfromtxt) oppure una ricerca sul web per ulteriori dettagli ed esempi
import numpy as np
help(np.genfromtxt)
np.genfromtxt("../Data/mixed_data.txt", skip_header=2 , dtype=None, encoding=None,
names=('gender', 'weight', 'age'))
Esempio usando gli array importati in precedenza:
info = "Misure esperimento Caduta Libera\n"
info +="Data: 13 Dicembre 2019\n"
info +="\n\n"
info +=" Punto Tempo(sec) Altezza(m) Incertezza(m)\n"
np.savetxt('../ShellPrograms/CLD.txt',
list(zip(dataPt, time, height, error)),
header=info, fmt="%12.1f")
Se si preferisce un file CSV:
np.savetxt('../ShellPrograms/CLD.csv',
list(zip(dataPt, time, height, error)),
header=info, fmt="%12.1f", delimiter=",")
La funzione zip, partendo da due o più oggetti iterabili (stringhe, liste, ntuple) di uguale lunghezza, restituisce un iteratore che combina
il contenuto degli oggetti di partenza. Per ottenere il risultato finale in un colpo solo, invece che un elemento alla volta, si passa l'iteratore generato da zip al comando list.
Esempi:
l1 = [1,2,3]
l2 = ['a','b','c']
list(zip(l1,l2))
l1 = 'pippo'
l2 = 'pluto'
list(zip(l1,l2))
l1 = [1,2,3,4,5]
l2 = 'pluto'
list(zip(l1,l2))
numeri quadrati.
…si possono trovare in: http://www.scipy.org/Numpy_Example_List