삐뚤빼뚤한 데이터 속에서 가장 그럴듯한 답 찾기
로봇이 보는 세상은 반듯하지 않다. 초음파 센서는 같은 벽을 재도 매번 조금씩 다른 거리를 말하고, 카메라는 같은 물체의 위치를 프레임마다 살짝 다르게 추정한다. 바퀴 엔코더는 바닥이 미끄러지면 실제보다 더 멀리 갔다고 착각한다. 로봇에게 데이터란 깔끔한 정답지가 아니라, 잉크가 번진 지도에 가깝다.
그렇다면 로봇은 어떻게 해야 할까? 측정값 하나하나가 조금씩 틀린다면, 어느 값을 믿어야 할까?
최소제곱법은 이 질문에 대한 가장 오래되고도 강력한 대답 중 하나다. 핵심은 간단하다.
"각 데이터가 조금씩 틀릴 수 있다고 인정하자. 대신 전체 데이터를 가장 덜 틀리게 설명하는 답을 찾자."
앞 장에서 최적화를 배울 때 우리는 비용 함수라는 말을 만났다. 비용 함수는 "이 답이 얼마나 나쁜가"를 숫자로 표현한 것이다. 최소제곱법에서는 보통 오차의 제곱합을 비용으로 사용한다. 즉, 모델이 예측한 값과 실제 측정값 사이의 차이를 구하고, 그 차이를 제곱해서 모두 더한다. 그런 다음 이 값이 가장 작아지는 모델을 찾는다.
수식으로 쓰면 다음과 같다.
$$ J = sum_{i=1}^{N} (y_i - hat{y}_i)^2 $$
여기서 $y_i$는 측정값, $hat{y}_i$는 모델이 예측한 값이다. 둘의 차이인 $y_i - hat{y}_i$를 오차, 또는 잔차 residual라고 부른다. 최소제곱법은 잔차들의 제곱합 $J$를 가장 작게 만드는 모델 파라미터를 찾는다.
복잡하게 들리지만, 로봇 예제로 보면 매우 자연스럽다. 로봇이 벽의 위치를 여러 번 측정했다. 측정값은 흔들린다. 그렇다고 첫 번째 값만 믿을 수도 없고, 마지막 값만 믿을 수도 없다. 모든 측정값을 함께 보았을 때 가장 그럴듯한 벽의 위치를 찾는 것, 이것이 최소제곱법의 감각이다.
왜 하필 제곱할까
오차를 그냥 더하면 문제가 생긴다. 어떤 데이터는 예측보다 위에 있고, 어떤 데이터는 예측보다 아래에 있다. 오차가 양수와 음수로 섞이면 서로 상쇄되어 버린다.
예를 들어 두 오차가 각각 +3, -3이라면 그냥 더했을 때는 0이다. 하지만 모델이 완벽한 것은 아니다. 두 데이터 모두 3만큼 틀렸다. 그런데 합이 0이라고 해서 "오차가 없다"고 말하면 곤란하다.
그래서 오차를 제곱한다.
$$ (+3)^2 = 9, quad (-3)^2 = 9 $$
제곱을 하면 양수와 음수가 모두 양수가 된다. 틀린 정도가 사라지지 않는다. 게다가 큰 오차는 더 크게 벌점을 받는다. 오차가 1이면 제곱도 1이지만, 오차가 10이면 제곱은 100이다. 즉, 최소제곱법은 "대부분의 데이터는 잘 맞추되, 지나치게 크게 벗어난 예측은 더 강하게 피하자"는 성격을 가진다.
제곱합을 쓰는 이유를 정리하면 다음과 같다.
양수 오차와 음수 오차가 서로 상쇄되지 않는다.
큰 오차에 더 큰 벌점을 준다.
수학적으로 다루기 편하고 계산이 안정적인 경우가 많다.
많은 로봇 추정 문제에서 기본 비용 함수로 사용하기 좋다.
물론 제곱합이 항상 완벽한 선택은 아니다. 아주 이상한 값, 즉 이상치 outlier가 섞이면 제곱 때문에 그 값의 영향이 지나치게 커질 수 있다. 하지만 기본적인 센서 노이즈가 비교적 고르게 흔들리는 상황에서는 최소제곱법이 매우 좋은 출발점이 된다.
직선을 맞춘다는 것
가장 쉬운 최소제곱 예제는 점 데이터에 직선을 맞추는 것이다. 로봇이 어떤 센서로 거리를 측정한다고 해 보자. 실제 거리 $x$가 커질수록 센서 출력 $y$도 거의 선형으로 증가한다. 하지만 노이즈 때문에 모든 점이 정확히 한 직선 위에 있지는 않다.
우리는 다음과 같은 모델을 가정한다.
$$ hat{y} = ax + b $$
여기서 $a$는 기울기, $b$는 절편이다. $a$는 입력 $x$가 1 증가할 때 출력이 얼마나 증가하는지를 나타내고, $b$는 $x=0$일 때의 기본값이다.
각 데이터 $(x_i, y_i)$에 대해 예측값은 다음과 같다.
$$ hat{y}_i = ax_i + b $$
오차는 측정값에서 예측값을 뺀 값이다.
$$ e_i = y_i - (ax_i + b) $$
최소제곱법은 다음 비용을 최소화하는 $a$, $b$를 찾는다.
$$ J(a, b) = sum_{i=1}^{N} left(y_i - (ax_i + b)right)^2 $$
말로 바꾸면 이렇다. "모든 점과 직선 사이의 세로 방향 차이를 구하고, 그 차이를 제곱해서 더했을 때 가장 작은 직선을 찾는다."
이제 Python으로 직접 해 보자.
import numpy as np
import matplotlib.pyplot as plt
np.random.seed(0)
# 실제 관계: y = 2.0 * x + 1.0
true_a = 2.0
true_b = 1.0
x = np.linspace(0, 10, 30)
noise = np.random.normal(0, 1.2, size=x.shape)
y = true_a * x + true_b + noise
plt.figure(figsize=(7, 5))
plt.scatter(x, y, label="measured data")
plt.xlabel("x")
plt.ylabel("y")
plt.axis("equal")
plt.grid(True)
plt.legend()
plt.show()그래프를 보면 점들이 대략 오른쪽 위로 올라가는 흐름을 갖지만, 완벽히 일직선은 아니다. 센서 데이터가 늘 그렇듯 조금씩 흔들린다. 이제 이 데이터에 가장 잘 맞는 직선을 찾아 보자.
행렬로 직선 모델 표현하기
직선 모델은 다음과 같다.
$$ hat{y} = ax + b $$
데이터가 여러 개 있을 때는 이것을 행렬 형태로 쓸 수 있다.
$$ begin{bmatrix} hat{y}_1 \ hat{y}_2 \ vdots \ hat{y}_N end{bmatrix}
begin{bmatrix} x_1 & 1 \ x_2 & 1 \ vdots & vdots \ x_N & 1 end{bmatrix} begin{bmatrix} a \ b end{bmatrix} $$
짧게 쓰면 다음과 같다.
$$ hat{y} = Atheta $$
여기서 $A$는 입력 데이터를 담은 행렬이고, $theta$는 우리가 찾고 싶은 파라미터 벡터다.
$$ A = begin{bmatrix} x_1 & 1 \ x_2 & 1 \ vdots & vdots \ x_N & 1 end{bmatrix}, quad theta = begin{bmatrix} a \ b end{bmatrix} $$
우리가 원하는 것은 $Atheta$가 측정값 $y$와 최대한 가까워지는 $theta$를 찾는 것이다.
$$ min_{theta} |Atheta - y|^2 $$
여기서 $|Atheta - y|^2$는 오차 벡터의 제곱합을 뜻한다. NumPy에는 이 문제를 풀어 주는 함수가 준비되어 있다. 바로 np.linalg.lstsq다. lstsq는 least squares, 즉 최소제곱의 줄임말이다.
A = np.column_stack([x, np.ones_like(x)])
theta, residuals, rank, singular_values = np.linalg.lstsq(A, y, rcond=None)
estimated_a, estimated_b = theta
print("estimated a:", estimated_a)
print("estimated b:", estimated_b)
print("residuals:", residuals)theta에는 추정된 기울기와 절편이 들어 있다. 이제 이 직선을 그려 보자.
y_hat = estimated_a * x + estimated_b
plt.figure(figsize=(7, 5))
plt.scatter(x, y, label="measured data")
plt.plot(x, y_hat, color="red", label="least squares line")
plt.xlabel("x")
plt.ylabel("y")
plt.grid(True)
plt.legend()
plt.show()빨간 선은 모든 점을 정확히 지나가지 않는다. 그럴 수밖에 없다. 점들이 노이즈 때문에 흩어져 있기 때문이다. 대신 빨간 선은 전체 점들을 보았을 때 가장 균형 잡힌 위치를 지난다. 어떤 점은 선보다 위에 있고, 어떤 점은 아래에 있다. 최소제곱법은 이 오차들의 제곱합이 가장 작아지도록 선을 고른다.
여기서 중요한 감각은 "정답 직선을 찾았다"가 아니다. 실제 관계가 완전히 선형이 아닐 수도 있고, 노이즈가 특이할 수도 있다. 더 정확히 말하면 "우리가 가정한 직선 모델 안에서, 데이터에 가장 잘 맞는 파라미터를 찾았다"가 맞다.
로봇공학에서는 이 차이가 중요하다. 최소제곱법은 마법처럼 현실의 모든 복잡함을 없애 주지 않는다. 대신 우리가 세운 모델 안에서 가장 그럴듯한 답을 찾아 준다.
직접 비용을 계산해 보기
np.linalg.lstsq를 쓰면 답이 바로 나오지만, 최소제곱법의 감각을 제대로 잡으려면 비용을 직접 계산해 보는 것이 좋다.
def compute_cost(a, b, x, y):
y_pred = a * x + b
errors = y - y_pred
cost = np.sum(errors ** 2)
return cost
cost_estimated = compute_cost(estimated_a, estimated_b, x, y)
cost_wrong_1 = compute_cost(1.0, 1.0, x, y)
cost_wrong_2 = compute_cost(3.0, -2.0, x, y)
print("cost estimated:", cost_estimated)
print("cost wrong 1:", cost_wrong_1)
print("cost wrong 2:", cost_wrong_2)좋은 직선일수록 비용이 작다. 물론 데이터에 노이즈가 있기 때문에 비용이 0이 되지는 않는다. 비용이 0이라는 것은 모든 점이 직선 위에 정확히 놓여 있다는 뜻인데, 실제 센서 데이터에서는 거의 기대하기 어렵다.
이번에는 기울기와 절편을 조금씩 바꾸면서 비용이 어떻게 변하는지 보자.
a_values = np.linspace(0, 4, 80)
b_values = np.linspace(-4, 6, 80)
cost_grid = np.zeros((len(b_values), len(a_values)))
for i, b_val in enumerate(b_values):
for j, a_val in enumerate(a_values):
cost_grid[i, j] = compute_cost(a_val, b_val, x, y)
plt.figure(figsize=(8, 6))
plt.contourf(a_values, b_values, cost_grid, levels=40, cmap="viridis")
plt.colorbar(label="cost")
plt.scatter([estimated_a], [estimated_b], color="red", label="least squares solution")
plt.xlabel("a")
plt.ylabel("b")
plt.grid(True)
plt.legend()
plt.show()이 그림에서 가로축은 기울기 $a$, 세로축은 절편 $b$다. 색은 비용을 나타낸다. 빨간 점은 최소제곱법이 찾은 답이다. 비용 지형을 내려다보면 그 주변이 가장 낮은 골짜기처럼 보인다.
이 장면은 최적화의 핵심 이미지를 보여 준다. 우리는 파라미터 공간 위에 펼쳐진 비용 지형을 보고 있다. 좋은 파라미터는 낮은 곳에 있고, 나쁜 파라미터는 높은 곳에 있다. 최소제곱법은 그 지형에서 가장 낮은 지점을 찾는다.
평균도 최소제곱이다
잠깐 더 쉬운 예를 보자. 로봇이 정지한 상태에서 같은 벽까지의 거리를 20번 측정했다.
np.random.seed(1)
true_distance = 5.0
measurements = true_distance + np.random.normal(0, 0.3, size=20)
print(measurements)
print("mean:", np.mean(measurements))우리는 보통 이런 상황에서 평균을 낸다. 그런데 평균은 단순한 습관이 아니다. 사실 평균은 다음 비용을 최소화하는 값이다.
$$ J(d) = sum_{i=1}^{N} (z_i - d)^2 $$
여기서 $z_i$는 각 측정값이고, $d$는 우리가 추정하려는 진짜 거리다. 즉, 모든 측정값과의 제곱 오차 합을 가장 작게 만드는 하나의 값이 평균이다.
def distance_cost(d, measurements):
return np.sum((measurements - d) ** 2)
d_candidates = np.linspace(3.5, 6.5, 200)
costs = np.array([distance_cost(d, measurements) for d in d_candidates])
best_d = d_candidates[np.argmin(costs)]
mean_d = np.mean(measurements)
plt.figure(figsize=(7, 5))
plt.plot(d_candidates, costs)
plt.axvline(mean_d, color="red", linestyle="--", label="mean")
plt.axvline(best_d, color="green", linestyle=":", label="grid best")
plt.xlabel("estimated distance")
plt.ylabel("sum of squared errors")
plt.grid(True)
plt.legend()
plt.show()
print("mean:", mean_d)
print("grid best:", best_d)평균은 최소제곱법의 가장 단순한 얼굴이다. 하나의 숫자를 추정할 때는 평균이 답이 되고, 직선을 추정할 때는 기울기와 절편이 답이 된다. 나중에 로봇의 위치, 자세, 지도, 카메라 파라미터를 추정할 때도 같은 철학이 반복된다. 측정값과 예측값의 차이를 만들고, 그 차이의 제곱합을 줄인다.
센서 보정: 측정값을 실제값으로 바꾸기
이제 로봇다운 예제로 들어가 보자. 어떤 거리 센서가 있다고 하자. 실제 거리가 1m, 2m, 3m일 때 센서가 정확히 1, 2, 3을 출력하면 좋겠지만, 현실에서는 그렇지 않다. 센서 출력이 실제보다 조금 크게 나오거나, 거리 증가에 대한 반응이 약간 다를 수 있다.
예를 들어 실제 거리와 센서 측정값 사이에 다음과 같은 관계가 있다고 가정해 보자.
$$ text{measured} approx alpha cdot text{true} + beta $$
여기서 $alpha$는 스케일 오차를 나타낸다. 실제 거리가 1m 늘어날 때 센서값이 얼마나 늘어나는지에 해당한다. $beta$는 편향, 즉 전체적으로 더해지는 오차다. 센서가 항상 0.2m 정도 크게 나온다면 $beta$가 양수일 가능성이 크다.
센서 보정의 목표는 이 관계를 추정한 뒤, 새 측정값을 실제값에 가깝게 바꾸는 것이다.
먼저 보정용 데이터를 만들어 보자.
np.random.seed(2)
true_distances = np.linspace(0.5, 8.0, 25)
# 센서의 숨겨진 특성
sensor_alpha = 1.08
sensor_beta = 0.25
sensor_noise = np.random.normal(0, 0.15, size=true_distances.shape)
measured_distances = sensor_alpha * true_distances + sensor_beta + sensor_noise
plt.figure(figsize=(7, 5))
plt.scatter(true_distances, measured_distances, label="calibration data")
plt.xlabel("true distance")
plt.ylabel("measured distance")
plt.grid(True)
plt.legend()
plt.show()이 데이터에서 우리는 실제 거리와 측정값을 둘 다 알고 있다. 이런 데이터를 보정 데이터 calibration data라고 생각하면 된다. 실험실에서 줄자로 실제 거리를 맞춰 놓고 센서값을 기록하는 상황이다.
이제 최소제곱법으로 $alpha$와 $beta$를 추정한다.
A = np.column_stack([true_distances, np.ones_like(true_distances)])
params, residuals, rank, singular_values = np.linalg.lstsq(A, measured_distances, rcond=None)
estimated_alpha, estimated_beta = params
print("estimated alpha:", estimated_alpha)
print("estimated beta:", estimated_beta)추정한 모델은 다음과 같다.
$$ text{measured} approx hat{alpha} cdot text{true} + hat{beta} $$
그래프로 확인해 보자.
predicted_measured = estimated_alpha * true_distances + estimated_beta
plt.figure(figsize=(7, 5))
plt.scatter(true_distances, measured_distances, label="calibration data")
plt.plot(true_distances, predicted_measured, color="red", label="estimated sensor model")
plt.xlabel("true distance")
plt.ylabel("measured distance")
plt.grid(True)
plt.legend()
plt.show()빨간 선은 센서가 실제 거리를 어떻게 측정값으로 바꾸는지 설명하는 모델이다. 그런데 로봇이 실제로 필요한 것은 반대 방향일 때가 많다. 센서가 5.6m라고 말했을 때, 실제 거리가 얼마인지 알고 싶다.
위 관계를 뒤집으면 된다.
$$ text{true} approx frac{text{measured} - hat{beta}}{hat{alpha}} $$
코드로 보정 함수를 만들 수 있다.
def calibrate_distance(measured):
return (measured - estimated_beta) / estimated_alpha
new_measured = np.array([1.5, 3.0, 5.5, 8.0])
calibrated = calibrate_distance(new_measured)
for m, c in zip(new_measured, calibrated):
print(f"measured={m:.2f} -> calibrated true distance={c:.2f}")센서 보정은 로봇에서 매우 흔한 작업이다. 카메라의 픽셀 위치를 실제 거리로 바꾸거나, 모터 명령과 실제 속도 사이의 관계를 추정하거나, 온도에 따라 변하는 센서 출력을 보정할 때도 비슷한 생각을 쓴다.
핵심은 단순하다. 실제값과 측정값 사이에 어떤 모델을 세우고, 보정 데이터를 사용해 그 모델의 파라미터를 최소제곱법으로 찾는다.
모델을 반대로 세우는 방법
방금 예제에서는 다음 모델을 추정했다.
$$ text{measured} approx alpha cdot text{true} + beta $$
그리고 마지막에 식을 뒤집어 실제값을 구했다. 하지만 처음부터 다음 모델을 세울 수도 있다.
$$ text{true} approx c cdot text{measured} + d $$
이 방식은 센서 입력을 넣으면 바로 보정된 값을 출력한다.
A_inv = np.column_stack([measured_distances, np.ones_like(measured_distances)])
params_inv, _, _, _ = np.linalg.lstsq(A_inv, true_distances, rcond=None)
c, d = params_inv
print("c:", c)
print("d:", d)
new_measured = np.array([1.5, 3.0, 5.5, 8.0])
calibrated_direct = c * new_measured + d
for m, c_value in zip(new_measured, calibrated_direct):
print(f"measured={m:.2f} -> calibrated={c_value:.2f}")두 방식은 비슷해 보이지만 완전히 같은 문제는 아니다. 첫 번째 방식은 실제 거리로 측정값을 예측하는 모델을 맞춘다. 두 번째 방식은 측정값으로 실제 거리를 예측하는 모델을 맞춘다. 어떤 쪽이 더 적절한지는 데이터가 어떻게 만들어졌는지, 어떤 값을 입력으로 사용할 것인지에 따라 달라진다.
로봇 소프트웨어를 작성할 때는 "내가 지금 무엇을 예측하려고 하는가?"를 먼저 물어야 한다. 측정값을 설명하려는 것인지, 측정값에서 실제 상태를 복원하려는 것인지에 따라 행렬 $A$의 구성이 달라진다.
로봇 위치 추정과 거리 측정
이제 더 흥미로운 문제로 가 보자. 2D 평면 위에 로봇이 있고, 주변에 위치를 알고 있는 비컨 beacon이 여러 개 있다고 하자. 비컨은 기준점 역할을 하는 장치다. 로봇은 각 비컨까지의 거리를 측정할 수 있다.
비컨 위치가 $(a_i, b_i)$이고 로봇 위치가 $(x, y)$라면, 로봇과 비컨 사이의 거리는 다음과 같다.
$$ r_i = sqrt{(x - a_i)^2 + (y - b_i)^2} $$
양변을 제곱하면 다음과 같다.
$$ (x - a_i)^2 + (y - b_i)^2 = r_i^2 $$
거리 측정값이 정확하다면 원이 만나는 지점에서 로봇 위치를 찾을 수 있다. 각 비컨을 중심으로 반지름 $r_i$인 원을 그리면, 여러 원이 한 점에서 만나는 식이다.
하지만 실제 센서에는 노이즈가 있다. 원들이 한 점에서 딱 만나지 않는다. 어떤 원은 조금 크게, 어떤 원은 조금 작게 그려진다. 이때도 최소제곱법의 사고방식이 등장한다.
"모든 거리 측정값을 최대한 잘 설명하는 로봇 위치는 어디인가?"
먼저 데이터를 만들어 보자.
np.random.seed(3)
# 비컨 위치
beacons = np.array([
[0.0, 0.0],
[8.0, 0.0],
[0.0, 6.0],
[8.0, 6.0],
[4.0, 9.0],
])
true_robot = np.array([3.2, 4.1])
true_ranges = np.linalg.norm(beacons - true_robot, axis=1)
range_noise = np.random.normal(0, 0.15, size=true_ranges.shape)
measured_ranges = true_ranges + range_noise
print("true ranges:", true_ranges)
print("measured ranges:", measured_ranges)거리 기반 위치 추정은 원래 식에 제곱이 들어가기 때문에 선형 문제가 아니다. $(x-a_i)^2$와 $(y-b_i)^2$ 때문에 $x^2$, $y^2$ 항이 생긴다. 그런데 기준 비컨 하나를 골라 다른 식에서 빼면 $x^2$와 $y^2$ 항이 사라져 선형 형태를 만들 수 있다.
비컨 1번과 비컨 $i$번의 식을 각각 쓰면 다음과 같다.
$$ (x-a_1)^2 + (y-b_1)^2 = r_1^2 $$
$$ (x-a_i)^2 + (y-b_i)^2 = r_i^2 $$
두 식을 빼면 $x^2$와 $y^2$가 사라지고, $x$와 $y$에 대한 1차식이 남는다.
정리하면 다음 형태가 된다.
$$ 2(a_i-a_1)x + 2(b_i-b_1)y
r_1^2 - r_i^2 + a_i^2 - a_1^2 + b_i^2 - b_1^2 $$
이것은 다음과 같은 선형 시스템으로 쓸 수 있다.
$$ A begin{bmatrix} x \ y end{bmatrix} = c $$
비컨이 3개 이상이면 식이 여러 개 생긴다. 노이즈가 없다면 모든 식이 한 점에서 만날 수 있지만, 노이즈가 있으면 보통 완벽히 맞는 $(x, y)$가 없다. 그래서 최소제곱으로 푼다.
def estimate_position_from_ranges(beacons, ranges):
ref = beacons[0]
ref_range = ranges[0]
A_rows = []
c_rows = []
a1, b1 = ref
for i in range(1, len(beacons)):
ai, bi = beacons[i]
ri = ranges[i]
A_rows.append([
2 * (ai - a1),
2 * (bi - b1)
])
c_rows.append(
ref_range2 - ri2
+ ai2 - a12
+ bi2 - b12
)
A = np.array(A_rows)
c = np.array(c_rows)
position, residuals, rank, singular_values = np.linalg.lstsq(A, c, rcond=None)
return position, residuals
estimated_robot, residuals = estimate_position_from_ranges(beacons, measured_ranges)
print("true robot position:", true_robot)
print("estimated robot position:", estimated_robot)
print("residuals:", residuals)결과는 실제 위치와 완전히 같지는 않을 것이다. 하지만 여러 개의 거리 측정값을 함께 고려하여 가장 그럴듯한 위치를 찾는다.
그래프로 확인해 보자.
plt.figure(figsize=(8, 7))
plt.scatter(beacons[:, 0], beacons[:, 1], marker="s", s=100, label="beacons")
plt.scatter(true_robot[0], true_robot[1], color="green", s=100, label="true robot")
plt.scatter(estimated_robot[0], estimated_robot[1], color="red", s=100, label="estimated robot")
ax = plt.gca()
for beacon, r in zip(beacons, measured_ranges):
circle = plt.Circle(beacon, r, fill=False, linestyle="--", alpha=0.4)
ax.add_patch(circle)
for i, beacon in enumerate(beacons):
plt.text(beacon[0] + 0.1, beacon[1] + 0.1, f"B{i}")
plt.axis("equal")
plt.xlabel("x")
plt.ylabel("y")
plt.grid(True)
plt.legend()
plt.show()점선 원들이 한 점에서 깔끔하게 만나지 않는 모습을 볼 수 있다. 이것이 현실의 센서 데이터다. 최소제곱법은 그 삐걱거리는 원들 사이에서 가장 납득할 만한 위치를 찾는다.
이 예제는 실제 로봇 위치 추정 문제의 축소판이다. GPS, UWB, Wi-Fi 거리 추정, 초음파 거리 측정 등에서 비슷한 아이디어가 쓰인다. 물론 실제 시스템은 더 복잡하다. 벽에 의한 반사, 신호 지연, 센서 편향, 시간 동기화 문제까지 고려해야 한다. 하지만 출발점은 같다. 측정값과 모델 예측값 사이의 오차를 만들고, 그 오차를 줄이는 위치를 찾는다.
거리 오차를 직접 최소화하면 어떨까
방금 우리는 식을 변형해서 선형 최소제곱 문제로 만들었다. 하지만 원래 거리 모델을 그대로 사용할 수도 있다.
로봇 위치 후보가 $p = (x, y)$일 때, 비컨까지의 예측 거리는 다음과 같다.
$$ hat{r}_i = sqrt{(x-a_i)^2 + (y-b_i)^2} $$
오차는 측정 거리와 예측 거리의 차이다.
$$ e_i = r_i - hat{r}_i $$
비용은 다음과 같다.
$$ J(x, y) = sum_{i=1}^{N} (r_i - hat{r}_i)^2 $$
이 비용을 가장 작게 만드는 $(x, y)$를 찾으면 된다. 이 문제는 $x$, $y$에 대해 비선형이다. 그래서 np.linalg.lstsq처럼 한 번에 딱 풀리지는 않는다. 다음 장에서 배울 경사하강법 같은 반복 최적화가 이런 문제에 자주 사용된다.
그래도 지금은 격자 탐색으로 비용 지형을 살짝 구경해 볼 수 있다.
def range_cost(position, beacons, measured_ranges):
predicted_ranges = np.linalg.norm(beacons - position, axis=1)
errors = measured_ranges - predicted_ranges
return np.sum(errors ** 2)
x_candidates = np.linspace(-1, 9, 120)
y_candidates = np.linspace(-1, 10, 120)
cost_map = np.zeros((len(y_candidates), len(x_candidates)))
for iy, y_val in enumerate(y_candidates):
for ix, x_val in enumerate(x_candidates):
pos = np.array([x_val, y_val])
cost_map[iy, ix] = range_cost(pos, beacons, measured_ranges)
best_index = np.unravel_index(np.argmin(cost_map), cost_map.shape)
best_y = y_candidates[best_index[0]]
best_x = x_candidates[best_index[1]]
grid_best = np.array([best_x, best_y])
print("grid best:", grid_best)
print("linear least squares estimate:", estimated_robot)
print("true robot:", true_robot)시각화해 보자.
plt.figure(figsize=(8, 7))
plt.contourf(x_candidates, y_candidates, cost_map, levels=50, cmap="viridis")
plt.colorbar(label="range cost")
plt.scatter(beacons[:, 0], beacons[:, 1], marker="s", s=90, color="white", edgecolor="black", label="beacons")
plt.scatter(true_robot[0], true_robot[1], color="lime", s=100, label="true robot")
plt.scatter(estimated_robot[0], estimated_robot[1], color="red", s=100, label="linear LS")
plt.scatter(grid_best[0], grid_best[1], color="cyan", s=100, label="grid best")
plt.axis("equal")
plt.xlabel("x")
plt.ylabel("y")
plt.grid(True)
plt.legend()
plt.show()색이 어두운 곳과 밝은 곳은 비용의 크기를 나타낸다. 비용이 가장 작은 골짜기 근처가 로봇 위치로 추정된다. 이 그림은 앞으로 배울 경사하강법의 무대이기도 하다. 로봇은 비용 지형 위에서 낮은 곳을 찾아 이동하는 셈이다. 단, 실제 로봇이 공간을 움직이는 것이 아니라 파라미터 후보가 수학적인 지형 위를 움직인다.
np.linalg.lstsq가 돌려주는 값들
지금까지 np.linalg.lstsq를 여러 번 사용했다. 이 함수는 보통 네 가지 값을 반환한다.
theta, residuals, rank, singular_values = np.linalg.lstsq(A, y, rcond=None)각 값의 의미는 다음과 같다.
theta는 최소제곱으로 찾은 파라미터다.residuals는 남은 오차의 제곱합에 해당하는 값이다.rank는 행렬A가 가진 독립적인 정보의 수를 나타낸다.singular_values는 행렬의 구조를 분석할 때 쓰이는 값으로, 문제가 안정적인지 판단하는 데 도움을 준다.
처음에는 theta만 잘 봐도 충분하다. 그러나 로봇 시스템이 커지면 나머지 값들도 중요해진다. 예를 들어 비컨들이 거의 일직선으로 놓여 있으면 로봇 위치를 안정적으로 추정하기 어려울 수 있다. 이때 행렬의 랭크나 특잇값을 보면 문제가 불안정하다는 신호를 얻을 수 있다.
좋은 센서 배치는 좋은 수학 문제를 만든다. 비컨이 모두 한쪽에 몰려 있으면 위치 추정이 흔들리고, 카메라가 특징점을 한 방향에서만 보면 자세 추정이 불안정해진다. 최소제곱법은 데이터를 잘 활용하지만, 데이터 자체가 나쁘면 답도 흔들릴 수밖에 없다.
가중 최소제곱: 더 믿을 센서에 더 큰 목소리 주기
모든 측정값을 똑같이 믿어도 될까? 그렇지 않은 경우가 많다. 가까운 비컨의 거리 측정은 정확하지만, 먼 비컨의 측정은 더 불안정할 수 있다. 카메라가 선명하게 본 특징점은 믿을 만하지만, 흐릿한 특징점은 덜 믿어야 한다.
이럴 때는 가중 최소제곱 weighted least squares을 생각할 수 있다. 기본 최소제곱은 모든 오차를 똑같이 더한다.
$$ J = sum_i e_i^2 $$
가중 최소제곱은 각 오차에 가중치 $w_i$를 곱한다.
$$ J = sum_i w_i e_i^2 $$
믿을 만한 측정값에는 큰 가중치를 주고, 불확실한 측정값에는 작은 가중치를 준다. 확률 장에서 배운 분산을 떠올려 보자. 분산이 작다는 것은 측정값이 덜 흔들린다는 뜻이다. 따라서 보통은 분산이 작은 센서에 더 큰 가중치를 준다.
간단한 직선 예제로 보자. 어떤 점들은 노이즈가 작고, 어떤 점들은 노이즈가 크다고 해 보자.
np.random.seed(4)
x = np.linspace(0, 10, 40)
true_y = 1.5 * x + 0.5
noise_std = np.where(x < 5, 0.3, 1.5)
y = true_y + np.random.normal(0, noise_std)
weights = 1.0 / (noise_std ** 2)
A = np.column_stack([x, np.ones_like(x)])
# 일반 최소제곱
params_plain, _, _, _ = np.linalg.lstsq(A, y, rcond=None)
# 가중 최소제곱
sqrt_w = np.sqrt(weights)
A_weighted = A * sqrt_w[:, None]
y_weighted = y * sqrt_w
params_weighted, _, _, _ = np.linalg.lstsq(A_weighted, y_weighted, rcond=None)
print("plain:", params_plain)
print("weighted:", params_weighted)가중 최소제곱에서는 행렬 A와 측정값 y에 sqrt_w를 곱했다. 왜 제곱근일까? 비용이 $sum_i w_i e_i^2$이기 때문이다. 오차 $e_i$에 $sqrt{w_i}$를 곱한 뒤 제곱하면 $w_i e_i^2$가 된다.
그래프로 비교해 보자.
y_plain = params_plain[0] * x + params_plain[1]
y_weighted_fit = params_weighted[0] * x + params_weighted[1]
plt.figure(figsize=(8, 5))
plt.scatter(x, y, c=noise_std, cmap="coolwarm", label="data")
plt.plot(x, true_y, "k--", label="true line")
plt.plot(x, y_plain, color="red", label="plain LS")
plt.plot(x, y_weighted_fit, color="green", label="weighted LS")
plt.xlabel("x")
plt.ylabel("y")
plt.grid(True)
plt.legend()
plt.show()색이 다른 점들은 노이즈 크기가 다르다는 뜻이다. 가중 최소제곱은 더 정확하다고 알려진 데이터에 더 많은 힘을 실어 준다. 이것은 센서 퓨전과도 직접 연결된다. 엔코더는 짧은 시간 동안은 정확하지만 오래 지나면 드리프트가 커지고, GPS는 절대 위치를 주지만 순간적으로 튈 수 있다. 각 센서의 불확실성을 알고 있다면, 최소제곱 문제 안에서 그 신뢰도를 가중치로 표현할 수 있다.
이상치가 있을 때의 조심스러움
최소제곱법은 큰 오차를 강하게 벌준다. 이것은 장점이지만 때로는 약점이 된다. 만약 데이터에 말도 안 되게 튄 값이 하나 들어오면 어떻게 될까?
np.random.seed(5)
x = np.linspace(0, 10, 30)
y = 2.0 * x + 1.0 + np.random.normal(0, 0.7, size=x.shape)
# 이상치 추가
y[25] += 15.0
A = np.column_stack([x, np.ones_like(x)])
params_outlier, _, _, _ = np.linalg.lstsq(A, y, rcond=None)
y_fit_outlier = params_outlier[0] * x + params_outlier[1]
plt.figure(figsize=(7, 5))
plt.scatter(x, y, label="data with outlier")
plt.plot(x, y_fit_outlier, color="red", label="least squares")
plt.xlabel("x")
plt.ylabel("y")
plt.grid(True)
plt.legend()
plt.show()튀는 점 하나가 직선을 끌어당기는 모습을 볼 수 있다. 오차를 제곱하기 때문에 큰 오차가 비용을 크게 흔든다.
로봇 센서에서도 이런 일이 자주 생긴다. LiDAR가 유리벽에서 이상한 반사를 받거나, 카메라가 그림자를 물체로 착각하거나, GPS가 건물 사이에서 갑자기 튈 수 있다. 이런 상황에서는 이상치를 제거하거나, 제곱 오차 대신 이상치에 덜 민감한 비용 함수를 사용하기도 한다.
하지만 여기서 너무 멀리 갈 필요는 없다. 지금 단계에서 기억할 것은 하나다. 최소제곱법은 기본 도구로 매우 강력하지만, 데이터에 이상치가 많으면 결과가 흔들릴 수 있다. 따라서 로봇 시스템에서는 "좋은 최적화 방법"만큼이나 "좋은 데이터 선별"이 중요하다.
로봇 문제에서 최소제곱법이 나타나는 곳
최소제곱법은 로봇공학의 여기저기에 숨어 있다. 이름이 직접 보이지 않아도, 안쪽에서는 오차의 제곱합을 줄이는 문제가 돌아가고 있는 경우가 많다.
대표적인 예는 다음과 같다.
센서 보정에서 실제값과 측정값 사이의 선형 관계를 추정한다.
거리 측정값을 이용해 로봇의 위치를 추정한다.
카메라 이미지의 특징점과 3D 점 사이의 재투영 오차를 줄인다.
SLAM에서 로봇 궤적과 지도 특징점의 오차를 동시에 줄인다.
로봇 팔의 관절각과 손끝 위치 사이의 모델 파라미터를 보정한다.
모든 예제가 같은 문장으로 요약된다.
"모델이 예측한 값과 센서가 측정한 값의 차이를 만들고, 그 차이가 전체적으로 가장 작아지는 상태나 파라미터를 찾는다."
이 문장은 로봇 수학에서 매우 강력한 주문처럼 작동한다. 위치 추정, 지도 작성, 센서 보정, 경로 최적화가 모두 이 문장 주변에 모인다.
작은 프로젝트: 거리 센서 보정 후 위치 추정하기
마지막으로 두 가지 아이디어를 연결해 보자. 먼저 거리 센서가 약간의 스케일 오차와 편향을 가진다고 하자. 우리는 보정 데이터로 센서 모델을 추정한다. 그런 다음 보정된 거리값을 사용해 로봇 위치를 추정한다.
import numpy as np
import matplotlib.pyplot as plt
np.random.seed(10)
# 1. 거리 센서 보정 데이터 만들기
calib_true = np.linspace(1.0, 10.0, 30)
alpha_real = 1.05
beta_real = 0.30
calib_measured = alpha_real * calib_true + beta_real + np.random.normal(0, 0.12, size=calib_true.shape)
# measured -> true 모델 직접 추정
A_calib = np.column_stack([calib_measured, np.ones_like(calib_measured)])
calib_params, _, _, _ = np.linalg.lstsq(A_calib, calib_true, rcond=None)
c, d = calib_params
print("calibration model: true ≈ c * measured + d")
print("c:", c)
print("d:", d)
def calibrate(measured):
return c * measured + d
# 2. 비컨과 로봇 위치 설정
beacons = np.array([
[0.0, 0.0],
[10.0, 0.0],
[0.0, 8.0],
[10.0, 8.0],
])
true_robot = np.array([6.0, 3.5])
true_ranges = np.linalg.norm(beacons - true_robot, axis=1)
# 센서가 편향된 거리 측정값을 출력한다고 가정
raw_ranges = alpha_real * true_ranges + beta_real + np.random.normal(0, 0.12, size=true_ranges.shape)
calibrated_ranges = calibrate(raw_ranges)
# 3. 보정 전후 위치 추정 비교
def estimate_position_from_ranges(beacons, ranges):
ref = beacons[0]
ref_range = ranges[0]
a1, b1 = ref
A_rows = []
c_rows = []
for i in range(1, len(beacons)):
ai, bi = beacons[i]
ri = ranges[i]
A_rows.append([2 * (ai - a1), 2 * (bi - b1)])
c_rows.append(
ref_range2 - ri2
+ ai2 - a12
+ bi2 - b12
)
A = np.array(A_rows)
c_vec = np.array(c_rows)
pos, _, _, _ = np.linalg.lstsq(A, c_vec, rcond=None)
return pos
estimated_raw = estimate_position_from_ranges(beacons, raw_ranges)
estimated_calibrated = estimate_position_from_ranges(beacons, calibrated_ranges)
print("true robot:", true_robot)
print("estimated with raw ranges:", estimated_raw)
print("estimated with calibrated ranges:", estimated_calibrated)
print("raw error:", np.linalg.norm(estimated_raw - true_robot))
print("calibrated error:", np.linalg.norm(estimated_calibrated - true_robot))이제 그림으로 확인하자.
plt.figure(figsize=(8, 7))
plt.scatter(beacons[:, 0], beacons[:, 1], marker="s", s=100, label="beacons")
plt.scatter(true_robot[0], true_robot[1], color="green", s=120, label="true robot")
plt.scatter(estimated_raw[0], estimated_raw[1], color="red", s=100, label="raw estimate")
plt.scatter(estimated_calibrated[0], estimated_calibrated[1], color="blue", s=100, label="calibrated estimate")
for i, beacon in enumerate(beacons):
plt.text(beacon[0] + 0.15, beacon[1] + 0.15, f"B{i}")
plt.axis("equal")
plt.xlabel("x")
plt.ylabel("y")
plt.grid(True)
plt.legend()
plt.show()보정 전 거리값으로 추정한 위치와 보정 후 거리값으로 추정한 위치가 다를 것이다. 노이즈가 있기 때문에 보정 후가 항상 극적으로 좋아진다고 단정할 수는 없지만, 센서에 일정한 스케일 오차와 편향이 있다면 보정은 위치 추정의 품질을 크게 높일 수 있다.
이 작은 프로젝트에는 로봇 수학의 중요한 흐름이 모두 들어 있다. 먼저 센서의 오차 구조를 모델로 표현한다. 그다음 최소제곱법으로 모델 파라미터를 추정한다. 마지막으로 보정된 데이터를 사용해 로봇의 상태를 다시 추정한다. 데이터가 모델을 만들고, 모델이 다시 데이터를 해석하는 순환이다.
최소제곱법을 사용할 때의 체크리스트
최소제곱법을 쓸 때는 코드를 작성하기 전에 문제를 문장으로 정리하는 습관이 좋다. 수식보다 먼저 다음 질문에 답해 보자.
내가 예측하려는 값은 무엇인가?
측정값과 예측값의 오차를 어떻게 정의할 것인가?
찾고 싶은 파라미터는 무엇인가?
모델은 선형인가, 비선형인가?
모든 측정값을 똑같이 믿을 것인가, 가중치를 줄 것인가?
이상치가 들어왔을 때 결과가 크게 흔들리지 않는가?
이 질문에 답할 수 있으면 행렬 $A$와 벡터 $y$를 만드는 일이 훨씬 쉬워진다. 최소제곱법에서 가장 어려운 부분은 np.linalg.lstsq를 호출하는 코드가 아니다. "무엇을 오차로 볼 것인가"를 정하는 일이다.
로봇에게 오차는 단순한 실패가 아니다. 오차는 방향을 알려 주는 신호다. 센서가 얼마나 빗나갔는지, 모델이 무엇을 놓치고 있는지, 현재 추정이 어디로 움직여야 하는지를 말해 준다. 최소제곱법은 그 오차들을 모아 하나의 비용으로 만들고, 더 나은 답을 향해 길을 열어 준다.
정리
최소제곱법은 측정값과 모델 예측값의 차이를 오차로 정의하고, 그 오차의 제곱합을 최소화하는 방법이다. 오차를 제곱하면 양수와 음수가 상쇄되지 않고, 큰 오차에 더 큰 벌점을 줄 수 있다.
직선 맞추기에서는 $y = ax + b$ 모델을 세우고, np.linalg.lstsq를 사용해 기울기와 절편을 추정했다. 센서 보정에서는 실제값과 측정값 사이의 선형 관계를 추정한 뒤, 새 측정값을 실제값에 가깝게 바꾸는 방법을 보았다. 거리 기반 위치 추정에서는 여러 비컨까지의 거리 측정값을 사용해 로봇의 2D 위치를 추정했다.
이 장에서 본 모든 예제의 중심에는 같은 구조가 있었다.
$$ text{오차} = text{측정값} - text{모델 예측값} $$
$$ text{비용} = text{오차 제곱의 합} $$
최소제곱법은 로봇이 불완전한 센서 데이터 속에서 가장 그럴듯한 답을 찾게 해 주는 기본 도구다. 다음 단계에서는 이 비용을 한 번에 풀기 어려울 때, 조금씩 답을 개선해 나가는 방법을 배운다. 그 주인공이 바로 경사하강법이다.
이 노트는 보도·비평·교육·연구 목적으로 작성되었습니다. 저작권 문의가 있으시면 에서 알려주세요.