Na dnešní lekci si do virtuálního prostředí nainstalujte následující balíčky:

$ python -m pip install --upgrade pip
$ python -m pip install notebook numpy scipy imageio matplotlib pillow

Mezitím, co se to bude instalovat, si stáhněte do adresáře static tyto soubory:

A až bude nainstalováno, spusťte si nový Notebook. (Viz lekce o Notebooku.)


NumPy #

NumPy je základní knihovna pro vědce a analytiky, kteří pracují s Pythonem. Definuje typ pro n-rozměrné homogenní pole (nejčastěji čísel) a API pro práci s takovým polem.

Téměř všechny knihovny, kde se objevují větší matice či tabulky, jsou buď postavené na NumPy, nebo podporují numpy.array: od pandas pro datovou analýzu a matplotlib pro grafy, přes scipy, kde najdete základní algoritmy pro interpolaci, integraci aj., astrofyzikální astropy, librosa pro analýzu hudby, až po integraci v knihovnách jako Pillow nebo Tensorflow.

Podobně jako „Djangonauti” kolem webového frameworku Django tvoří vědci a datoví analytici podskupinu pythonní komunity s vlastními konferencemi (PyData), organizacemi (NumFocus, Continuum Analytics) a knihovnami jako NumPy, Pandas, SciPy, Matplotlib či Astropy. Potřeby této komunity se samozřejmě odrážejí i v Pythonu samotném (např. ... a @, které si ukážeme dále, byly do jazyka přidány pro ulehčení výpočtů) a naopak (na rozdíl od specializovaných jazyků jako R nebo Matlab se tu stále indexuje od nuly). Většina těchto knihoven ale má jednu zvláštnost, kterou ve zbytku pythonního světa tolik nevidíme: důraz na použití v interaktivním režimu.

Nejednoznačnost a zkratky #

Čísla můžeme buď prozkoumávat, hrát si s nimi, zjišťovat zajímavé souvislosti; anebo můžeme připravovat programy, které nějaké výpočty provedou automaticky. Na obojí se používají podobné nástroje. Automaticky pouštěné skripty musí být samozřejmě robustní. Nástroje ke zkoumání dat ale bývají přívětivé k vědcům, často na úkor robustnosti nebo „dobrých programátorských mravů”. Například některé funkce tak trochu „hádají”, co uživatel chtěl, a v tutoriálech se setkáte se zkratkami jako import numpy as np či dokonce from numpy import *.

Toto je kurz programovací, kde nám záleží více na znovupoužitelném kódu než na jednom konkrétním výsledku. Budeme proto preferovat explicitní a jednoznačné operace. Ty jsou v použitých knihovnách vždy vedle zkratek k dispozici a popsány v dokumentaci.

Jak s polem pracovat? Nejprve si ho vytvoříme, nejlépe ze seznamu čísel:

In [1]:
import numpy
array = numpy.array([[1, 2, 3], [4, 5, 6], [7, 8, 9]])
array
Out[1]:
array([[1, 2, 3],
       [4, 5, 6],
       [7, 8, 9]])

Každé pole má dvě základní, neměnné vlastnosti: tvar (shape), neboli velikost, a datový typ (dtype).

In [2]:
array.shape
Out[2]:
(3, 3)
In [3]:
array.dtype
Out[3]:
dtype('int64')

Tvar je n-tice, kde n je počet dimenzí pole; shape=(3, 3) dtype='int64' znamená pole 3×3 celých čísel. Dimenzí může být libovolně mnoho, např. trojrozměrnou matici můžeme vytvořit z trojnásobně vnořených seznamů a NumPy ji „rozumně” vypíše:

In [4]:
cube = numpy.array([[[1, 2], [3, 4]], [[5, 6], [7, 8]]])
print(cube.shape)
cube
(2, 2, 2)
Out[4]:
array([[[1, 2],
        [3, 4]],

       [[5, 6],
        [7, 8]]])

Základní operace #

Základní vlastnost NumPy polí je to, že aritmetické operace se skalárními hodnotami (např. čísly) fungují po prvcích.

In [5]:
array
Out[5]:
array([[1, 2, 3],
       [4, 5, 6],
       [7, 8, 9]])
In [6]:
array - 1
Out[6]:
array([[0, 1, 2],
       [3, 4, 5],
       [6, 7, 8]])
In [7]:
array / 2
Out[7]:
array([[0.5, 1. , 1.5],
       [2. , 2.5, 3. ],
       [3.5, 4. , 4.5]])
In [8]:
-(array % 4)
Out[8]:
array([[-1, -2, -3],
       [ 0, -1, -2],
       [-3,  0, -1]])

Kromě aritmetických operací takto funguje i porovnávání. Následujícím výrazem dostanete pravdivostní tabulku, která má True na místech, kde pro příslušný prvek platí podmínka:

In [9]:
array > 5
Out[9]:
array([[False, False, False],
       [False, False,  True],
       [ True,  True,  True]])

Protože Python neumožňuje předefinovat chování operátorů and a or, logické spojení operací se tradičně dělá přes bitové operátory & (a) a | (nebo). Ty mají ale neintuitivní prioritu, proto se jednotlivé výrazy hodí uzavřít do závorek:

In [10]:
(array > 3) & (array < 7)
Out[10]:
array([[False, False, False],
       [ True,  True,  True],
       [False, False, False]])

Operace s jinými poli pracují po prvcích. Obě pole musí být stejně velké.

In [11]:
array + array
Out[11]:
array([[ 2,  4,  6],
       [ 8, 10, 12],
       [14, 16, 18]])
In [12]:
array * numpy.array([[0, 1, 0], [1, 0, 1], [0, 1, 0]])
Out[12]:
array([[0, 2, 0],
       [4, 0, 6],
       [0, 8, 0]])

Sekvence (jako seznamy) jsou před operací převedeny na pole.

In [13]:
array * [[0, 1, 0], [1, 0, 1], [0, 1, 0]]
Out[13]:
array([[0, 2, 0],
       [4, 0, 6],
       [0, 8, 0]])

Indexování #

NumPy pole jde různými způsoby indexovat. Základní způsoby známe už ze samotného Pythonu – pole se chová jako sekvence menších polí:

In [14]:
array
Out[14]:
array([[1, 2, 3],
       [4, 5, 6],
       [7, 8, 9]])
In [15]:
array[0]
Out[15]:
array([1, 2, 3])
In [16]:
array[0:-1]
Out[16]:
array([[1, 2, 3],
       [4, 5, 6]])
In [17]:
array[0][1]
Out[17]:
2

Na rozdíl od Pythonu se ale dá vícerozměrné pole indexovat přímo n-ticí. Toto je dokonce preferovaný způsob – přehlednější a mnohem rychlejší:

In [18]:
array[0, 1]
Out[18]:
2
In [19]:
array[0:-1, 1:]
Out[19]:
array([[2, 3],
       [5, 6]])

Chceme-li vybrat určitý sloupec, řekneme si „kompletním intervalem“ (:) o všechny řádky:

In [20]:
array[:, 2]
Out[20]:
array([3, 6, 9])

Další způsob specifický pro NumPy je indexování pravdivostní tabulkou.

Když potřebujete z matice vybrat prvky, pro které platí nějaká podmínka, napřed si vytvořte pole hodnot True/False:

In [21]:
array > 4
Out[21]:
array([[False, False, False],
       [False,  True,  True],
       [ True,  True,  True]])

To pak použijte jako index:

In [22]:
array[array > 4]
Out[22]:
array([5, 6, 7, 8, 9])

Výsledek je jednorozměrný – informace o původních pozicích prvků se ztratí.

Pro úplnost uvedu ještě dva způsoby indexování. První je seznamem indexů, kterým můžete vybírat, přehazovat nebo i duplikovat konkrétní řádky:

In [23]:
array[[0, 2, 1, 1]]  # Řádky 0, 2, 1 a 1
Out[23]:
array([[1, 2, 3],
       [7, 8, 9],
       [4, 5, 6],
       [4, 5, 6]])
In [24]:
array[:, [2, 2, 0, 0]]  # Podobně pro sloupce
Out[24]:
array([[3, 3, 1, 1],
       [6, 6, 4, 4],
       [9, 9, 7, 7]])

Druhý je indexování pomocí n-tice „řádkových souřadnic“ a n-tice odpovídajících „sloupcových souřadnic“:

In [25]:
array[(0, 1, 2), (1, 2, 0)]  # Vybere prvky (0, 1), (1, 2), (2, 0)
Out[25]:
array([2, 6, 7])

Trochu specifické je indexování vícerozměrných polí. U nich se často využije „kompletní interval“ (:):

In [26]:
cube
Out[26]:
array([[[1, 2],
        [3, 4]],

       [[5, 6],
        [7, 8]]])
In [27]:
cube[0, :, :]  # První „vrstva“ - to samé jako cube[0]
Out[27]:
array([[1, 2],
       [3, 4]])
In [28]:
cube[:, 0, :]  # První řádky - to samé jako cube[:, 0]
Out[28]:
array([[1, 2],
       [5, 6]])
In [29]:
cube[:, :, 0]  # První sloupce
Out[29]:
array([[1, 3],
       [5, 7]])

Má-li pole hodně rozměrů, je psaní spousty :, zdlouhavé a nepřehledné. Existuje proto speciální hodnota ... (Ellipsis), která doplní tolik „kompletních intervalů“ (:), aby souhlasil počet dimenzí:

In [30]:
cube[..., 0]  # První sloupce – ekvivalent [:, :, 0]
Out[30]:
array([[1, 3],
       [5, 7]])

Broadcasting a změny #

Už jsme si ukázali, že aritmetické operace se skalárními hodnotami se provede pro všechny prvky, zatímco operace mezi dvěma stejně velkými poli se provede po prvcích:

In [31]:
array
Out[31]:
array([[1, 2, 3],
       [4, 5, 6],
       [7, 8, 9]])
In [32]:
array * 3
Out[32]:
array([[ 3,  6,  9],
       [12, 15, 18],
       [21, 24, 27]])
In [33]:
array * [[0, 1, 0], [1, 0, 1], [0, 1, 0]]
Out[33]:
array([[0, 2, 0],
       [4, 0, 6],
       [0, 8, 0]])

Jak je to ale s různě velkými poli?

Nemá-li sekvence, se kterou pracujeme, dost dimenzí, poslední dimenze se „rozšíří“, jako bychom pracovali v každém sloupci se skalární hodnotou. Tomuto „rozšiřování” se obecně říká broadcasting.

In [34]:
array
Out[34]:
array([[1, 2, 3],
       [4, 5, 6],
       [7, 8, 9]])
In [35]:
array * [0, 1, 10]   # vynásobí 1. sloupec nulou, 2. jedničkou, 3. deseti
Out[35]:
array([[ 0,  2, 30],
       [ 0,  5, 60],
       [ 0,  8, 90]])

Podobné rozšiřování nastane, má-li některá dimenze velikost 1:

In [36]:
array * [[0], [1], [10]]  # vynásobí 1. *řádek* nulou, atd.
Out[36]:
array([[ 0,  0,  0],
       [ 4,  5,  6],
       [70, 80, 90]])

Jednotlivé hodnoty v poli lze měnit:

In [37]:
array[0, 0] = 123
array
Out[37]:
array([[123,   2,   3],
       [  4,   5,   6],
       [  7,   8,   9]])

...a i na měnění se vztahuje broadcasting:

In [38]:
array[0] = 123
array
Out[38]:
array([[123, 123, 123],
       [  4,   5,   6],
       [  7,   8,   9]])
In [39]:
array[:] = 123
array
Out[39]:
array([[123, 123, 123],
       [123, 123, 123],
       [123, 123, 123]])

Obecně platí, že lze-li něčím vybírat prvky, lze tím i pole měnit:

In [40]:
array[(1, 2, 0), (0, 2, 1)] = 0
array
Out[40]:
array([[123,   0, 123],
       [  0, 123, 123],
       [123, 123,   0]])

Další způsob, jak pole měnit, je rozšířeným přiřazením.

In [41]:
array *= 2
array
Out[41]:
array([[246,   0, 246],
       [  0, 246, 246],
       [246, 246,   0]])

Tato operace není totéž co array = array * 2, ačkoli má stejný výsledek.

array *= 2 změní existující pole, zatímco array = array * 2 vytvoří nové pole, které pak přiřadí do původní proměnné.

Pozor na to, že není možné měnit typ pole:

In [42]:
try:
    array /= 2
except Exception as e:
    print("Chyba!!", type(e), e)
Chyba!! <class 'TypeError'> No loop matching the specified signature and casting was found for ufunc true_divide

Tvoření matic, část 2 #

Časté druhy matic se dají vytvořit pomocí pomocných funkcí. Výsledky se dají použít přímo nebo naplnit vypočítanými daty:

In [43]:
numpy.zeros((4, 4))
Out[43]:
array([[0., 0., 0., 0.],
       [0., 0., 0., 0.],
       [0., 0., 0., 0.],
       [0., 0., 0., 0.]])
In [44]:
numpy.ones((4, 4))
Out[44]:
array([[1., 1., 1., 1.],
       [1., 1., 1., 1.],
       [1., 1., 1., 1.],
       [1., 1., 1., 1.]])
In [45]:
numpy.full((4, 4), 12.34)
Out[45]:
array([[12.34, 12.34, 12.34, 12.34],
       [12.34, 12.34, 12.34, 12.34],
       [12.34, 12.34, 12.34, 12.34],
       [12.34, 12.34, 12.34, 12.34]])
In [46]:
numpy.eye(4)  # Jednotková matice (je čtvercová – n×n)
Out[46]:
array([[1., 0., 0., 0.],
       [0., 1., 0., 0.],
       [0., 0., 1., 0.],
       [0., 0., 0., 1.]])
In [47]:
numpy.diag([1, 2, 3, 4])  # Diagonální matice
Out[47]:
array([[1, 0, 0, 0],
       [0, 2, 0, 0],
       [0, 0, 3, 0],
       [0, 0, 0, 4]])

U těchto funkcí lze obecně použít argument dtype, kterým specifikujeme datový typ:

In [48]:
int_zeros = numpy.zeros((4, 4), dtype='int8')
print(int_zeros.dtype)
int_zeros
int8
Out[48]:
array([[0, 0, 0, 0],
       [0, 0, 0, 0],
       [0, 0, 0, 0],
       [0, 0, 0, 0]], dtype=int8)

Další funkce tvoří jednorozměrné matice. Základní je arange, která bere stejné argumenty jako range v Pythonu:

In [49]:
numpy.arange(50)  # Celočíselné – argumenty jako range() v Pythonu
Out[49]:
array([ 0,  1,  2,  3,  4,  5,  6,  7,  8,  9, 10, 11, 12, 13, 14, 15, 16,
       17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33,
       34, 35, 36, 37, 38, 39, 40, 41, 42, 43, 44, 45, 46, 47, 48, 49])

Navíc umí pracovat s reálnými čísly (float). Pozor ale na to, že reálná čísla jsou nepřesná! arange k začátku sekvence postupně přičítá „krok”, takže chyba narůstá celkem rychle:

In [50]:
numpy.arange(0.0, 50.0, 0.3)[-1]
Out[50]:
49.799999999999997

V krajních případech takto dokonce můžeme dostat pole jiné velikosti, než jsme zamýšleli. Proto arange používejte jen pro celá čísla; pro reálná je tu linspace, která bere začátek a konec intervalu, plus počet prvků:

In [51]:
numpy.linspace(0, 50, num=11)  # vždy 11 prvků
Out[51]:
array([ 0.,  5., 10., 15., 20., 25., 30., 35., 40., 45., 50.])

Reshape #

Ačkoli indexování polí v NumPy je dost mocné, v paměti jsou jednotlivé hodnoty reprezentovány jako (metadata a) jednorozměrné pole, známé z jazyka C (ačkoli samotné rozmístění prvků může být jiné než po řádcích, jak jsme zvyklí u C).

Je jednoduché změnit tvar pole, nezmění-li se tím celkový počet prvků:

In [52]:
array = numpy.arange(12)
array
Out[52]:
array([ 0,  1,  2,  3,  4,  5,  6,  7,  8,  9, 10, 11])
In [53]:
reshaped = array.reshape((3, 4))
reshaped
Out[53]:
array([[ 0,  1,  2,  3],
       [ 4,  5,  6,  7],
       [ 8,  9, 10, 11]])

Pozor na to, že reshape sice vrací nový objekt, ale může (ale nemusí!) to být jen nový pohled na existující data. Změny v pohledu se projeví i v původním poli:

In [54]:
reshaped[2, 2] = -99
reshaped
Out[54]:
array([[  0,   1,   2,   3],
       [  4,   5,   6,   7],
       [  8,   9, -99,  11]])
In [55]:
array
Out[55]:
array([  0,   1,   2,   3,   4,   5,   6,   7,   8,   9, -99,  11])

Podobně tvoří pohledy i jednoduché indexování:

In [56]:
a_slice = array[2:4]
a_slice[:] = -99, -99

array
Out[56]:
array([  0,   1, -99, -99,   4,   5,   6,   7,   8,   9, -99,  11])

Potřebujete-li kopii dat, použijte metodu copy:

In [57]:
array.reshape((3, 4)).copy()
Out[57]:
array([[  0,   1, -99, -99],
       [  4,   5,   6,   7],
       [  8,   9, -99,  11]])

Podobně jako reshape funguje transpozice, což je tak častá operace, že má jednopísmennou zkratku – atribut T. (Tohle hodně napomáhá tomu, že zápis maticových výpočtů v NumPy se podobá odpovídajícím matematickým vzorcům.)

In [58]:
reshaped.T
Out[58]:
array([[  0,   4,   8],
       [  1,   5,   9],
       [-99,   6, -99],
       [-99,   7,  11]])

Až budete NumPy zkoušet, doporučuji si u nových funkcí najít, zda tvoří nová pole, vracejí pohled nebo modifikují původní pole. U některých funkcí najdete pojmenovaný argument inplace (modifikuje původní pole), případně out, („naplní“ jiné, existující pole).

Datové typy #

Jak už bylo řečeno, matice v NumPy mají určené datové typy. Ty jdou nastavit ve většině funkcí, které matice tvoří, argumentem dtype:

In [59]:
numpy.zeros(4, dtype=int)
Out[59]:
array([0, 0, 0, 0])
In [60]:
numpy.zeros(4, dtype=float)
Out[60]:
array([0., 0., 0., 0.])
In [61]:
numpy.zeros(4, dtype=bool)
Out[61]:
array([False, False, False, False])

Nejobecnější typ je object (jehož použitím ale přicházíme o většinu výhod NumPy).

In [62]:
numpy.array([0, 1.3, "foobar"], dtype=object)
Out[62]:
array([0, 1.3, 'foobar'], dtype=object)

Kromě pythonních typů bere dtype i řetězcové specifikace typu:

In [63]:
numpy.array([1, 8, 500], dtype='int8')  # 8bitové celé číslo
Out[63]:
array([  1,   8, -12], dtype=int8)

Znáte-li modul array ze standardní knihovny, můžete jeho specifikaci použít i tady:

In [64]:
numpy.zeros(4, dtype='<I')
Out[64]:
array([0, 0, 0, 0], dtype=uint32)

Navíc dtype umí řetězcové a bytestring typy. Tyto mají danou maximální velikost a nesmí obsahovat \0 (resp. znakem '\0' jsou ukončeny):

In [65]:
numpy.full(4, 'abcdef', dtype=('U', 10))  # Unicode
Out[65]:
array(['abcdef', 'abcdef', 'abcdef', 'abcdef'], dtype='<U10')
In [66]:
numpy.full(4, 'abcdef', dtype=('a', 3))  # "ASCII"
Out[66]:
array([b'abc', b'abc', b'abc', b'abc'], dtype='|S3')

Typy v NumPy můžou být poměrně složité; např. existují i složené datové typy (records). Ty nebudeme používat, ale je dobré o nich aspoň tušit:

In [67]:
record_type = numpy.dtype([('a', int), ('b', float), ('c', ('U', 3))])
numpy.array([(1, 2, 'abc')] * 4, record_type)
Out[67]:
array([(1, 2., 'abc'), (1, 2., 'abc'), (1, 2., 'abc'), (1, 2., 'abc')],
      dtype=[('a', '<i8'), ('b', '<f8'), ('c', '<U3')])

Maticové násobení #

Kromě základních aritmetických operací se u vícerozměrných polí často setkáme s maticovým násobením. Předpokládám, že jako bakaláři jste se s ním už setkali a tušíte co dělá – jestli ne, tuto sekci ignorujte.

V Pythonu 3.5 byl na výzvu vědecké komunity do jazyka přidán operátor @ (mATrix multiplication), který je vyhrazen pro maticové násobení. V samotném Pythonu ani ve standardní knihovně není typ, který ho podporuje, ale matice v NumPy tuto operaci samozřejmě umí.

In [68]:
array1 = numpy.array([[1, 2, 3], [4, 5, 6], [7, 8, 9]])
array2 = numpy.array([[1, 0, 0], [0, 2, 0], [0, 0, 3]])
array1 @ array2
Out[68]:
array([[ 1,  4,  9],
       [ 4, 10, 18],
       [ 7, 16, 27]])

Ve starších verzích Pythonu je potřeba používat metodu nebo funkci dot.

Důvod přidání operátoru @ byl prostý – zjednodušení zápisu maticových operací. Jako příklad uvedený v návrhu je uveden tento vzorec pro testování hypotéz v lineárním regresním modelu:

$ S=(H\beta-r)^T(HVH^T)^{-1} (H\beta-r) $

V NumPy se dá přepsat jako:

from numpy import dot
from numpy.linalg import inv, solve

Pomocí funkce dot:

S = dot((dot(H, beta) - r).T,
       dot(inv(dot(dot(H, V), H.T)), dot(H, beta) - r))

Pomocí metoody dot:

S = (H.dot(beta) - r).T.dot(inv(H.dot(V).dot(H.T))).dot(H.dot(beta) - r)

Pomocí operátoru @:

S = (H @ beta - r).T @ inv(H @ V @ H.T) @ (H @ beta - r)

Poslední varianta nápadně připomíná původní vzorec; u prvních dvou se člověk snadno ztratí ve změti závorek.

Booleovské hodnoty polí #

Použijeme-li pole v příkazu if, NumPy nám vynadá. Standardní pythonní seznam je „pravdivý“ pokud obsahuje nějaké prvky, ale u pole, které má fixní velikost, je tahle informace téměř vždy zbytečná.

In [69]:
try:
    if numpy.eye(3):
        pass
except ValueError as e:
    print("Chyba!", type(e), e)
Chyba! <class 'ValueError'> The truth value of an array with more than one element is ambiguous. Use a.any() or a.all()

Musíme říct přesně, co chceme.

In [70]:
if numpy.eye(3).any():
    print('Alespoň jeden prvek je nenulový')
if numpy.eye(3).all():
    print('Všechny prvky jsou nenulové')
if numpy.eye(3).size:
    print('Pole obsahuje nějaké prvky')
Alespoň jeden prvek je nenulový
Pole obsahuje nějaké prvky

Z historických důvodů existují dvě výjimky: pole s právě jedním prvkem má pravdivostní hodnotu podle daného prvku a prázdné pole je „nepravdivé“:

In [71]:
if numpy.ones((1, 1, 1, 1)):
    print('Ano')
if numpy.zeros((1, 1, 1, 1)):
    print('Ne')
Ano
In [72]:
if numpy.ones((0, 0)):
    print('Ano')
DeprecationWarning: The truth value of an empty array is ambiguous. Returning False, but in future this will result in an error. Use `array.size > 0` to check that an array is not empty.

Další operace #

Modul numpy obsahuje spoustu základních funkcí, které pracují s maticemi; mimo jiné většinu funkcií z pythonního modulu math. Oproti math zvládají funkce z NumPy broadcasting.

In [73]:
array = numpy.linspace(0, numpy.pi, num=1000)
array[:10]
Out[73]:
array([0.        , 0.00314474, 0.00628947, 0.00943421, 0.01257895,
       0.01572369, 0.01886842, 0.02201316, 0.0251579 , 0.02830264])
In [74]:
sine = numpy.sin(array)
sine[:10]
Out[74]:
array([0.        , 0.00314473, 0.00628943, 0.00943407, 0.01257862,
       0.01572304, 0.0188673 , 0.02201138, 0.02515525, 0.02829886])

Další operace doporučuji hledat buď v Notebooku přes tab, v dokumentaci, nebo obecně na Internetu (kde najdete i případné knihovny, které implementují operace, které v NumPy nejsou).

Příklady použití #

Dost teorie. Tahle n-rozměrná pole se používají na spoustu zajímavých věcí. Podívejme se na některé příklady.

Matematika a grafy #

Jak se používají matice, jistě znáte z matematiky a cílem tohoto kurzu není vás to naučit. Ukážu ale pár ochutnávek.

Použijeme knihovnu Matplotlib, která vykresluje grafy. Jak ji použít dohledáte v dokumentaci nebo – často efektivněji – v galerii příkladů.

Matplotlib nemá automatickou integraci s Jupyter Notebookem, proto ji je potřeba po importu zapnout:

In [75]:
from matplotlib import pyplot

# Zapnutí integrace s notebookem – `%` je "magický" příkaz IPythonu, podobně jako `!` pro shell
%matplotlib inline

A teď můžeme nakreslit třeba graf funkce:

In [76]:
x = numpy.linspace(0, numpy.pi * 4)   # definiční obor
y = numpy.sin(x)                      # odpovídající hodnoty funkce

pyplot.plot(x, y)
Out[76]:
[<matplotlib.lines.Line2D at 0x7fab9ebd62d0>]
No description has been provided for this image
In [77]:
s = numpy.linspace(-10, 10, num=100)

# meshgrid vrátí dvě 100×100 matice:
# - jednu, kde v každém řádku jsou čísla od -10 do 10,
# - druhou, kde v každém sloupci jsou čísla od -10 do 10.
x, y = numpy.meshgrid(s, s)

# vyhodnotíme (x**2 + y**2) pro každý prvek
z = x ** 2 + y ** 2

# výsledek vykreslíme jako obrázek
pyplot.imshow(z)

# přidáme legendu
pyplot.colorbar()
Out[77]:
<matplotlib.colorbar.Colorbar at 0x7fab9e6e56d0>
No description has been provided for this image
In [78]:
# Ta samá data můžeme vykreslit i ve 3D
from mpl_toolkits.mplot3d import Axes3D

fig = pyplot.figure()
axes = fig.gca(projection='3d')

surf = axes.plot_surface(x, y, z, cmap='viridis')
No description has been provided for this image

Obrázky #

Typický barevný obrázek není nic než matice $m \times n \times 3$ čísel: $m \times n$ pixelů na šířku a výšku a 3 kanály pro červenou, zelenou a modrou barvu.

Knihovna pillow (nástupce knihovny PIL, který se stále importuje jako PIL) obsahuje nástroje na práci s obrázky, např. „nakresli čáru“ nebo „převeď na černobílý obrázek“ nebo „načti PNG“. Není postavena přímo na NumPy, ale umí obrázky převádět z a na NumPy pole, pokud máme NumPy nainstalované.

Knihovna imageio slouží k načítání a ukládání různých (njen obrázkových) formátů do/z NumPy matic.

Nás bude zajímat funkce imageio.imread, která pomocí Pillow/PIL načte obrázek jako 3D matici 8bitových čísel. Já načtu obrázek hada, vy najděte na internetu jakýkoli RGB obrázek a načtěte si ten.

Použitý obrázek je stažený z Wikimedia Commons a je pod licencí CC-BY-SA 3.0. Autor je uživatel Mokele na anglické Wikipedii.

In [79]:
import imageio
img = imageio.imread('static/python.jpg')
img
Out[79]:
Array([[[172, 173, 165],
        [173, 174, 166],
        [173, 174, 168],
        ...,
        [172, 171, 167],
        [173, 172, 168],
        [173, 172, 168]],

       [[175, 176, 168],
        [174, 175, 167],
        [172, 173, 167],
        ...,
        [172, 171, 167],
        [173, 172, 168],
        [174, 173, 169]],

       [[176, 177, 169],
        [174, 175, 167],
        [171, 172, 166],
        ...,
        [173, 172, 168],
        [173, 172, 168],
        [172, 171, 167]],

       ...,

       [[209, 211, 206],
        [210, 212, 207],
        [211, 213, 208],
        ...,
        [202, 203, 197],
        [200, 201, 195],
        [200, 201, 195]],

       [[208, 210, 205],
        [209, 211, 206],
        [209, 211, 206],
        ...,
        [202, 203, 197],
        [202, 203, 197],
        [203, 204, 198]],

       [[207, 209, 204],
        [209, 211, 206],
        [210, 212, 207],
        ...,
        [201, 202, 196],
        [201, 202, 196],
        [201, 202, 196]]], dtype=uint8)

Pomocí nám už známé knihovny matplotlib takovou matici můžeme zobrazit:

In [80]:
pyplot.imshow(img)
Out[80]:
<matplotlib.image.AxesImage at 0x7fabd75fcc50>
No description has been provided for this image

Podívejme se teď na strukturu matice:

In [81]:
img.shape
Out[81]:
(887, 1037, 3)

První rozměr jsou řádky (y); můj obrázek je 887 pixelů vysoký. Druhý jsou sloupce (x); tento obrázek má na šířku 1037 px. Třetí rozměr jsou barevné kanály.

Pomocí indexování se můžeme na jednotlivé barevné kanály dostat: je to poslední index, takže řádky a sloupce nahradíme buď dvěma kompletními intervaly (:, :) nebo vynechávkou (...). Červený kanál tedy bude [..., 1], modrý [..., -1].

Zobrazení chceme černobílé; na to má matplotlib pojmenovaný argument cmap. Výchozí způsob obarvování je vhodný spíše pro grafy funkcí.

In [82]:
blue_channel = img[..., -1]
pyplot.imshow(blue_channel, cmap='gray')
Out[82]:
<matplotlib.image.AxesImage at 0x7fabd75bfc50>
No description has been provided for this image

Zajímavé využití obrázku jako matice je steganografie: ukrytí informace v obrazových datech.

Načteme jiný obrázek stejné velikosti, tentokrát černobílý (s módem L). Informace v něm schováme do posledního bitu modrého kanálu.

In [83]:
secret = imageio.imread('static/secret.png', pilmode='L')

img[..., -1] = (img[..., -1] & 0b11111110) + (secret.astype(bool))

Obrázek vypadá na první pohled stejně...

In [84]:
pyplot.imshow(img)
Out[84]:
<matplotlib.image.AxesImage at 0x7fabd75813d0>
No description has been provided for this image

... ale v posledím modrém bitu se skrývá tajná informace.

In [85]:
pyplot.imshow(img[..., -1] & 1, cmap='gray')
Out[85]:
<matplotlib.image.AxesImage at 0x7fab9c631d90>
No description has been provided for this image

Výsledek je dobré uložit v bezztrátovém formátu (PNG), aby se informace neztratila:

In [86]:
imageio.imsave('python.png', img)

Zvuk #

Jako pole lze reprezentovat i zvukový záznam. Mám záznam, na kterém zkouším zpívat; pomocí funkce scipy.io.wavfile ho můžu načíst jako NumPy pole:

In [87]:
import scipy.io.wavfile
sample_rate, sound = scipy.io.wavfile.read('static/sample.wav')
print(sample_rate)
sound
48000
Out[87]:
array([[ -58,  -88],
       [ -65,  -49],
       [  56,  -18],
       ...,
       [   1,  231],
       [ -85,  234],
       [-118,  212]], dtype=int16)

Zvuk je stereo, má dvě stopy; jednu z nich si vykreslím:

In [88]:
channel = sound[..., 1]
pyplot.plot(channel)
Out[88]:
[<matplotlib.lines.Line2D at 0x7fab9c051b90>]
No description has been provided for this image

Případně můžu využít možností Jupyter Notebooku a HTML a zvuk si přehrát:

In [89]:
from IPython.display import Audio
Audio(data=channel, rate=sample_rate)
print('(Zkuste si to sami; tento print vymažte)')
(Zkuste si to sami; tento print vymažte)

Podívám se na detail první „noty”, kde je patrná vlna s nějakou frekvencí:

In [90]:
segment = channel[50000:55000]
pyplot.plot(segment)
Out[90]:
[<matplotlib.lines.Line2D at 0x7fab96685510>]
No description has been provided for this image
In [91]:
from IPython.display import Audio
Audio(data=segment, rate=sample_rate)
Out[91]:

Jaká to je frekvence? Znáte-li analýzu signálů, tušíte, že na podobné otázky odpovídá Fourierova transformace. NumPy obsahuje diskrétní Fourierovu transformaci v modulu numpy.fft spolu s funkcí, která spočítá odpovídající frekvence v Hz:

In [92]:
spectrum = numpy.fft.fft(segment)
freqs = numpy.fft.fftfreq(len(spectrum), 1/sample_rate)
pyplot.xlabel('Frekvence (Hz)')
pyplot.plot(freqs, numpy.abs(spectrum))
Out[92]:
[<matplotlib.lines.Line2D at 0x7fab965a4ad0>]
No description has been provided for this image

V tomto grafu hledám maximum. Můžu se zaměřit na prvních pár hodnot spektra:

In [93]:
pyplot.xlabel('Frekvence (Hz)')
pyplot.plot(freqs[:100], numpy.abs(spectrum[:100]))
Out[93]:
[<matplotlib.lines.Line2D at 0x7fab96514910>]
No description has been provided for this image

Maximum je někde kolem 120 Hz; abych to zjistil přesně, použiji funkci argmax:

In [94]:
amax = numpy.argmax(abs(spectrum))
amax
Out[94]:
13

... a najdu odpovídající frekvenci:

In [95]:
freqs[amax]
Out[95]:
124.80000000000001

Což je podle seznamu not skoro H$_2$ (123,5 Hz).