Επίλυση γραμμικών εξισώσεων σε

Σχετικά με την έννοια Επίλυση γραμμικών εξισώσεων

Θυμάσαι όταν έμαθες για πρώτη φορά άλγεβρα, στο Λύκειο;

Συνήθως μας έδιναν ένα σύστημα εξισώσεων και μας ζητούσαν να λύσουμε ως προς x και y.

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

Δύο εξισώσεις, δύο άγνωστοι: λίγη αντικατάσταση και σύντομα βλέπεις ότι x = 1, y = 2.

Κοίτα ξανά το αριστερό μέλος του προηγούμενου παραδείγματος. Μοιάζει με έναν πίνακα [4 -3; 2 7] που πολλαπλασιάζει έναν άγνωστο διάνυσμα [x, y] για να δώσει το διάνυσμα [-2, 16].

Στον συμβολισμό της Γραμμικής Άλγεβρας: A x = b, ακολουθώντας τη συνήθη σύμβαση των κεφαλαίων γραμμάτων για τους πίνακες και των πεζών για τα διανύσματα.

A x = 0 και ο μηδενόχωρος

Μια σημαντική ειδική περίπτωση είναι όταν όλες οι εξισώσεις μας έχουν μηδέν στο δεξί μέλος.

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

Στο παραπάνω παράδειγμα, η μόνη λύση είναι όταν x = y = 0, που συνήθως δεν παρουσιάζει ενδιαφέρον.

Μη τετριμμένες λύσεις για A x = 0 υπάρχουν μόνο όταν ο A έχει μηδενική ορίζουσα.

julia> using LinearAlgebra

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

julia> det(A)  # not zero!
34.0

Ας δούμε, αντίθετα, αυτές τις εξισώσεις:

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

Η δεύτερη εξίσωση είναι απλώς η διπλάσια της πρώτης, χωρίς να προσθέτει νέα πληροφορία. Οποιεσδήποτε τιμές όπου x = 0.75y θα αποτελούν λύση.

Αυτές οι τιμές (x, y) βρίσκονται πάνω σε μια ευθεία στον δισδιάστατο χώρο. Υπάρχει άπειρος αριθμός λύσεων και σχηματίζουν τον μηδενόχωρο του πίνακα.

Οι γραμμές του πίνακα δεν είναι γραμμικά ανεξάρτητες σε αυτή την περίπτωση, και ο βαθμός του πίνακα είναι μικρότερος από τον αριθμό γραμμών ή στηλών.

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

Η συνάρτηση nullspace() επιστρέφει ένα μόνο σημείο, επιλεγμένο ώστε να είναι κανονικοποιημένο σε μοναδιαίο διάνυσμα, αλλά οποιοδήποτε βαθμωτό πολλαπλάσιο αυτού του διανύσματος ανήκει επίσης στον μηδενόχωρο.

Επίλυση του A x = b

Τώρα φαντάσου ότι έχεις 1000 εξισώσεις με 1000 αγνώστους; Αν αυτό ακούγεται ανόητο, θυμήσου ότι οι ασύρματοι ψηφιακοί αισθητήρες είναι πλέον φθηνοί και ευέλικτοι (μετρούν παραμόρφωση, ταχύτητα ανέμου, επιτάχυνση σε 3 άξονες, ό,τι να 'ναι...). Το να τα έχεις 1000 από αυτά να παρακολουθούν μια σύγχρονη κατασκευή όπως μια κρεμαστή γέφυρα είναι απολύτως λογικό. Χρειάζεται ένα σύστημα ελέγχου που να ερμηνεύει τη ροή δεδομένων, έτοιμο να σημάνει συναγερμό αν τα πράγματα γίνουν ανησυχητικά.

Και πάλι, A x = b, όπου έχουμε το A (από τη μελέτη του μηχανικού) και το b (μετρημένες τιμές από τους αισθητήρες), αλλά πρέπει να βρούμε το x.

Μια πρώτη σκέψη, για όποιον δεν είναι εξοικειωμένος με αυτά, είναι να "διαιρέσει με το A" για να το μετακινήσει στο δεξί μέλος.

Στην πραγματικότητα, το A έχει αντίστροφο (τα περισσότερα, αλλά όχι όλα τα τετραγωνικά μητρώα έχουν: η ορίζουσα πρέπει να είναι μη μηδενική), και ο υπολογισμός λειτουργεί (αργά!).

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

Δυστυχώς, ο υπολογισμός του αντιστρόφου είναι αργός για μικρούς πίνακες και παγερά αργός για μεγάλους. Όχι και το ιδανικό, αν η γέφυρά σου καταρρεύσει κατά τη διάρκεια του υπολογισμού!

Ευτυχώς, άλλοι αλγόριθμοι είναι δραματικά ταχύτεροι (σε αυτή την περίπτωση, η απαλοιφή Gauss).

Η Julia (αντιγράφοντας τη Matlab) χρησιμοποιεί απλώς μια ανάστροφη κάθετο για τον επιλυτή (τεχνικά, "αριστερή διαίρεση").

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

Αυτό είναι ένα απλό παράδειγμα, αλλά υπάρχουν μερικές λεπτομέρειες που πρέπει να προσέξεις.

• Οι συντελεστές πρέπει να ευθυγραμμίζονται σε στήλες, οπότε βάλε μηδενικά όπου χρειάζεται για τους όρους που λείπουν.
• Οι εξισώσεις (και επομένως οι γραμμές του πίνακα) πρέπει να είναι γραμμικά ανεξάρτητες. Αν η γραμμή 3 είναι απλώς το άθροισμα των γραμμών 1 και 2, δεν προσθέτει νέα πληροφορία και το πρόβλημα είναι υποπροσδιορισμένο.

Η εντολή rank() θα σου δώσει τον αριθμό των γραμμικά ανεξάρτητων γραμμών/στηλών, οπότε στόχευσε αυτό να ισούται με τον αριθμό των άγνωστων μεταβλητών.

julia> rank(A)
2

Ορθογώνιοι πίνακες

Κάποια στιγμή πριν από την εφηβεία, σου έμαθαν μάλλον ότι χρειάζεσαι ένα σύστημα N εξισώσεων για να λύσεις N αγνώστους.

Όπως με τα περισσότερα πράγματα στα μαθηματικά, η πραγματικότητα είναι λίγο πιο σύνθετη (όχι μόνο τη λεπτομέρεια της γραμμικής ανεξαρτησίας που αναφέρθηκε στην προηγούμενη ενότητα).

"Πολύ λίγες" γραμμές

Με N-1 εξισώσεις (γραμμές στον πίνακα A), το πρόβλημα είναι υποπροσδιορισμένο. Τα παρατάμε από απόγνωση;

Εξαρτάται!

Οι N-1 εξισώσεις περιέχουν ακόμα πολλή πληροφορία. Γεωμετρικά, χρειαζόμαστε N για να προσδιορίσουμε τη λύση ως σημείο στον N-διάστατο χώρο, αλλά οι N-1 θα μας πουν ότι η λύση πρέπει να βρίσκεται κάπου σε μια ευθεία.

Αυτό μας αφήνει με άπειρο αριθμό λύσεων, αλλά και άπειρο χώρο παραμέτρων με καμία λύση.

Η ευθεία των λύσεών σου περνάει από μια περιοχή του χώρου παραμέτρων που θα έπρεπε να ανησυχήσει έναν αξιόπιστο μηχανικό; Ή που θα μπορούσε να δικαιολογήσει το κλείσιμο της κατασκευής για το κοινό (οι κρεμαστές γέφυρες μπορούν να κουνιούνται έντονα με δυνατούς πλάγιους ανέμους); Ίσως αξίζει να κάνεις μερικούς ακόμη υπολογισμούς για να το διαπιστώσεις!

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

Ο επιλυτής με την ανάστροφη κάθετο μας δίνει μια λύση! Ανεπίσημα, αυτή είναι προφανώς η λύση με τη μικρότερη νόρμα (πιο κοντά στην αρχή, γεωμετρικά).

Για να πάρουμε τις άλλες λύσεις, χρειαζόμαστε τον nullspace του A2.

Προσθέτοντας οποιοδήποτε βαθμωτό πολλαπλάσιο του μηδενόχωρου στην αρχική λύση θα δώσει άλλη μια έγκυρη λύση.

# 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

Παρομοίως, οι N-2 εξισώσεις θα περιορίσουν τις (άπειρες) λύσεις σε ένα επίπεδο, και οι ίδιες αρχές ισχύουν.

"Πάρα πολλές" γραμμές

Η αντίθετη κατάσταση προκύπτει όταν υπάρχουν περισσότερες εξισώσεις από αγνώστους, αλλά παρ' όλα αυτά είναι ακόμα γραμμικά ανεξάρτητες.

Στην πραγματική μηχανική, αυτό είναι απολύτως φυσιολογικό και θεωρείται καλό!

Οι αισθητήρες έχουν περιορισμένη ακρίβεια, το διάνυσμα b έχει περιορισμένη ακρίβεια, και υπάρχει θόρυβος στον υπολογισμό. Τώρα η λύση είναι μια προσαρμογή ελαχίστων τετραγώνων στα θορυβώδη δεδομένα.

Η απλούστερη τεχνική χρησιμοποιεί τον pseudoinverse του πίνακα, τον οποίο η Julia υλοποιεί ως τη συνάρτηση pinv().

Αυτός μπορεί να χρησιμοποιηθεί σχεδόν όπως ο αντίστροφος ενός (μη ιδιάζοντος) τετραγωνικού πίνακα, δίνοντας μια εκτίμηση ελαχίστων τετραγώνων των μεταβλητών αντί για ακριβή λύση.

# 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

Ιδιοτιμές και Ιδιοδιανύσματα

Η προηγούμενη ενότητα περιέγραψε λύσεις για A x = b, που ισοδυναμούν με τη γνωστή άλγεβρα.

Αυτή η ενότητα αφορά λύσεις της A x = λ x, όπου το λ είναι βαθμωτό μέγεθος.

Αυτή είναι μια ξεχωριστή έννοια της Γραμμικής Άλγεβρας, αλλά μπορούμε να την ερμηνεύσουμε γεωμετρικά:

Για κατάλληλες τιμές των x και λ, ο τετραγωνικός πίνακας A κλιμακώνει το μη μηδενικό διάνυσμα x κατά έναν παράγοντα λ σε μήκος, χωρίς να αλλάζει τη διεύθυνσή του (εκτός από την περίπτωση αρνητικού λ που την αντιστρέφει).

Αυτό ακούγεται εξειδικευμένο, αλλά αποδεικνύεται παραλόγως χρήσιμο!

Ορολογία: Οι έγκυρες τιμές του λ είναι οι ιδιοτιμές του A, και οι αντίστοιχες τιμές του x είναι τα ιδιοδιανύσματα του A. Δυστυχώς, πρέπει απλώς να συμβιβαστούμε με λέξεις που αλλάζουν γλώσσα στη μέση.

Στους μαθητές συνήθως διδάσκεται πώς να υπολογίζουν ιδιοτιμές/ιδιοδιανύσματα για πίνακες 2×2 με το χέρι, αλλά η χρήση υπολογιστή είναι πολύ πιο εύκολη (αν και ακόμα αρκετά αργή, με κλιμάκωση O(n^3) για έναν πίνακα n×n στη γενική περίπτωση).

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

Γενικά, ένας πίνακας n×n θα έχει n ιδιοτιμές, αν και όχι πάντα διακριτές τιμές. Σκέψου τις ως τις ρίζες ενός πολυωνύμου n-οστού βαθμού (που ονομάζεται χαρακτηριστικό πολυώνυμο), οι οποίες μπορεί να είναι επαναλαμβανόμενες και είναι συχνά μιγαδικές ακόμα και για έναν πραγματικό πίνακα.

Κάθε ιδιοδιάνυσμα αντιπροσωπεύει μια διεύθυνση, και οποιοδήποτε βαθμωτό πολλαπλάσιό του είναι επίσης έγκυρο ιδιοδιάνυσμα. Για ευκολία σε περαιτέρω υπολογισμούς, η Julia επιστρέφει μοναδιαία διανύσματα με νόρμα ίση με 1.

Εφαρμογές

Το να εξηγήσεις γιατί τα ιδιοδιανύσματα είναι σημαντικά είναι ιδανικά δουλειά για ένα βιβλίο 500 σελίδων, όχι για λίγες σύντομες παραγράφους. Αυτό το θέμα διαπερνά τόσο μεγάλο μέρος των σύγχρονων εφαρμοσμένων μαθηματικών.

Σε ένα γενικό επίπεδο, τα ιδιοδιανύσματα αντιπροσωπεύουν τους "πιο σημαντικούς" άξονες (διευθύνσεις) σε ένα σύνολο δεδομένων. Οι ιδιοτιμές υποδεικνύουν τη "σχετική σημασία" κάθε άξονα (με την επιφύλαξη ορισμένων παραδοχών για την κανονικοποίηση των εισόδων).

Το τι σημαίνει αυτό στην πραγματικότητα εξαρτάται από την εφαρμογή.

Ανάλυση Κυρίων Συνιστωσών

Σε ένα τυπικό πρόβλημα επιστήμης δεδομένων, μπορεί να έχουμε 100 ή περισσότερα "χαρακτηριστικά", αποθηκευμένα ως στήλες δεδομένων. Σχεδόν αναπόφευκτα, θα υπάρχει θόρυβος, πλεονασμός και ανεπιθύμητες συσχετίσεις.

Πρέπει να κάνουμε μείωση διαστάσεων, και η PCA είναι ένας τρόπος να βάλουμε τάξη στο χάος:

• Υπολόγισε τον πίνακα συνδιακύμανσης των δεδομένων.
• Τώρα έχουμε έναν τετραγωνικό πίνακα, οπότε στη συνέχεια υπολόγισε τις ιδιοτιμές και τα ιδιοδιανύσματα.
• Ταξινόμησε τα ιδιοδιανύσματα κατά φθίνουσα σειρά των ιδιοτιμών τους (πιο συγκεκριμένα, της απόλυτης τιμής |λ|).
• Τα πρώτα k ιδιοδιανύσματα είναι τώρα οι principal components σου, όπου το k είναι σημαντικά μικρότερο από τον αρχικό αριθμό χαρακτηριστικών.
• Προβάλλε το πλήρες σύνολο δεδομένων σε αυτούς τους k άξονες και άρχισε να ψάχνεις για ενδιαφέροντα μοτίβα.

Σε αυτή τη μορφή, η PCA έχει χρησιμοποιηθεί από ομάδες τόσο διαφορετικές όσο επιστήμονες της ιατρικής που προσπαθούν να βελτιώσουν μοριακές ιδιότητες στον σχεδιασμό φαρμάκων, και πολιτικοί ακτιβιστές που προσπαθούν να κατανοήσουν τις προτιμήσεις των ψηφοφόρων από δεδομένα δημοσκοπήσεων. Πιθανώς επίσης από ομάδες μάρκετινγκ που προσπαθούν να σου πουλήσουν πράγματα, αλλά κάθε τεχνολογία μπορεί να χρησιμοποιηθεί για καλό ή για κακό.

Η PCA δεν αποτελεί μέρος μιας ελάχιστης εγκατάστασης της Julia, αλλά (έξω από το Exercism) το πακέτο MultivariateStats περιέχει ό,τι χρειάζεσαι.

Επεξεργασία εικόνας

Πηγαίνοντας την ιδέα της PCA ένα βήμα παραπέρα, οι ψηφιακές εικόνες είναι απλώς πίνακες τιμών pixel, και μπορούμε να υπολογίσουμε κύριες συνιστώσες για αυτές.

Αυτό χρησιμοποιείται εδώ και δεκαετίες στη συμπίεση εικόνων, όπου η PCA καθοδηγεί ως προς το τι είναι πιο σημαντικό να διατηρηθεί όταν μειώνεται το μέγεθος του αρχείου.

Όλο και περισσότερο, η PCA αποτελεί ζωτικό μέρος ταξινομητών εικόνων, όπως η αναγνώριση προσώπων. Την επόμενη φορά που θα περπατήσεις σε ένα αεροδρόμιο, ο Μεγάλος Αδερφός δεν παρακολουθεί απλώς, χρησιμοποιεί επίσης τη Γραμμική Άλγεβρα για να καταλάβει τι βλέπει!

Μηχανολογία

Λιγότερο αμφιλεγόμενα, οι κύριοι άξονες ενός μηχανικού εξαρτήματος είναι τα ιδιοδιανύσματα του τανυστή ροπής αδράνειας I (ένας πίνακας με ελαφρώς διαφορετικό όνομα).

Κάθε τροχός του αυτοκινήτου σου έχει πιθανώς ένα ή περισσότερα αντίβαρα, ρυθμισμένα ώστε να μηδενίζουν τα εκτός διαγωνίου στοιχεία αυτού του τανυστή I. Ο μηχανικός του συνεργείου σίγουρα αποφεύγει να κάνει τα μαθηματικά (σε αντίθεση με τους σχεδιαστές αεροσκαφών και τους μηχανικούς πυραύλων), αλλά το "wobble" είναι απλώς μια καθημερινή λέξη για να περιγράψει τους μικτούς όρους του τανυστή, που σε αυτή την περίπτωση θα έκαναν τη διαδρομή σου λιγότερο άνετη και θα αύξαναν τη μηχανική φθορά στα ρουλεμάν.