Lineáris egyenletek megoldása ebben a kurzusban:

A(z) Lineáris egyenletek megoldása fogalomról

Emlékszel, amikor először tanultál algebrát, még a középiskolában?

Általában adtak egy egyenletrendszert, és azt kérték, hogy oldjuk meg x-re és y-ra.

4x - 3y = -2
2x + 7y = 16

Két egyenlet, két ismeretlen: egy kis behelyettesítés, és máris láthatod, hogy x = 1, y = 2.

Nézd meg újra az előző példa bal oldalát. Úgy néz ki, mint egy [4 -3; 2 7] mátrix, amely megszoroz egy ismeretlen [x, y] vektort, és eredményül a [-2, 16] vektort adja.

Lineáris algebrai jelöléssel: A x = b, követve a szokásos konvenciót, miszerint a mátrixokat nagybetűk, a vektorokat kisbetűk jelölik.

A x = 0 és a nulltér

Egy fontos különleges eset, amikor minden egyenletünk jobb oldalán nulla áll.

4x - 3y = 0
2x + 7y = 0

A fenti példában az egyetlen megoldás az, amikor x = y = 0, ami általában nem túl érdekes.

A A x = 0 egyenletnek csak akkor vannak nemtriviális megoldásai, ha A determinánsa nulla.

julia> using LinearAlgebra

julia> A = [4 -3; 2 7]
2×2 Matrix{Int64}:
 4  -3
 2   7

julia> det(A)  # not zero!
34.0

Vegyük inkább ezeket az egyenleteket:

4x - 3y = 0
8x + 6y = 0

A második egyenlet csak a kétszerese az elsőnek, nem ad új információt. Bármely olyan érték megoldás lesz, amelyre x = 0.75y.

Ezek az (x, y) értékek egy egyenes mentén fekszenek a kétdimenziós térben. Végtelen sok megoldás van, és ezek alkotják a mátrix nullterét.

Ebben az esetben a mátrix sorai lineárisan nem függetlenek, és a mátrix rangja kisebb, mint a sorok vagy oszlopok száma.

julia> A = [4 -3; 8 -6]
2×2 Matrix{Int64}:
 4  -3
 8  -6

julia> det(A)
0.0

julia> size(A)  # total number of rows and columns
(2, 2)

julia> rank(A)  # number of linearly-independent rows (or columns)
1

julia> nullspace(A)
2×1 Matrix{Float64}:
 -0.6
 -0.8

julia> nullspace(A) |> norm  # Julia returns unit vectors when possible
1.0

A nullspace() függvény egyetlen pontot ad vissza, amelyet normalizáltak egységvektorrá, de ennek a vektornak bármely skalárszorosa is a nulltérben van.

Az A x = b megoldása

Most pedig képzeld el, hogy 1000 egyenleted van 1000 ismeretlennel? Ha ez hülyeségnek hangzik, ne feledd, hogy a vezeték nélküli digitális érzékelők ma már olcsók és sokoldalúak (mérnek nyúlást, szélsebességet, háromtengelyes gyorsulást, bármit...). Teljesen ésszerű, ha 1000 ilyen érzékelő figyel egy modern szerkezetet, például egy függőhidat. Kell lennie egy vezérlőrendszernek, amely értelmezi az adatfolyamot, és készen áll riasztást adni, ha a dolgok aggasztóvá válnak.

Megint csak A x = b, ahol A (a mérnöki tervből) és b (az érzékelők mért értékei) adott, de meg kell keresnünk x-et.

Aki még nem ismeri ezt, annak első gondolata, hogy valahogy »mindkét oldalt elosztja A-val«, hogy átvigye a jobb oldalra.

Valójában A-nak van inverze (a legtöbb négyzetes mátrixnak van, de nem mindegyiknek: a determináns nem lehet nulla), és a számítás működik (lassan!).

julia> A = [4 -3; 2 7]
2×2 Matrix{Int64}:
 4  -3
 2   7

julia> b = [-2, 16]
2-element Vector{Int64}:
 -2
 16

# the inverse matrix
julia> inv(A)
2×2 Matrix{Float64}:
  0.205882   0.0882353
 -0.0588235  0.117647

julia> inv(A) * b
2-element Vector{Float64}:
 0.9999999999999999
 2.0

Sajnos az inverz kiszámítása lassú a kis mátrixoknál, a nagyoknál pedig dermesztően lassú. Nem éppen ideális, ha a híd összeomlik a számítás közben!

Szerencsére más algoritmusok drámaian gyorsabbak (ebben az esetben a Gauss-elimináció).

A Julia (a Matlabot utánozva) egyszerűen egy visszaperjelet használ a megoldóhoz (technikailag »bal osztás«).

julia> x = A \ b
2-element Vector{Float64}:
 1.0
 2.0

Ez egy triviális példa, de van néhány részlet, amire érdemes figyelni.

• Az együtthatóknak oszlopokba kell rendeződniük, ezért a hiányzó tagok helyére szúrj be nullákat.
• Az egyenleteknek (és így a mátrix sorainak) lineárisan függetleneknek kell lenniük. Ha a 3. sor csak az 1. és 2. sor összege, nem ad új információt, és a feladat alulhatározott.

A rank() parancs megadja a lineárisan független sorok/oszlopok számát, ezért törekedj arra, hogy ez megegyezzen az ismeretlen változók számával.

julia> rank(A)
2

Téglalap alakú mátrixok

Valamikor kiskamaszként valószínűleg azt tanultad, hogy N ismeretlen megoldásához N egyenletből álló rendszer kell.

Ahogy a matematikában a legtöbb dologgal, a valóság itt is egy kicsit árnyaltabb (és nem csupán a lineáris függetlenség részletével, amelyet az előző szakaszban említettünk).

„Túl kevés” sor

N-1 egyenlet (az A mátrix N-1 sora) esetén a feladat alulhatározott. Feladjuk kétségbeesésünkben?

Attól függ!

Az N-1 egyenlet még mindig sok információt tartalmaz. Geometriai értelemben N egyenlet kell ahhoz, hogy a megoldást egy pontként határozzuk meg az N-dimenziós térben, de N-1 megmondja, hogy a megoldásnak valahol egy egyenesen kell lennie.

Így végtelen sok megoldásunk marad, de a paramétertérnek is végtelen sok olyan része, ahol nincs megoldás.

Átmegy a megoldásaid egyenese a paramétertérnek olyan tartományán, amely miatt egy jó hírű mérnöknek aggódnia kellene? Vagy amely indokolttá tehetné, hogy a szerkezetet lezárják a nyilvánosság elől (a függőhidak élénkek tudnak lenni erős oldalszélben)? Lehet, hogy megér néhány további számítást, hogy ezt kiderítsük!

julia> A3 = [4 -3 1; 2 2 2; 3 -4 -3]
3×3 Matrix{Int64}:
 4  -3   1
 2   2   2
 3  -4  -3

julia> b3 = [1, 12, -14]
3-element Vector{Int64}:
   1
  12
 -14

# fully-determined solution
julia> A3 \ b3
3-element Vector{Float64}:
 1.0
 2.0
 3.0

# remove row 3 from A3 and B3
julia> A2 = A3[1:2, :]
2×3 Matrix{Int64}:
 4  -3  1
 2   2  2

julia> b2 = b3[1:2]
2-element Vector{Int64}:
  1
 12

# an under-determined solution
julia> z1 = A2 \ b2
3-element Vector{Float64}:
 1.594594594594594
 2.445945945945947
 1.9594594594594592

A visszaperjeles megoldó ad nekünk egy megoldást! Nem hivatalosan ez láthatóan a legkisebb normájú megoldás (geometriailag a legközelebb az origóhoz).

A többi megoldás eléréséhez az A2 nullspace függvényére van szükségünk.

Ha a nulltér bármely skalárszorosát hozzáadjuk az eredeti megoldáshoz, egy másik érvényes megoldást kapunk.

# get the nummspace of A2
julia> N = nullspace(A2)
3×1 Matrix{Float64}:
 -0.46499055497527725
 -0.3487429162314578
  0.813733471206735

# add a random multiple of the nullspace to the earlier solution
julia> z = z1 + rand() * N
3×1 Matrix{Float64}:
 1.274331860650258
 2.205748895487695
 2.5199192438620472

# our random z is a valid solution, recovering b2 = [1, 12]
julia> A2 * z
2×1 Matrix{Float64}:
  0.9999999999999947
 12.0

Hasonlóképpen az N-2 egyenlet egy síkra korlátozza a (végtelen sok) megoldást, és ugyanezek az elvek érvényesek.

„Túl sok” sor

Az ellenkező helyzet akkor áll elő, amikor több egyenlet van, mint ismeretlen, mégis valahogy lineárisan függetlenek maradnak.

A valódi mérnöki gyakorlatban ez teljesen normális, sőt jó dolognak számít!

Az érzékelők pontossága véges, a b vektor pontossága véges, és zaj van a számításban. Ilyenkor a megoldás a zajos adatokra illesztett legkisebb négyzetek szerinti illesztés.

A legegyszerűbb technika a mátrix pseudoinverse függvényét használja, amelyet a Julia a pinv() függvényként valósít meg.

Ez nagyon hasonlóan használható, mint egy (nem szinguláris) négyzetes mátrix inverze: a változók legkisebb négyzetek szerinti becslését adja pontos megoldás helyett.

# create 5-row A and b for 2 variables
julia> A5 = [4 -3; 2 7; -1 2; 1 -1; 3 1]
5×2 Matrix{Int64}:
  4  -3
  2   7
 -1   2
  1  -1
  3   1

julia> rank(A5)
2

# add some random-normal noise to b5
julia> b5n = [-2, 16, 3, -1, 5] + randn(5) * 0.1
5-element Vector{Float64}:
 -1.9337349223581917
 16.024913008059457
  3.0087646177373992
 -0.8675748721176717
  4.914568047308805

# solve to get x ≈ 1, y ≈ 2, using the pseudoinverse
julia> pinv(A5) * b5n
2-element Vector{Float64}:
 1.006117942769543
 1.9962973764508272

Sajátértékek és sajátvektorok

Az előző szakasz a A x = b megoldásait írta le, ami megfelel a jól ismert algebrának.

Ez a szakasz a A x = λ x megoldásairól szól, ahol λ egy skalár.

Ez egy jellegzetesen lineáris algebrai fogalom, de geometriailag is értelmezhetjük:

Az x és λ megfelelő értékei esetén az A négyzetes mátrix a nullvektortól különböző x vektort λ-szorosára nyújtja hosszában, miközben az irányát nem változtatja meg (legfeljebb egy negatív λ fordítja meg).

Ez elsőre nagyon speciálisnak hangzik, de kiderül, hogy nevetségesen hasznos!

Terminológia: A λ érvényes értékei az A sajátértékei, a megfelelő x értékek pedig az A sajátvektorai. Sajnos együtt kell élnünk azokkal a szavakkal, amelyek félúton nyelvet váltanak.

Az iskolában általában megtanítják, hogyan számítsuk ki kézzel a 2×2-es mátrixok sajátértékeit/sajátvektorait, de számítógéppel sokkal könnyebb (bár még mindig meglehetősen lassú, és általános esetben egy n×n-es mátrixnál O(n^3)-ként skálázódik).

julia> A = rand(-9:9, 2, 2)
2×2 Matrix{Int64}:
 9  5
 3  2

julia> F = eigen(A)
Eigen{Float64, Float64, Matrix{Float64}, Vector{Float64}}
values:
2-element Vector{Float64}:
  0.27984674554472466
 10.720153254455274
vectors:
2×2 Matrix{Float64}:
 -0.497417  0.945605
  0.867511  0.325317

# eigenvalues
julia> F.values
2-element Vector{Float64}:
  0.27984674554472466
 10.720153254455274

# each column is an eigenvector, normalized to a unit vector
julia> F.vectors
2×2 Matrix{Float64}:
 -0.497417  0.945605
  0.867511  0.325317

Általánosságban egy n×n-es mátrixnak n sajátértéke van, bár nem mindig különböző értékek. Tekintsd őket egy n-edfokú polinom gyökereiként (ezt karakterisztikus polinomnak nevezik), amelyek ismétlődhetnek, és gyakran komplexek, még valós értékű mátrix esetében is.

Minden sajátvektor egy irányt képvisel, és bármely skalárszorosa is érvényes sajátvektor. A további számítások megkönnyítése érdekében a Julia 1 normájú egységvektorokat ad vissza.

Alkalmazások

Annak megmagyarázása, hogy miért fontosak a sajátvektorok, ideális esetben egy 500 oldalas tankönyv feladata, nem néhány rövid bekezdésé. Ez a téma a modern alkalmazott matematika olyan sok területét áthatja.

A legmagasabb szinten a sajátvektorok egy adathalmaz »legfontosabb« tengelyeit (irányait) képviselik. A sajátértékek az egyes tengelyek »relatív fontosságát« jelzik (néhány, a bemenetek normalizálására vonatkozó feltevés függvényében).

Hogy ez valójában mit jelent, az az alkalmazástól függ.

Főkomponens-analízis

Egy tipikus adattudományi feladatban lehet akár 100 vagy több »jellemzőnk«, amelyeket adatoszlopokként tárolunk. Szinte elkerülhetetlenül lesz zaj, redundancia és nem kívánatos korreláció.

Dimenziócsökkentésre van szükségünk, és a PCA az egyik módja annak, hogy rendet teremtsünk a káoszban:

• Számítsd ki az adatok kovarianciamátrixát.
• Most már van egy négyzetes mátrixunk, úgyhogy következőként számítsuk ki a sajátértékeket és sajátvektorokat.
• Rendezd a sajátvektorokat sajátértékeik csökkenő sorrendjébe (pontosabban az |λ| abszolút érték szerint).
• A legfelső k sajátvektor most a principal components (főkomponensek), ahol a k lényegesen kisebb az eredeti jellemzők számánál.
• Vetítsd a teljes adathalmazt ezekre a k tengelyekre, és kezdj érdekes mintázatokat keresni.

Ebben a formájában a PCA-t olyan különböző csoportok használták, mint az orvostudósok, akik a gyógyszertervezésben molekulák tulajdonságait próbálják javítani, és a politikai kampányolók, akik közvélemény-kutatási adatokból próbálják megérteni a választói preferenciákat. Valószínűleg marketingcsoportok is, amelyek dolgokat próbálnak eladni neked, de bármely technológia használható jóra és rosszra is.

A PCA nem része egy minimális Julia-telepítésnek, de (az Exercism-ön kívül) a MultivariateStats csomag mindent tartalmaz, amire szükséged van.

Képfeldolgozás

A PCA ötletét továbbgondolva a digitális képek csupán pixelértékek mátrixai, és ki tudjuk számítani a főkomponenseiket.

Ezt évtizedek óta használják a képtömörítésben, ahol a PCA útmutatóként szolgál ahhoz, hogy mit érdemes leginkább megőrizni a fájlméret csökkentésekor.

A PCA egyre inkább létfontosságú része az olyan kép-osztályozóknak, mint az arcfelismerés. Amikor legközelebb átsétálsz egy repülőtéren, a Nagy Testvér nem csupán figyel: lineáris algebrát is használ, hogy megértse, amit lát!

Gépészmérnöki tudomány

Kevésbé vitatható módon egy mechanikai alkatrész főtengelyei a tehetetlenségi nyomaték tenzorának I sajátvektorai (egy mátrix, csak kicsit más néven).

Az autód minden kerekén valószínűleg van egy vagy több kiegyensúlyozó súly, amelyeket úgy állítanak be, hogy nullázzák e tenzor I nem átlós elemeit. A műhely szerelője kétségtelenül elkerüli a matekozást (ellentétben a repülőgép-tervezőkkel és a rakétamérnökökkel), de a »billegés« csupán hétköznapi szó a tenzor kereszttagjaira, amelyek ebben az esetben kényelmetlenebbé tennék az utazást, és növelnék a csapágyak mechanikai kopását.