Nauka #22 — Algebra liniowa i data science

Algebra liniowa to język współczesnej informatyki. PageRank, kompresja obrazów, systemy rekomendacji, symulacje fizyczne, rendering 3D i sieci neuronowe to w swej istocie ciągi operacji macierzowych. Mimo to algebra liniowa jest często nauczana jako czysta manipulacja symbolami — wiersze i kolumny bez żadnego znaczenia geometrycznego. Ten wpis eksploruje ją w odwrotny sposób: zaczynamy od obrazu, a dopiero potem wyprowadzamy wzór.

1. Transformacje macierzowe — co macierze robią

Macierz 2×2 robi dokładnie jedną rzecz: przekształca wektory na płaszczyźnie. Może rozciągać, obracać, odbijać, ścinać lub zapadać. Wyznacznik — pojedyncza liczba — mówi, jak zmieniają się pola powierzchni. Gdy det = 0, macierz zapada płaszczyznę do linii lub punktu; informacja jest wtedy trwale utracona.

Transformacja:  T(v) = Mv
  M = [[a, b], [c, d]]   odwzorowuje kwadrat jednostkowy na równoległobok.

Wyznacznik:  det(M) = ad − bc
  |det| = współczynnik skalowania pola
  det > 0: orientacja zachowana  |  det < 0: odbicie
  det = 0: macierz osobliwa (rząd niepełny) — odwzorowuje w niższy wymiar

Wektory własne:  Mv = λv
  (M − λI)v = 0   ⟹  det(M − λI) = 0  (wielomian charakterystyczny)
  Wartości własne λ₁,λ₂ mówią o współczynnikach rozciągania w specjalnych kierunkach.

Dla 2×2:  λ = (tr(M)/2) ± √((tr(M)/2)² − det(M))
  Rzeczywiste wartości własne: M ma ortogonalną bazę własną (M symetryczna)
  Zespolone wartości własne: M ma składową rotacji

Ślad = λ₁ + λ₂ = a + d
  (suma wartości własnych = suma elementów przekątnej)

Twierdzenie spektralne (M symetryczna, M = Mᵀ):
  M = QΛQᵀ,   Q ortogonalna,  Λ = diag(λ₁, λ₂)
  Każdą macierz symetryczną można zdiagonalizować — podstawa PCA.

Symulacja transformacji macierzowych uruchamia 8 gotowych przekształceń — obrót, skalowanie, ścinanie, odbicie, rzut — oraz tryb swobodnej edycji. Przełącz okrąg jednostkowy, strzałki wektorów własnych i deformację siatki, by zobaczyć wyznacznik jako pole powierzchni, a ślad jako sumę rozciągnięć. Zwróć uwagę na przypadki tożsamości i macierzy nilpotentnej: jedna nic nie robi, druga zapada w dwóch krokach.

2. Regresja liniowa — dopasowywanie linii do zaszumionych danych

Metoda najmniejszych kwadratów (OLS) to najprostszy model w statystyce i jeden z najczęściej stosowanych. Mając n punktów danych, szuka się linii minimalizującej sumę kwadratów odległości pionowych. Rozwiązanie ma piękną postać zamkniętą — i jest szczególnym przypadkiem równań normalnych w algebrze liniowej.

Cel:  minimalizuj SSE = Σ(yᵢ − ŷᵢ)² = ‖y − Xβ‖²

Równania normalne:  (XᵀX)β = Xᵀy
  Rozwiązanie:  β̂ = (XᵀX)⁻¹Xᵀy  (jednoznaczne, jeśli X ma pełny rząd kolumnowy)

Dla regresji prostej (jedna zmienna objaśniająca + wyraz wolny):
  β̂₁ = Σ(xᵢ−x̄)(yᵢ−ȳ) / Σ(xᵢ−x̄)²  = Cov(x,y)/Var(x)
  β̂₀ = ȳ − β̂₁x̄

Współczynnik determinacji:
  R² = 1 − SSE/SST = 1 − Σ(yᵢ−ŷᵢ)²/Σ(yᵢ−ȳ)²
  R² = r²  (dla regresji prostej, r = współczynnik korelacji Pearsona)

Korelacja Pearsona r:  r = Cov(x,y)/[σₓ σᵧ] ∈ [−1, 1]
  |r| = 0: brak zależności liniowej
  |r| = 1: idealna zależność liniowa

Interpretacja geometryczna:
  β̂ = (XᵀX)⁻¹Xᵀy to rzut ortogonalny y na col(X).
  Wektor reszt e = y − Xβ̂ ⊥ col(X)  →  Xᵀe = 0  (zawsze!)

Twierdzenie Gaussa-Markowa:
  OLS jest BLUE (Best Linear Unbiased Estimator) gdy
  błędy są nieskorelowane, homoskedastyczne, o wartości oczekiwanej zero.
  Naruszenie → heteroskedastyczność, użyj WLS lub GLS.

Kliknij na wykres punktowy, aby dodać punkty i obserwować natychmiastową aktualizację linii regresji. Linie reszt (odległości pionowe) pokazują dokładnie, co minimalizuje OLS. Przeciągnij gotowy scenariusz z wartością odstającą, by zobaczyć dźwignię — jak jeden ekstremalny punkt może obrócić całą linię — i zrozumieć, dlaczego samo R² nie wystarcza do oceny modelu (kwartet Anscombe'a).

3. Wnioskowanie bayesowskie — aktualizowanie przekonań na podstawie danych

Twierdzenie Bayesa mieści się w kilku słowach: posterior to prior razy wiarygodność, znormalizowane. Tam, gdzie statystyka częstościowa pyta „Jakie jest prawdopodobieństwo tych danych przy założeniu hipotezy zerowej?”, statystyka bayesowska pyta „Jakie jest prawdopodobieństwo hipotezy przy tych danych?” — znacznie bardziej naturalne pytanie w większości problemów naukowych.

Twierdzenie Bayesa:  P(H|E) = P(E|H)·P(H) / P(E)
  H = hipoteza,  E = dowód (zaobserwowane dane)
  P(H)   = prawdopodobieństwo a priori (prior)
  P(E|H) = wiarygodność (likelihood)
  P(E)   = prawdopodobieństwo brzegowe (stała normalizująca)
  P(H|E) = prawdopodobieństwo a posteriori (posterior)

Postać szans (unika obliczania P(E)):
  Szanse a posteriori = Szanse a priori × Współczynnik wiarygodności (czynnik Bayesa)
  BF = P(E|H₁)/P(E|H₀)   (BF > 10 → „silny” dowód)

Priory sprzężone (postać zamknięta posteriora):
  Beta-dwumianowy:
    Prior: θ ~ Beta(α, β)
    Wiarygodność: k sukcesów w n próbach → rozkład dwumianowy(n, θ)
    Posterior: θ | dane ~ Beta(α+k, β+n−k)

Przykład testu medycznego — znaczenie częstości bazowej:
  Częstość choroby P(D) = 0,001 (1 na 1000)
  Czułość testu P(+|D) = 0,99
  Swoistość testu P(−|¬D) = 0,99
  P(D|+) = 0,99×0,001 / [0,99×0,001 + 0,01×0,999] ≈ 0,09
  → 91% wyników fałszywie dodatnich mimo testu o dokładności 99%!

Symulacja odtwarza scenariusz testu medycznego, pozwalając na dostosowanie częstości choroby, czułości i swoistości. Wynik — że pozytywny test na rzadką chorobę zwykle jest błędny — to najważniejszy i najbardziej sprzeczny z intuicją rezultat w probabilistyce stosowanej. Przełącz na aktualizację sekwencyjną i dodawaj obserwacje jedna po drugiej, by zobaczyć, jak posterior zaostrza się z płaskiego prioru.

4. Klastrowanie K-Means — podział bez etykiet

Nadzorowane uczenie maszynowe wymaga danych z etykietami — kosztownych w pozyskaniu. Klastrowanie jest nienadzorowane: mając nieoznaczone punkty, szuka naturalnych grupowań. K-Means to najprostszy i najczęściej stosowany algorytm klastrowania, na przemian przypisujący punkty do najbliższego centroidu i przesuwający każdy centroid do średniej jego klastra.

Algorytm Lloyda (K-Means):
  Zainicjuj k centroidów μ₁,…,μₖ  (losowo lub K-Means++)
  Powtarzaj aż do zbieżności:
    1. Przypisz:   zᵢ = argmin_j ‖xᵢ − μⱼ‖²  (najbliższy centroid)
    2. Zaktualizuj:   μⱼ = (1/|Cⱼ|) Σ_{i∈Cⱼ} xᵢ  (średnia klastra)
  Funkcja celu:  J = Σᵢ ‖xᵢ − μ_{z_i}‖²  (suma kwadratów odległości wewnątrz klastra)

  Zbieżność zagwarantowana (J maleje monotonicznie)
  Ale nie jest to optimum globalne — czuła na inicjalizację!

K-Means++ (Arthur i Vassilvitskii, 2007):
  1. Wybierz μ₁ losowo z rozkładu jednostajnego.
  2. Dla j = 2,…,k:
     Wybierz μⱼ z prawdopodobieństwem ∝ D(x)² = min_i ‖x − μᵢ‖²
  Oczekiwane J ≤ 8(ln k + 2) · J*   (gwarancja aproksymacji)

Wybór k — metoda łokcia:
  Narysuj wykres J od k.  Szukaj „łokcia”, gdzie ∂J/∂k gwałtownie maleje.

Współczynnik sylwetki:   s(i) = (b(i) − a(i)) / max(a(i), b(i))
  a(i) = średnia odległość wewnątrzklastrowa
  b(i) = średnia odległość do najbliższego innego klastra
  s ∈ [−1, 1]:  1 = idealne klastrowanie, 0 = granica, −1 = błędny klaster

Diagram Voronoi:   każda komórka = {x : ‖x−μⱼ‖ ≤ ‖x−μₗ‖ ∀l≠j}
  Granice K-Means są zawsze liniowe (teselacja Voronoi).

Symulacja renderuje pełny podział Voronoi piksel po pikselu przy użyciu ImageData — każdy piksel jest kolorowany zgodnie z najbliższym centroidem. Przejdź krok po kroku przez algorytm Lloyda, by zobaczyć zbieżność, przełącz na K-Means++, aby porównać jakość inicjalizacji, i skorzystaj z wykresu łokciowego, by odkryć „właściwe” k dla trzech różnych syntetycznych zbiorów danych.

5. Metoda gradientu prostego — nawigacja po krajobrazie funkcji straty

Niemal każdy współczesny model uczenia maszynowego jest trenowany metodą gradientu prostego: wielokrotnie oblicza się gradient funkcji straty względem parametrów, a następnie wykonuje mały krok w dół zbocza. Warianty — SGD, Momentum, RMSprop, Adam — różnią się jedynie sposobem wykorzystania historii gradientu do adaptacji wielkości kroku.

Zwykły SGD:
  θ ← θ − α∇L(θ)   (α = współczynnik uczenia)
  Zbieżność dla wypukłego L: O(1/√T) dla niegładkiego, O(1/T) dla gładkiego

Momentum (Polyak, 1964):
  v ← βv − α∇L(θ)
  θ ← θ + v
  β ≈ 0,9.  Przyspiesza wzdłuż spójnych kierunków gradientu.
  Wariant Nesterova: oblicza gradient w θ + βv (wyprzedzenie)

RMSprop (Hinton, 2012, niepublikowane):
  s ← ρs + (1−ρ)(∇L)²
  θ ← θ − (α/√(s+ε))∇L
  Adaptuje współczynnik uczenia dla każdego parametru — przydatne dla celów niestacjonarnych

Adam (Kingma i Ba, 2014):
  m ← β₁m + (1−β₁)∇L        (pierwszy moment)
  v ← β₂v + (1−β₂)(∇L)²     (drugi moment)
  m̂ = m/(1−β₁ᵗ),  v̂ = v/(1−β₂ᵗ)   (korekta obciążenia)
  θ ← θ − α·m̂/(√v̂+ε)
  Wartości domyślne: α=1e-3, β₁=0,9, β₂=0,999, ε=1e-8
  Zbieżność ≈ O(1/√T) dla celów niewypukłych

Funkcja Rosenbrocka (bananowa):
  f(x,y) = (1−x)² + 100(y−x²)²
  Globalne minimum w (1,1). Notorycznie trudna: strome zbocze z łagodnym nachyleniem.
  Adam ≈ 800 iteracji do zbieżności wobec SGD ≈ 10 000+

Symulacja renderuje krajobraz funkcji straty 60×60 jako kolorowaną powierzchnię 3D (algorytm malarza) i śledzi trajektorię optymalizatora w czasie rzeczywistym. Uruchom wszystkie cztery optymalizatory jednocześnie na bananowej funkcji Rosenbrocka i wielominimalnej funkcji Himmelblaua. Wizualna różnica między zygzakiem SGD, przeregulowaniem Momentum i czystym zejściem Adama wyjaśnia, dlaczego Adam stał się uniwersalnym domyślnym wyborem w głębokim uczeniu.

6. Konwolucyjne sieci neuronowe — ekstrakcja cech z założenia

W pełni połączona sieć neuronowa zastosowana do obrazu 224×224 wymagałaby ponad 50 000 wag na neuron — obliczeniowo niewykonalne. Kluczowa obserwacja stojąca za CNN: sąsiadujące piksele są skorelowane, a ta sama cecha (krawędź, tekstura) może pojawić się w dowolnym miejscu obrazu. Zastąpienie pełnych połączeń lokalnymi filtrami i współdzieleniem wag redukuje liczbę parametrów o rzędy wielkości, jednocześnie wbudowując niezmienność względem przesunięcia.

Dyskretna konwolucja 2D:
  (f * g)[i,j] = Σₘ Σₙ f[m,n]·g[i−m, j−n]
  W CNN g jest wyuczonym filtrem (jądrem) o rozmiarze k×k.

Rozmiar wyjściowy (bez wypełnienia):
  H_out = (H_in − k) / stride + 1

Obliczanie mapy cech (dla jednego filtra W, biasu b):
  Z[i,j] = (X * W)[i,j] + b  (przed aktywacją)
  A[i,j] = ReLU(Z[i,j]) = max(0, Z[i,j])

Popularne filtry (zaprojektowane ręcznie):
  Krawędź pozioma:  [[-1,-1,-1],[0,0,0],[1,1,1]]
  Krawędź pionowa:    [[-1,0,1],[-1,0,1],[-1,0,1]]
  Wyostrzenie:          [[0,-1,0],[-1,5,-1],[0,-1,0]]
  Rozmycie Gaussa:    1/16 · [[1,2,1],[2,4,2],[1,2,1]]

Max pooling (2×2, stride 2):
  Zmniejsza wymiary przestrzenne o połowę, zachowuje maksymalną aktywację.
  Zapewnia niezmienność względem drobnych przesunięć.

Propagacja wsteczna przez konwolucję:
  ∂L/∂W = wejście * ∂L/∂Z  (korelacja, nie konwolucja)
  ∂L/∂X = pełna konwolucja ∂L/∂Z z odwróconym W

Liczba parametrów (jedna warstwa konwolucyjna, k=3, C_in kanałów, C_out filtrów):
  (k² · C_in + 1) · C_out  =  (9 · C_in + 1) · C_out
  vs. w pełni połączona: H·W·C_in · H·W·C_out  (rzędy wielkości więcej)

Wybierz jeden z pięciu obrazów wejściowych 8×8 i przełączaj między czterema ręcznie zaprojektowanymi filtrami. Nawiguj po zakładkach warstw — Wejście → Konwolucja → ReLU → MaxPool → głębsze warstwy — i obserwuj, jak animowana pozycja skanowania wydobywa cechy krok po kroku. Panel architektury śledzi wymiary przestrzenne i liczbę parametrów, czyniąc obliczeniową przewagę współdzielenia wag konkretną i namacalną.

Jedność algebry liniowej i data science

Wszystkie sześć symulacji to różne oblicza tej samej struktury geometrycznej. Regresja liniowa to rzut. K-Means dzieli przestrzeń diagramem Voronoi, którego granice są symetralnymi zdefiniowanymi przez metrykę euklidesową. Metoda gradientu prostego podąża wzdłuż gradientu formy kwadratowej — która w pobliżu minimum jest zawsze zdominowana przez macierz Hessego i jej wartości własne. Warstwa konwolucyjna to ustrukturyzowane rzadkie mnożenie macierzowe. Nawet twierdzenie Bayesa, zastosowane do gaussowskiego prioru i gaussowskiej wiarygodności, daje gaussowski posterior, którego średnia jest ważoną średnią — operacją liniową.

Seria o uczeniu maszynowym jest kontynuowana w tematach drzew decyzyjnych (podziały nieliniowe) i uczenia ze wzmocnieniem w Nauce #9 — Uczenie maszynowe, a podstawy statystyczne omówiono w Nauce #15 — Prawdopodobieństwo i statystyka. Stronę przetwarzania sygnałów w data science — FFT, filtry, spektrogramy — omawia Nauka #14 — Przetwarzanie sygnałów.