Tudományos csomagok (Python Scientific stack)¶
| csomag | |
|---|---|
| NumPy | Hatékony N-dimenziós tömb |
| SciPy | Numerikus számítások |
| Matplotlib | Grafikonok és rajzok |
| IPython (Jupyter) | Interaktív notebook |
| SymPy | Szimbolikus számítások |
| Pandas | Adatbányászat |
Numpy¶
import numpy as np
Vektorizáció¶
Ciklusok helyett automatikus elemenkénti műveletvégzés és optimalizálás a háttérben.
- Egyszerűbb, átlátható kód
- Gyorsabb
a=np.arange(0,2000)
b=np.arange(2000,0,-1)
%%timeit
ans=0
for i in range(len(a)):
ans=ans+a[i]*b[i]
ans
%%timeit
ans=0
for x,y in zip(a,b):
ans=ans+x*y
ans
%%timeit
sum((x*y for x,y in zip(a,b)))
%%timeit
np.dot(a,b)
Copy¶
view(): ugyanahhoz az adathoz fér hozzá és az adtok megváltoztatása az eredeti tömböt is megváltoztatjacopy(): másolatot készít a tömbről, megváltoztatása nem befolyásolja az eredetit
a=np.arange(12).reshape(3,4)
b=a
print(b is a)
c = a.view()
print(c is a)
c = c.reshape((2, 6)) # a alakja nem változik
a.shape
c[0, 4] = 1234
a
A slicing is view-t hoz létre. Viszont a tömbökkel indexelés nem!
s = a[ : , 1:3]
s[:] = 10
a
e=a[[1,2]]
print(e)
e[:]=11
print(e)
a
Copy¶
Teljes másolat jön létre.
d = a.copy()
print(d is a)
d[0,0] = 9999
a
Sokszor megéri egy slicera copyt hívni, ha tömb többi részére később már nem lesz szükségünk. Így memóriát szabadíthatunk fel.
a = np.arange(int(1e8))
b = a[:100].copy()
del a # az a áltál foglalt memória felszabadul
Ha csak b = a[:100] szerepelne, akkor a del a parancs kiadása után sem szabadulna fel a memória, hiszen b még hivatkozik rá.
np.random¶
Numpynak van egy gazdag random számmokkal operáló alcsomagja.
A numpy.random.rand például float64 számokat randomol a [0, 1) intervallumról egyenletesen.
np.random.rand(2, 3, 2)
Más eloszlások:
np.random.uniform(1, 2, (2, 2))
np.random.standard_normal(10)
np.random.normal(10, 1, size=(1,10))
Diszkrét eloszlások:
np.random.choice(["A", "2", "3", "4", "5", "6", "7", "8", "9",
"10", "J", "Q", "K"], 5, replace=True)
choice a valszóínűségeket is megadhatjuk:
np.random.choice(range(1, 7), 10,
p=[0.1, 0.1, 0.1, 0.1, 0.1, 0.5])
print(np.random.permutation(["A", "2", "3", "4", "5", "6",
"7", "8", "9", "10", "J", "Q", "K"]))
A permutation az első dimenzió szerint permutál.
print(np.random.permutation(np.arange(9).reshape((3, 3))))
Lineáris algebra¶
linalg csomag:
https://numpy.org/doc/stable/reference/routines.linalg.html
A mátrix szorzásra három jelölés is létezik:
A = np.array([[1, 2], [3, 4]])
A.dot(A)
A*A
A @ A
np.matmul(A,A)
Tudunk invertálni, de mint mindig, figyelni kell a kerkítési hibákra.
A_inv = np.linalg.inv(A)
print(A_inv)
A_inv.dot(A)
np.round(A_inv.dot(A),5)
u = np.array([[0.0, -1.0], [1.0, 0.0]])
np.trace(u) # mátrix nyoma
Lineáris egyenletrendszerek, sajátértékek
y = np.array([[5.], [7.]])
x=np.linalg.solve(A, y)
x
A @ x
np.linalg.eig(u)
import matplotlib.pyplot as plt
import imageio #kép beolvasáshoz
Legyen $A\in\mathbb{R}^{n\times m}$ egy tetszőleges mátrix. Ekkor léteznek $U\in\mathbb{R}^{n\times n}$, $D\in\mathbb{R}^{n\times m}$ és $V\in\mathbb{R}^{m\times m}$ mátrxiok, melyekre
$$ A = UDV^*, $$ahol $U$ és $V$ unitér mátrixok, tehát $U^*U = UU^* = I_n$ és $V^*V=VV^*=I_m$. $D$ pedig diakonális mátrix, azaz $d_{ij}=0$ ha $i\ne j$. A csillag operátor konjugált transzponáltja, tehát $(V^*)_{ij} = \overline V_{ji}$, de mivel mi valós mátrixokkal dolgozunk, elegendő transzponálásra gondolni. $D$ diagonális elemei nemnegatívak és csökkenő sorrenben jönnek: $d_{ii} = \sigma_i$, $\sigma_1\geq\sigma_2\geq\ldots\geq\sigma_r > \sigma_{r+1} = \ldots = \sigma_{\min(n,m)} = 0$, ahol $r$ az $A$ mátrix rangja. Ezeket a $\sigma$ értékeket az $A$ mátrix szinguláris értékeinek hívjuk.
Alacsony rangú aproximáció¶
Legyen $k\in\mathbb{N}$ egy természetes szám, ahol $k\leq\text{rank}(A)\leq\min\{n, m\}$. Egy olyan $A_k$ mátrixok keresünk, melyre $\text{rank}(A_k) = k$ és amelyik az $A$ mátrix legjobb közelítése a $k$ rangú mátrixok között. Tehát a következő minimalizálási feladatot szeretnénk megoldani:
$$ \left|\left| A - B \right|\right|_F \to \min !\qquad \mbox{ ahol }\quad B\in\mathbb{R}^{n\times m}, \ \text{rank}(B) = k. $$Itt $\left|\left| X \right|\right|_F$ a Frobenius normát jelöli, ami az $X$ eleminek négyzetösszegének a gyöke.
A feladat megoldását mekapjuk az SVD-felbontás segítségével. Ha $A = UDV^*$, akkor megtartjuk $D$ átlójának első $k$ elemét és a többi szinguláris értéket pedig 0-ra állítjuk. Legyen $D_k$ az így kapott diagonális mátrix. Ezután újra kiszámítjuk a $UD_kV^*$ szorzatot. A nullára állított értékek miatt elég $U$ első k oszlopát és $V$ első $k$ sorát megtartani.
Összefoglalva, az $A_k := U_kD_kV_k^*$ mátrix van legközelebb az $A$ mátrixhoz a Frobenius norma szerint a $k$ rangú mátrixok között, ahol $U_k$ és $V_k$ az első $k$ oszlopla és sora $U$-nak és $V$-nek.
# A rangja 2. Ha A[2,1] értéke 4 lenne, akkor a rang csak 1 lenne.
A = np.array([[1, 2, 0, 0, 2, 3, -1, -2]]).reshape((4, 2))
print(A)
n, m = A.shape
rank_A = np.linalg.matrix_rank(A)
print(f"sorok száma: {n}, oszlopok száma: {m}")
print(f"A rangja: {rank_A}")
U, d, V_T = np.linalg.svd(A, full_matrices=True)
print(f"szinguláris értékek: {len(d)}")
D = np.concatenate((np.diag(d), np.zeros((n - len(d), m))), axis=0) # szinguláris értékekből diagonális mátrix
A_restored = np.matmul(U, np.matmul(D, V_T))
A_restored = np.dot(U, np.dot(D, V_T))
A_restored = (U @ D) @ V_T
print(A_restored)
B = np.array([[1, 2, 0, 0, 2, 4, -1, -2]]).reshape((4, 2))
print(B)
n, m = B.shape
rank_B = np.linalg.matrix_rank(B)
print(f"sorok száma: {n}, os zlopok száma: {m}")
print(f"B rangja: {rank_B}")
distance = np.linalg.norm(A - B, "fro")
print(f'Frobenius távolság A és B között: {distance}')
U_1 = U[:, :1]
D_1 = D[:1, :1]
V_T_1 = V_T[:1, :]
A_1 = np.matmul(U_1, np.matmul(D_1, V_T_1))
print(f"A legjobb 1 rangú aproximációj A-nak \n {A_1}")
dist_from_A_1 = np.linalg.norm(A - A_1, "fro")
print("A Frobenius norma szerintei eltérés: {}.".format(dist_from_A_1))
Képek közelítése¶
A képeket a korábban látott módon egy $n\times m\times 3$-as tömbben tároljuk.
A fenti közelítést mindhárom színre elvégezzük.
import imageio
image = imageio.imread('dragon.jpg')
image
image = image / 255
row, col, _ = image.shape
print("felbontás: {} * {}".format(row, col))
import matplotlib.pyplot as plt
fig = plt.figure(figsize=(15, 10))
imgplot = plt.imshow(image)
plt.show()
image_red = image[:, :, 0]
image_green = image[:, :, 1]
image_blue = image[:, :, 2]
U_r, d_r, V_r = np.linalg.svd(image_red, full_matrices=True)
U_g, d_g, V_g = np.linalg.svd(image_green, full_matrices=True)
U_b, d_b, V_b = np.linalg.svd(image_blue, full_matrices=True)
k = 50
U_r_k = U_r[:, :k]
V_r_k = V_r[:k, :]
U_g_k = U_g[:, :k]
V_g_k = V_g[:k, :]
U_b_k = U_b[:, :k]
V_b_k = V_b[:k, :]
d_r_k = d_r[:k]
d_g_k = d_g[:k]
d_b_k = d_b[:k]
image_red_approx = np.matmul(U_r_k, np.matmul(np.diag(d_r_k), V_r_k))
image_green_approx = np.matmul(U_g_k, np.matmul(np.diag(d_g_k), V_g_k))
image_blue_approx = np.matmul(U_b_k, np.matmul(np.diag(d_b_k), V_b_k))
image_reconstructed = np.zeros((row, col, 3))
image_reconstructed[:, :, 0] = image_red_approx
image_reconstructed[:, :, 1] = image_green_approx
image_reconstructed[:, :, 2] = image_blue_approx
fig = plt.figure(figsize=(15, 10))
plt.imshow(image_reconstructed)
plt.show()