선형 방정식 풀이 에서

선형 방정식 풀이 소개

고등학교 때 대수학을 처음 배웠던 게 기억나요?

보통 연립방정식이 주어지고, 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과 영공간

중요한 특수한 경우는 모든 방정식의 우변이 0일 때예요.

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

위 예제에서는 유일한 해가 x = y = 0일 때뿐인데, 보통은 흥미롭지 않아요.

A x = 0의 비자명해는 A의 행렬식이 0일 때만 존재해요.

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) 값들은 2차원 공간에서 하나의 직선 위에 있어요. 해가 무한히 많고, 이 해들이 행렬의 영공간을 이뤄요.

이 경우 행렬의 행들은 선형 독립이 아니고, 행렬의 계수는 행 또는 열의 개수보다 작아요.

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로 나눠서 x를 우변에 남기면 되지 않을까 생각해요.

사실 A에는 역행렬이 있고(모든 정사각 행렬은 아니지만 대부분 있어요. 행렬식이 0이 아니어야 하죠), 계산은 (느리게!) 동작해요.

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

안타깝게도 역행렬 계산은 작은 행렬에서는 느리고, 큰 행렬에서는 빙하처럼 느려요. 계산 중에 다리가 무너진다면 이상적이지 않겠죠!

다행히 다른 알고리즘은 훨씬 빨라요(이 경우에는 가우스 소거법).

Julia는 (Matlab을 따라) 풀이에 백슬래시를 사용해요(엄밀히는 '왼쪽 나눗셈').

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

이건 아주 간단한 예제지만, 주의해야 할 몇 가지 세부 사항이 있어요.

• 계수들은 열에 맞춰 정렬되어야 하므로, 빠진 항에는 필요에 따라 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

백슬래시 풀이법은 해를 하나 내놓아요! 비공식적으로 듣기로는, 이것이 노름이 가장 작은 해(기하학적으로 원점에 가장 가까운 해)라고 해요.

다른 해를 구하려면 A2의 nullspace가 필요해요.

원래 해에 영공간의 임의의 스칼라 배수를 더하면 또 다른 유효한 해가 돼요.

# 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는 0이 아닌 벡터 x의 길이를 λ배로 조정하면서, 방향은 바꾸지 않아요(음수 λ가 방향을 뒤집는 경우를 제외하고).

별로 대중적이지 않은 것처럼 들리지만, 알고 보면 말도 안 되게 유용해요!

용어: λ의 유효한 값은 A의 _고윳값_이고, 대응하는 x 값은 A의 _고유벡터_예요. 안타깝게도 중간에 언어가 바뀌는 단어들은 그냥 받아들여야 해요.

학생들은 보통 2×2 행렬의 고윳값/고유벡터를 손으로 계산하는 법을 배우지만, 컴퓨터를 쓰는 게 훨씬 쉬워요(그래도 꽤 느리고, 일반적인 경우 n×n 행렬에서 O(n^3)으로 증가해요).

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 아이디어를 한 단계 더 밀고 나가면, 디지털 이미지는 그저 픽셀 값의 행렬이고, 이에 대한 주성분을 계산할 수 있어요.

이것은 수십 년 동안 이미지 압축에 사용되어 왔어요. 파일 크기를 줄일 때 무엇을 가장 중요하게 보존해야 하는지 PCA가 길잡이가 되어 줘요.

점점 더 PCA는 얼굴 인식 같은 이미지 분류기의 핵심 부분이 되고 있어요. 다음에 공항을 걸어 지나갈 때, 빅 브라더는 그저 지켜보는 것만이 아니라 선형대수를 이용해 보이는 것을 이해해요!

기계 공학

덜 논쟁적으로 말하면, 기계 부품의 주축은 관성 모멘트 텐서 I (조금 다른 이름의 행렬)의 고유벡터예요.

자동차의 각 바퀴에는 아마 이 텐서 I의 비대각 원소를 0으로 만들기 위해 조정된 균형추가 하나 이상 있을 거예요. 자동차 정비사는 분명히 수학을 하는 걸 피하겠지만(항공기 설계자나 로켓 엔지니어와는 달리), '흔들림'은 이 텐서의 교차항을 설명하는 일상적인 말일 뿐이에요. 이 경우 교차항은 여정을 덜 편안하게 만들고 베어링의 기계적 마모를 늘려요.