비어 있는 시각 채우기: 보간과 스플라인

어제 6호기의 배터리를 지켜본 다음 날 아침, 민재 팀장이 사고 보고서 양식을 내밉니다. "4호기가 산마을 가는 길에 돌풍을 맞았어요. 다행히 무사히 돌아오고 있는데, 항공청 보고서에 돌풍이 분 순간의 위치를 적어야 해요. 기체 센서로는 12.3초에 돌풍이 왔대요. 그런데 GPS 기록은 0초, 5초, 10초, 15초… 이렇게 5초마다만 있어요."

서윤은 10초와 15초의 기록을 직선으로 이어 12.3초의 값을 읽습니다. 항로에서 옆으로 밀려난 거리가 10초에 0.1 m, 15초에 −0.8 m이니 12.3초에는 약 −0.31 m입니다. "거의 항로 위에 있었네요."

오 박사가 그래프를 들여다봅니다. "점 사이를 직선으로 잇는 건 '그 5초 동안 똑바로 날았다'는 가정이야. 이 기체는 바람에 좌우로 흔들리면서 날았어. 10초와 15초 사이에 더 크게 밀렸다가 돌아왔을 수도 있지. 이번에는 회귀처럼 큰 흐름만 그리면 안 돼. 기록된 점은 정확히 지나면서, 점 사이를 그럴듯하게 채우는 곡선이 필요해."

이번 교시는 그런 곡선을 만드는 보간(interpolation) 을 다룹니다. 점을 모두 지나는 다항식 하나로 시작해, 점이 많아지면 그것이 왜 위험한지(룽게 현상) 보고, 실무의 표준인 3차 스플라인을 직접 짭니다. 4호기가 돌아오면 비행 기록 장치(블랙박스)의 0.1초 간격 기록과 견줘 누가 맞았는지 확인합니다.

10-1. 보간과 회귀는 무엇이 다른가 (5분)

4호기의 기록입니다. 값은 계획 항로에서 옆으로 밀려난 거리(m)이고, 오른쪽이 +입니다. 블랙박스 기록은 아직 없으니, 그것을 흉내 내는 함수 truth를 만들어 두고 관제실은 모른다고 칩니다.

>>> def truth(t):                      # 블랙박스 기록을 흉내 낸 함수 (관제실은 아직 모른다)
...     return 3 * np.sin(0.35 * t) + 0.12 * t
>>> t = np.arange(0, 31, 5)            # GPS 기록 시각 (초)
>>> d = np.round(truth(t), 1)          # 항로에서 옆으로 밀려난 거리 (m)
>>> t
array([ 0,  5, 10, 15, 20, 25, 30])
>>> d
array([ 0. ,  3.6,  0.1, -0.8,  4.4,  4.9,  1. ])

9교시의 회귀와 이번 교시의 보간은 둘 다 "점들로 곡선 만들기"지만 목적이 다릅니다.

GPS 기록이 잡음투성이라면 보간은 잡음까지 정확히 지나는 곡선을 그리므로 회귀가 맞습니다. 여기서는 위치 오차가 수 cm인 정밀 GPS라고 보고 보간을 씁니다. 또 하나, 보간은 점 사이(내삽) 를 채우는 일입니다. 마지막 기록 뒤를 내다보는 외삽은 9교시에서 봤듯이 어떤 곡선이든 위험합니다.

10-2. 점을 모두 지나는 다항식: 라그랑주와 뉴턴 (12분)

점 두 개를 지나는 직선은 하나뿐입니다. 점 세 개(가로 위치가 서로 다른)를 지나는 2차식도 하나뿐입니다. 일반적으로 점 $n$개를 모두 지나는 $n-1$차 이하의 다항식은 딱 하나입니다. 이것을 보간 다항식(interpolating polynomial) 이라고 합니다. 하나뿐이니 어떤 방법으로 구하든 같은 다항식이 나옵니다. 방법마다 다른 것은 계산의 편리함입니다.

라그랑주 형태(Lagrange form). 점마다 "그 점에서는 1, 다른 점에서는 0"인 다항식 $L_i(x)$를 만들고, 높이 $y_i$를 곱해 더합니다.

$\prod$(대문자 파이)는 "모두 곱한다"는 기호입니다. 10초, 15초, 20초의 세 점으로 12.3초의 값을 손으로 구해 봅시다.

세 무게의 합은 정확히 1입니다(0.4158 + 0.7084 − 0.1242). 값은 $0.1 \times 0.4158 - 0.8 \times 0.7084 + 4.4 \times (-0.1242) \approx -1.0716$ m입니다. 직선으로 읽은 −0.31 m보다 세 배 넘게 밀려나 있습니다. 같은 계산을 함수로 적으면 이렇습니다.

>>> def lagrange(xs, ys, x):
...     x = np.asarray(x, dtype=float)
...     total = np.zeros_like(x)
...     for i in range(len(xs)):
...         L = np.ones_like(x)
...         for j in range(len(xs)):
...             if j != i:
...                 L = L * (x - xs[j]) / (xs[i] - xs[j])
...         total = total + ys[i] * L
...     return total
>>> round(float(lagrange(t[2:5], d[2:5], 12.3)), 4)
-1.0716
>>> round(float(lagrange(t, d, 12.3)), 4)
-1.2479

점 7개를 모두 쓴 6차 보간 다항식은 −1.2479 m를 내놓습니다.

뉴턴 형태(Newton form). 라그랑주 형태는 점을 하나 더하면 $L_i$를 모두 새로 계산해야 합니다. 뉴턴 형태는 항을 하나씩 덧붙이는 꼴이라, 점을 더할 때 앞의 계산을 그대로 씁니다.

계수 $c_k$는 분할 차분(divided differences) 표로 구합니다. 첫 열은 $y$ 값, 다음 열은 "이웃한 두 값의 차 ÷ 가로 거리"입니다. $f[x_0, x_1] = \dfrac{y_1 - y_0}{x_1 - x_0}$는 기울기이고, $f[x_0, x_1, x_2] = \dfrac{f[x_1, x_2] - f[x_0, x_1]}{x_2 - x_0}$는 기울기가 바뀌는 빠르기입니다.

굵은 숫자가 계수입니다. 표 아래에 점을 한 줄 더하면 계수가 하나 늘 뿐, 위의 계수는 바뀌지 않습니다. 코드로는 표를 한 열씩 덮어쓰며 만듭니다.

>>> def divided_diff(xs, ys):
...     coef = np.array(ys, dtype=float)
...     for k in range(1, len(xs)):
...         coef[k:] = (coef[k:] - coef[k-1:-1]) / (xs[k:] - xs[:-k])
...     return coef
>>> def newton_eval(xs, coef, x):
...     p = coef[-1]
...     for k in range(len(coef) - 2, -1, -1):
...         p = p * (x - xs[k]) + coef[k]
...     return p
>>> x3 = np.array([10, 15, 20])
>>> c3 = divided_diff(x3, [0.1, -0.8, 4.4])
>>> np.round(c3, 4), round(newton_eval(x3, c3, 12.3), 4)
(array([ 0.1  , -0.18 ,  0.122]), np.float64(-1.0716))
>>> x4 = np.array([10, 15, 20, 5])
>>> c4 = divided_diff(x4, [0.1, -0.8, 4.4, 3.6])
>>> np.round(c4, 4), round(newton_eval(x4, c4, 12.3), 4)
(array([ 0.1   , -0.18  ,  0.122 ,  0.0047]), np.float64(-0.8485))

newton_eval은 1교시부터 써 온 점화식처럼 안쪽 괄호부터 계산합니다(호너 방법). 그런데 점을 하나 더했더니 값이 −1.07에서 −0.85로 크게 움직였습니다. 점이 늘면 다항식이 좋아지기만 할까요? 다음 절의 주제입니다.

np.polyfit(t, d, len(t) - 1)로도 보간 다항식을 얻을 수 있습니다. 차수가 점 개수 − 1이면 최소제곱의 잔차가 0이 되기 때문입니다. 다만 9교시에서 본 반데르몽드 행렬의 조건수 문제가 그대로 따라옵니다.

>>> c6 = np.polyfit(t, d, 6)
>>> f"{c6[-1]:.1e}", round(float(np.polyval(c6, 12.3)), 4)
('1.5e-12', -1.2479)
>>> f"{np.linalg.cond(np.vander(t)):.1e}"
'2.1e+09'

$t = 0$에서 값이 0이므로 상수항은 정확히 0이어야 하는데 $1.5 \times 10^{-12}$가 나왔습니다. 조건수 $2 \times 10^9$만큼 반올림 오차가 커진 흔적입니다. 여기서는 해가 없지만, 점이 더 많아지면 이 오차가 눈에 보이게 커집니다.

10-3. 점을 늘리면 좋아질까: 룽게 현상과 체비쇼프 점 (8분)

1901년 수학자 룽게(Runge)가 놀라운 예를 보였습니다. 종 모양의 매끄러운 함수 $f(x) = \dfrac{1}{1 + 25x^2}$을 $[-1, 1]$에서 같은 간격의 점들로 보간하면, 점을 늘릴수록 가운데는 잘 맞지만 양 끝에서 다항식이 점점 크게 출렁입니다. 이것을 룽게 현상(Runge's phenomenon) 이라고 합니다.

>>> f = lambda x: 1 / (1 + 25 * x**2)
>>> xx = np.linspace(-1, 1, 401)
>>> for n in [5, 9, 13, 15]:
...     xe = np.linspace(-1, 1, n)
...     k = np.arange(n)
...     xc = np.cos((2 * k + 1) * np.pi / (2 * n))
...     err_e = np.abs(lagrange(xe, f(xe), xx) - f(xx)).max()
...     err_c = np.abs(lagrange(xc, f(xc), xx) - f(xx)).max()
...     print(n, round(err_e, 3), round(err_c, 3))
5 0.438 0.402
9 1.045 0.171
13 3.658 0.069
15 7.192 0.047

가운데 열이 등간격 점의 가장 큰 오차입니다. 점 15개에서는 높이가 1도 안 되는 함수를 보간했는데 오차가 7입니다. 이유는 보간 오차의 식에 있습니다.

$\xi$(크시)는 구간 안 어딘가의 점입니다. 뒤쪽의 곱 $\prod (x - x_i)$가 점들의 배치로 정해지는 부분인데, 등간격이면 이 곱이 양 끝 근처에서 가운데보다 수백 배 큽니다($n = 13$에서 약 400배). 끝 쪽 점들은 한쪽에만 이웃이 있어 붙잡아 주는 힘이 약하기 때문입니다.

해결책은 점을 양 끝에 몰아 두는 것입니다. 오른쪽 열의 체비쇼프 점(Chebyshev nodes) $x_k = \cos\dfrac{(2k+1)\pi}{2n}$는 반원 위에 같은 간격으로 찍은 점을 바닥으로 내린 위치입니다. 이렇게 하면 곱 $\prod(x - x_i)$의 가장 큰 값이 $2^{1-n}$으로 가능한 한 작아지고, 점을 늘릴수록 오차가 꾸준히 줄어듭니다.

하지만 GPS처럼 기록 시각을 고를 수 없는 데이터에는 체비쇼프 점을 쓸 수 없습니다. 이때의 교훈은 하나입니다. 점이 많은 긴 구간을 고차 다항식 하나로 덮지 마라. 대신 구간을 잘게 나눠 조각마다 낮은 차수를 씁니다. 다음 두 절이 그 방법입니다. 실습대(10-6)의 '룽게 예제'에서 점이 하나씩 늘며 끝이 출렁이는 모습을 직접 볼 수 있습니다.

10-4. 가장 단순한 답: 구간별 선형 보간 (5분)

이웃한 두 점을 직선으로 잇습니다. 서윤이 처음 한 방법입니다. $t_i \le t \le t_{i+1}$이면

12.3초는 10초와 15초 사이이니 $0.1 + \dfrac{-0.8 - 0.1}{5} \times 2.3 = -0.314$ m입니다. NumPy에서는 np.interp(구할 위치, x 기록, y 기록) 한 줄입니다. 여러 위치를 한꺼번에 넣을 수 있습니다.

>>> np.interp(12.3, t, d)
np.float64(-0.31400000000000006)
>>> np.round(np.interp([2.5, 12.3, 27], t, d), 3)
array([ 1.8  , -0.314,  3.34 ])

구간별 선형 보간은 절대 출렁이지 않습니다. 값이 언제나 이웃한 두 기록 사이에 있습니다. 오차는 구간 폭 $h$의 제곱에 비례해서(대략 $\frac{h^2}{8}\max\lvert f''\rvert$), 점 간격을 반으로 줄이면 오차가 1/4쯤 됩니다. 다항 보간과 달리 점을 늘리면 확실히 좋아집니다.

약점은 이음매에서 꺾인다는 것입니다. 위치를 선형으로 이으면 기울기(속도)가 기록 시각마다 갑자기 바뀝니다. 순간적으로 속도가 바뀌는 드론은 없으니, 속도나 가속도가 필요한 분석(경로 계획, 카메라 움직임)에는 모자랍니다. 또 곡선이 휘는 곳에서는 안쪽을 잘라 먹습니다. 12.3초처럼 두 기록 사이에서 더 크게 밀렸다가 돌아온 경우가 바로 그렇습니다.

10-5. 매끄럽게 잇기: 3차 스플라인 (12분)

두 방법의 좋은 점만 모읍니다. 구간마다 따로 낮은 차수(3차식)를 쓰되(출렁임 방지), 이음매에서 값·기울기·휘어짐(2계 도함수) 이 모두 이어지게 합니다(꺾임 방지). 이것이 3차 스플라인(cubic spline) 입니다. 이름은 옛날 설계사들이 배의 곡면을 그릴 때 쓰던 휘는 나무 자(spline)에서 왔습니다. 핀으로 고정한 점들 사이에서 나무 자가 가장 덜 휘는 모양을 잡는데, 3차 스플라인이 바로 그 모양입니다.

점이 $n$개면 구간이 $n-1$개, 3차식마다 계수가 4개라 모르는 수가 $4(n-1)$개입니다. 조건은 양 끝값 맞추기 $2(n-1)$개, 안쪽 이음매 $n-2$곳에서 기울기와 휘어짐 잇기 $2(n-2)$개로 $4n - 6$개입니다. 두 개가 모자라서 양 끝에 조건을 하나씩 더합니다. 끝에서 휘어짐을 0으로 두는 것이 자연 스플라인(natural spline) 입니다($M_0 = M_{n-1} = 0$). 핀 밖으로 나무 자가 곧게 뻗은 모양입니다.

계산은 이음매의 휘어짐 $M_i = S''(t_i)$를 먼저 구하면 깔끔해집니다. 3차식의 2계 도함수는 직선이므로, 구간 안의 $S''$는 $M_i$와 $M_{i+1}$을 잇는 직선입니다. 이것을 두 번 적분하고 양 끝값을 맞추면 구간 $[t_i, t_{i+1}]$의 식이 나옵니다($h_i = t_{i+1} - t_i$, $a = t_{i+1} - t$, $b = t - t_i$).

남은 조건은 이음매에서 기울기가 같다는 것입니다. 이것을 정리하면 안쪽 점마다 식이 하나씩 나옵니다.

식마다 이웃한 $M$ 세 개만 나오니 행렬은 대각선과 그 양옆만 채워진 삼중대각 행렬(tridiagonal matrix) 입니다. 7교시의 가우스 소거로 바로 풀 수 있고, 가운데 계수 $2(h_{i-1} + h_i)$가 양옆의 합 $h_{i-1} + h_i$보다 커서 8교시의 대각 우세도 만족합니다. 피벗팅이 필요 없고, 가우스-자이델 반복도 수렴합니다.

>>> n, h = len(t), np.diff(t)
>>> A = np.zeros((n, n)); b = np.zeros(n)
>>> A[0, 0] = A[-1, -1] = 1                # 자연 스플라인: 양 끝 M = 0
>>> for i in range(1, n - 1):
...     A[i, i-1:i+2] = h[i-1], 2 * (h[i-1] + h[i]), h[i]
...     b[i] = 6 * ((d[i+1] - d[i]) / h[i] - (d[i] - d[i-1]) / h[i-1])
>>> A[:3]
array([[ 1.,  0.,  0.,  0.,  0.,  0.,  0.],
       [ 5., 20.,  5.,  0.,  0.,  0.,  0.],
       [ 0.,  5., 20.,  5.,  0.,  0.,  0.]])
>>> M = np.linalg.solve(A, b)
>>> np.round(M, 4)
array([ 0.    , -0.4688,  0.1712,  0.408 , -0.3392, -0.1792,  0.    ])

$M$을 구했으면 위의 $S(t)$ 식에 넣기만 하면 됩니다. np.searchsorted로 각 시각이 몇 번째 구간에 있는지 찾습니다. 전체를 함수 하나로 묶습니다(실습대의 NumPy 칸도 이 함수를 씁니다).

>>> def natural_spline(x, y, xq):
...     n, h = len(x), np.diff(x)
...     A = np.zeros((n, n)); b = np.zeros(n)
...     A[0, 0] = A[-1, -1] = 1
...     for i in range(1, n - 1):
...         A[i, i-1:i+2] = h[i-1], 2 * (h[i-1] + h[i]), h[i]
...         b[i] = 6 * ((y[i+1] - y[i]) / h[i] - (y[i] - y[i-1]) / h[i-1])
...     M = np.linalg.solve(A, b)
...     i = np.clip(np.searchsorted(x, xq) - 1, 0, n - 2)
...     a, c = x[i+1] - xq, xq - x[i]
...     return ((M[i] * a**3 + M[i+1] * c**3) / (6 * h[i])
...             + (y[i] / h[i] - M[i] * h[i] / 6) * a + (y[i+1] / h[i] - M[i+1] * h[i] / 6) * c)
>>> round(float(natural_spline(t, d, 12.3)), 4)
-1.2034

정말 이음매에서 꺾이지 않는지, 10초 기록의 왼쪽과 오른쪽 기울기(옆으로 밀리는 속도, m/s)를 아주 작은 걸음으로 재 봅니다.

>>> e = 1e-6
>>> for name, g in [("선형", lambda q: np.interp(q, t, d)), ("스플라인", lambda q: natural_spline(t, d, q))]:
...     print(name, round((g(10) - g(10 - e)) / e, 3), round((g(10 + e) - g(10)) / e, 3))
선형 -0.7 -0.18
스플라인 -0.805 -0.805

구간별 선형은 10초에서 속도가 −0.7 m/s에서 −0.18 m/s로 순간 이동하지만, 스플라인은 −0.805 m/s 그대로 매끄럽게 지나갑니다. 연립방정식은 점 $n$개에 $n \times n$이지만 삼중대각이라, 전용 풀이(토마스 알고리즘, 소거 한 번에 이웃 한 줄만 건드림)를 쓰면 계산량이 $n^3$이 아니라 $n$에 비례합니다. GPS 기록이 수만 개여도 순식간입니다.

10-6. 세 방법 견주기 (3분)

오후에 4호기가 돌아왔습니다. 블랙박스의 0.1초 간격 기록(여기서는 truth)과 견줘 봅니다. 12.3초의 값과, 0~30초 전체에서 가장 크게 틀린 곳을 봅니다.

>>> tt = np.linspace(0, 30, 301)
>>> fits = {"다항": lagrange(t, d, tt), "선형": np.interp(tt, t, d), "스플라인": natural_spline(t, d, tt)}
>>> for name, v in fits.items():
...     err = np.abs(v - truth(tt))
...     print(name, round(float(v[123]), 3), round(float(err.max()), 3), round(tt[err.argmax()], 1))
다항 -1.248 1.223 28.4
선형 -0.314 1.049 22.5
스플라인 -1.203 0.55 28.0
>>> round(float(truth(12.3)), 3)
-1.278

12.3초에 4호기는 항로에서 왼쪽으로 1.28 m 밀려나 있었습니다. 선형 보간(−0.31 m)은 1 m 가까이 틀렸고, 다항 보간과 스플라인은 5 cm 안팎으로 맞혔습니다. 그런데 구간 전체로 보면 다항 보간이 가장 나쁩니다(28.4초, 끝 구간에서 1.22 m). 룽게 현상의 작은 판입니다. 스플라인은 가장 큰 오차가 0.55 m로 가장 작고, 그것도 끝 구간(28초)입니다. 자연 스플라인이 끝의 휘어짐을 억지로 0으로 두었기 때문입니다.

아래 실습대에서 세 방법을 겹쳐 볼 수 있습니다. 남색이 다항 보간, 청록이 구간별 선형, 자주색이 3차 스플라인이고 '강조할 방법'으로 고른 곡선이 굵게 그려집니다. 회색 점선은 블랙박스 기록(참값), 호박색 세로선은 값을 읽는 탐침입니다. 오른쪽 표에 탐침 위치에서 각 방법의 값과 참값과의 차이가 나옵니다. 검은 점(GPS 기록)을 끌거나, '점 더하기'로 기록을 촘촘히 해 보세요. '룽게 예제'를 고르면 회색 곡선이 1/(1+25x²)로 바뀌고, 점 개수 슬라이더와 '점 늘리기 재생'으로 점을 5개에서 15개까지 하나씩 늘리며 다항 보간의 끝이 출렁이는 모습을 볼 수 있습니다. 이때 표의 마지막 열은 구간 전체의 최대 오차입니다.

책의 기호 읽는 법

실무에서는

np.interp(xq, x, y)                                  # 구간별 선형. x는 오름차순이어야 한다
scipy.interpolate.CubicSpline(x, y, bc_type="natural")  # 3차 스플라인. 기본값은 'not-a-knot' 끝 조건
scipy.interpolate.make_interp_spline(x, y, k=3)      # B-스플라인 꼴의 보간 (차수 k를 고른다)
scipy.interpolate.PchipInterpolator(x, y)            # 모양 보존: 단조로운 데이터는 단조롭게, 넘치지 않게
scipy.interpolate.BarycentricInterpolator(x, y)      # 다항 보간을 안정적으로 (체비쇼프 점과 함께)

CubicSpline의 기본 끝 조건(not-a-knot)은 양 끝 두 구간을 같은 3차식으로 만들어, 자연 스플라인보다 끝 오차가 작은 경우가 많습니다. 배터리 %나 누적 거리처럼 넘치면 안 되는 값(100%를 넘거나 거꾸로 줄면 안 됨)에는 스플라인이 살짝 넘칠 수 있으니 PchipInterpolator를 씁니다. 2차원 경로는 동쪽·북쪽 좌표를 시간에 대해 따로 보간합니다. 기록에 잡음이 크면 보간 대신 9교시의 회귀나 매끄럽게 하는 스플라인(make_smoothing_spline)을 씁니다.

자주 하는 실수 (2분)

np.interp에 거꾸로 된 순서 넣기. np.interp는 x가 오름차순이라고 믿고 확인하지 않습니다. 순서가 뒤집히면 오류 없이 엉뚱한 값을 냅니다. 먼저 np.argsort로 정렬하세요.

>>> np.interp(12.3, t[::-1], d[::-1])
np.float64(0.0)

점이 많은 데이터에 고차 다항식 하나. 점 30개를 29차 다항식으로 보간하면 룽게 현상과 조건수 문제를 함께 만납니다. 구간별 방법(스플라인)을 쓰세요.

보간 곡선으로 기록 밖 예측하기. 다항식과 스플라인 모두 마지막 점 밖에서는 근거가 없습니다(9교시의 외삽).

잡음 낀 데이터를 보간하기. 잡음까지 정확히 지나느라 곡선이 들쭉날쭉해집니다. 잡음이 있으면 회귀입니다.

연습 문제 (3분)

1. 27.5초의 위치를 세 방법으로 구하고 참값(truth(27.5))과 견줘 보세요. 끝 구간에서는 어느 방법이 가장 나은가요?
2. 룽게 함수를 등간격 점 11, 21, 41, 81개로 구간별 선형 보간하고 최대 오차를 구해 보세요. 점을 늘릴 때 다항 보간과 무엇이 다른가요?
3. 10-5의 삼중대각 연립방정식 A M = b를 8교시의 가우스-자이델 반복으로 20번 돌려 풀고, np.linalg.solve의 답과 같은지 확인해 보세요.

풀이

>>> [round(float(v), 3) for v in (lagrange(t, d, 27.5), np.interp(27.5, t, d), natural_spline(t, d, 27.5), truth(27.5))]
[1.666, 2.95, 3.23, 2.703]
>>> for m in [11, 21, 41, 81]:
...     xe = np.linspace(-1, 1, m)
...     print(m, round(np.abs(np.interp(xx, xe, f(xe)) - f(xx)).max(), 4))
11 0.0674
21 0.0418
41 0.014
81 0.0037
>>> Mg = np.zeros(len(t))
>>> for sweep in range(20):
...     for i in range(1, len(t) - 1):
...         Mg[i] = (b[i] - A[i, i-1] * Mg[i-1] - A[i, i+1] * Mg[i+1]) / A[i, i]
>>> bool(np.abs(Mg - M).max() < 1e-10)
True

1번은 참값 2.703 m에 구간별 선형(2.95)이 가장 가깝습니다. 스플라인(3.23)은 0.53 m, 다항 보간(1.666)은 1 m 넘게 빗나갑니다. 끝 구간에서는 자연 스플라인의 '끝에서 휘어짐 0'이라는 가정이 실제 경로(끝에서도 휘고 있음)와 맞지 않아 오차가 커집니다. 끝 조건이 중요한 이유이고, CubicSpline의 기본값이 자연 스플라인이 아닌 이유입니다. 어느 방법이든 끝 구간은 가운데보다 덜 믿으세요.

2번은 점을 늘릴수록 오차가 꾸준히 줄어듭니다. 41 → 81개에서 간격이 반이 되자 오차가 약 1/4(0.014 → 0.0037)이 되어, 오차가 $h^2$에 비례한다는 것을 확인할 수 있습니다. 다항 보간은 점을 늘릴수록 오히려 나빠졌습니다. 구간별 방법은 "점을 더 모으면 더 좋아진다"는 당연한 기대를 지켜 줍니다.

3번은 가운데 계수 20이 양옆의 합 10보다 커서(대각 우세) 반복이 빠르게 수렴합니다. 스플라인 연립방정식은 크기가 커도 늘 이렇게 풀기 쉬운 모양입니다.

보고서를 보낸 뒤, 민재 팀장이 하나를 더 묻습니다. "항공청에서 그날 4호기가 실제로 날아간 거리와 쓴 에너지(Wh)도 적으래요. 그런데 기록에는 속도와 전력만 있어요." 속도를 시간에 대해 쌓으면 거리, 전력을 쌓으면 에너지입니다. 다음 교시에는 기록된 값의 넓이를 구하는 수치 적분을 배웁니다. 이번에 배운 보간이 그대로 쓰입니다. 점 사이를 직선으로 이으면 사다리꼴 공식, 곡선으로 이으면 심프슨 공식이 됩니다.