하이브리드 모델과 대리모델

결합 방식

왜 필요한가. 24.1장에서 본 것처럼 물리 모델의 오차에는 파라미터가 틀린 오차와 식의 구조가 틀린(빠진 항이 있는) 오차가 섞여 있다. 어느 쪽이냐에 따라 데이터로 고쳐야 할 대상이 다르다. 결합 방식은 데이터가 모델의 어디를 고치게 할지를 정하는 선택이다.

정의. 물리 모델을 $f(\mathbf{x}; \boldsymbol{\theta})$, 학습되는 함수를 $g_\phi$라 하자. $\mathbf{x}$는 입력(조건), $\boldsymbol{\theta}$는 물리 파라미터, $\phi$는 학습 파라미터다. 자주 쓰는 결합 방식은 다섯 가지다.

직렬 결합에서 학습된 함수가 미분방정식 안의 한 항이 되면 universal differential equation(24.3장)이 된다.

직관. 잔차 학습은 물리 모델을 그대로 두고 그 위에 보정을 덧칠한다. 파라미터 학습은 물리 모델의 식은 믿고 계수만 데이터에 맞춘다. 보정이 학습 범위 밖에서 어떻게 행동하는지가 둘을 가른다. 잔차 모델이 GP나 신경망이면 학습 범위 밖에서 보정은 대략 상수(GP의 사전 평균)나 임의의 직선(ReLU 신경망, 7.2장)이 된다. 반면 파라미터를 고친 물리 모델은 식의 구조를 따라 외삽한다. 구조가 맞다면 이것이 가장 큰 장점이고, 구조가 틀렸다면 가장 큰 위험이다.

예로 화학 반응의 속도 상수 $k$가 온도 $T$에 따라 변하는 Arrhenius 식을 쓴다.

$A$는 빈도 인자, $E_a$는 활성화 에너지(J/mol), $R$은 기체 상수다. 참 과정은 $\ln A = 18$, $E_a = 60$ kJ/mol이지만, 실험자는 문헌값 $\ln A = 16$, $E_a = 55$ kJ/mol을 가지고 있다. 식의 구조는 맞고 파라미터만 틀린 경우다. 300–340 K에서 15번 측정하고(잡음 표준편차 0.05), 같은 범위의 내삽 시험점과 380–420 K의 외삽 시험점에서 네 모델을 비교한다. GP 회귀는 13.7장의 식을 그대로 쓴다.

import numpy as np

R_J = 8.314                                          # 기체 상수, J/(mol K)
def ln_rate(temp_k, ln_a, ea_j):                     # Arrhenius 식: ln k = ln A - Ea / (R T)
    return ln_a - ea_j / (R_J * temp_k)

def gp(x_tr, y_tr, x_te, ell, sf=1.0, sn=0.05):      # GP 회귀(13.7장): 사후 평균·분산, 평균은 데이터 평균
    k = lambda a, b: sf**2 * np.exp(-0.5 * (a[:, None] - b[None, :])**2 / ell**2)
    L = np.linalg.cholesky(k(x_tr, x_tr) + sn**2 * np.eye(x_tr.size))
    v = np.linalg.solve(L, k(x_tr, x_te))            # (N, M)
    alpha = np.linalg.solve(L.T, np.linalg.solve(L, y_tr - y_tr.mean()))
    return y_tr.mean() + k(x_te, x_tr) @ alpha, sf**2 - np.sum(v**2, axis=0)

truth = lambda t: ln_rate(t, 18.0, 60e3)             # 참 과정(실험자는 모른다)
physics = lambda t: ln_rate(t, 16.0, 55e3)           # 문헌값으로 만든 물리 모델: 식의 구조는 맞고 파라미터가 틀림
rng = np.random.default_rng(0)
temp_tr_k = np.sort(rng.uniform(300, 340, 15))       # (15,) 실험한 온도
ln_k_tr = truth(temp_tr_k) + 0.05 * rng.normal(size=15)
temp_in_k, temp_out_k = np.linspace(305, 335, 50), np.linspace(380, 420, 50)   # 내삽·외삽 시험점
print(f"학습 범위 {temp_tr_k.min():.0f}–{temp_tr_k.max():.0f} K, ln k {ln_k_tr.min():.2f} … {ln_k_tr.max():.2f}")

학습 범위 300–337 K, ln k -6.08 … -3.38

X = np.column_stack([np.ones(15), -1 / (R_J * temp_tr_k)])   # (15, 2) ln k는 (ln A, Ea)에 선형
ln_a_fit, ea_fit = np.linalg.lstsq(X, ln_k_tr, rcond=None)[0]
resid_tr = ln_k_tr - physics(temp_tr_k)                       # 물리 모델이 설명하지 못한 부분
models = {
    "데이터 전용 GP      ": lambda t: gp(temp_tr_k, ln_k_tr, t, ell=15.0)[0],
    "물리 전용(문헌값)   ": physics,
    "잔차 학습 하이브리드": lambda t: physics(t) + gp(temp_tr_k, resid_tr, t, ell=15.0)[0],
    "파라미터 학습       ": lambda t: ln_rate(t, ln_a_fit, ea_fit),
}
rmse = lambda f, t: np.sqrt(np.mean((f(t) - truth(t))**2))
for name, f in models.items():
    print(f"{name}  내삽 RMSE {rmse(f, temp_in_k):.3f}   외삽 RMSE {rmse(f, temp_out_k):.3f}")
print(f"추정 Ea = {ea_fit / 1e3:.1f} kJ/mol (참값 60, 문헌값 55)")

데이터 전용 GP        내삽 RMSE 0.033   외삽 RMSE 4.477
물리 전용(문헌값)     내삽 RMSE 0.130   외삽 RMSE 0.497
잔차 학습 하이브리드  내삽 RMSE 0.031   외삽 RMSE 0.371
파라미터 학습         내삽 RMSE 0.001   외삽 RMSE 0.004
추정 Ea = 60.1 kJ/mol (참값 60, 문헌값 55)

결과는 결합 방식의 성격을 그대로 드러낸다.

• 데이터 전용 GP는 내삽에서 좋지만(0.033) 외삽에서 무너진다(4.477). 학습 범위를 벗어나면 사후 평균이 데이터 평균으로 돌아가기 때문이다.
• 물리 전용 모델은 내삽에서 가장 나쁘지만(0.130) 외삽 오차가 0.497에 그친다. 식이 맞아서 추세는 따라간다.
• 잔차 학습 하이브리드는 내삽에서 GP만큼 좋지만(0.031) 외삽에서는 물리 모델의 오차를 거의 그대로 물려받는다(0.371). 학습된 잔차는 300–340 K의 잔차를 기억할 뿐, 잔차가 온도에 따라 커진다는 사실은 모른다.
• 파라미터 학습은 두 영역 모두에서 가장 좋다. $E_a$를 60.1 kJ/mol로 되찾았다. 이 시드에서 내삽 오차 0.001은 운이 좋은 값이고, 시드 20개 평균은 내삽 0.017, 외삽 0.073이다(실습).

수식으로 본 혼동. 실제 연구에서는 파라미터도 틀리고 구조도 조금 틀린 경우가 대부분이다. Kennedy와 O'Hagan(2001)은 이를

로 썼다. $\boldsymbol{\theta}$는 보정(calibration)할 파라미터, $\delta(\mathbf{x})$는 모델 불일치(model discrepancy) 함수, $\varepsilon_i$는 측정 잡음이다. $\delta$를 GP로 두면 잔차 학습과 파라미터 학습을 함께 하는 셈이다. 문제는 $\boldsymbol{\theta}$를 바꿔서 줄일 수 있는 오차와 $\delta$로 흡수할 수 있는 오차가 겹친다는 점이다. 유연한 $\delta$는 틀린 $\boldsymbol{\theta}$가 만든 오차까지 흡수하고, 그러면 데이터는 $\boldsymbol{\theta}$에 대해 거의 아무것도 말해 주지 않는다. 그래서 $\delta$에는 물리적으로 그럴듯한 크기와 매끄러움을 사전분포로 주고, $\boldsymbol{\theta}$의 사후분포가 사전분포에서 얼마나 움직였는지를 함께 보고한다. 이 식별 문제는 24.6장과 24.7장에서 다시 다룬다.

연구에서의 의미. 결합 방식은 세 질문으로 고른다. 물리 모델의 구조를 얼마나 믿는가(믿으면 파라미터 학습, 의심하면 잔차·직렬 결합). 데이터가 얼마나 많은가(적으면 학습 부분을 작게). 예측이 학습 범위 밖으로 얼마나 나가야 하는가(멀리 나갈수록 물리 구조에 기댄다). 셋을 적어 두지 않으면 결합 방식은 결과가 좋게 나온 것으로 사후에 정해진다.

대리모델

왜 필요한가. 설계 최적화, 불확실성 전파, 베이지안 보정은 시뮬레이터를 수천에서 수백만 번 부른다. 한 번에 몇 시간 걸리는 유한요소·전산유체 시뮬레이션으로는 불가능하다. 대리모델(surrogate model, emulator)은 시뮬레이터의 입력–출력 관계를 적은 실행으로 학습해, 밀리초 단위로 근사값을 내는 모델이다.

정의. 시뮬레이터를 $\mathbf{u} = \mathcal{S}(\mathbf{x})$라 하자. $\mathbf{x}$는 설계·조건 변수(보통 수 개에서 수십 개), $\mathbf{u}$는 출력(스칼라 하나일 수도, 수십만 격자점의 장일 수도 있다)이다. 대리모델 $\hat{\mathcal{S}}$는 실행 결과 $\{(\mathbf{x}_i, \mathbf{u}_i)\}_{i=1}^N$으로 학습한다. 대표적인 세 계열은 다음과 같다.

축소 모델. 출력이 $n$개 격자점의 장이면 $n$이 수십만이라 출력 하나하나를 학습하기 어렵다. 시뮬레이션 $M$번의 결과를 열로 쌓은 snapshot 행렬 $\mathbf{U} = [\mathbf{u}_1, \dots, \mathbf{u}_M] \in \mathbb{R}^{n \times M}$의 SVD $\mathbf{U} = \mathbf{\Phi}\mathbf{\Sigma}\mathbf{W}^\top$에서 앞 $r$개의 왼쪽 특이벡터 $\mathbf{\Phi}_r$을 공간 모드로 쓰면, 각 장은 $\mathbf{u} \approx \mathbf{\Phi}_r \mathbf{a}$로 $r$개의 계수 $\mathbf{a}$만으로 요약된다. 이것이 POD(proper orthogonal decomposition)이고, 계산은 7.6장의 SVD·PCA와 같다(평균을 뺀 뒤 SVD를 하는 것이 보통이다). 모드 수 $r$은 보존한 "에너지" 비율

로 정한다. 그 뒤 입력 $\mathbf{x}$에서 계수 $\mathbf{a}$로 가는 사상만 GP나 신경망으로 학습하면 된다. 지배 방정식을 모드 위로 투영(Galerkin projection)해 작은 ODE를 얻는 방법도 있다.

검증. 대리모델의 검증은 일반 회귀 모델의 검증(16.6장)보다 까다롭다. 세 가지를 지킨다.

• 시험 실행은 학습과 다른 입력에서 새로 돌린 시뮬레이션이어야 한다. 같은 실행의 시간 단계를 나눠 학습·시험에 쓰면 누수다.
• 대리모델이 실제로 쓰일 영역에서 평가한다. 최적화는 설계 공간의 모서리로 가는 경향이 있으므로, 균일하게 뽑은 시험점의 평균 오차만으로는 부족하다. 최적화가 찾은 해는 원래 시뮬레이터로 다시 확인한다.
• 속도 이득은 학습 비용을 포함해 보고한다. 학습용 실행 500번이 필요한 대리모델이 최적화에서 실행 300번을 아꼈다면 손해다.

다중 충실도

왜 필요한가. 같은 현상을 여러 정확도로 계산할 수 있는 경우가 많다. 성긴 격자와 촘촘한 격자, 단순화한 물리와 전체 물리, 시뮬레이션과 실험이 그렇다. 정확한 것은 비싸고 싼 것은 편향되어 있다. 다중 충실도(multi-fidelity) 모델은 싼 저충실도 실행 여러 번과 비싼 고충실도 실행 몇 번을 결합한다.

정의. Kennedy와 O'Hagan(2000)의 자기회귀 구조가 가장 널리 쓰인다.

$f_L$은 저충실도 모델, $f_H$는 고충실도 모델, $\rho$는 척도 계수, $\delta$는 두 충실도의 차이를 나타내는 함수다. $f_L$과 $\delta$에 각각 GP 사전분포를 두면, 저충실도 데이터는 $f_L$의 모양을, 소수의 고충실도 데이터는 $\rho$와 $\delta$를 정한다. $\delta$가 $f_H$보다 단순할수록(더 매끄러울수록) 고충실도 실행을 아낄 수 있다.

코드. 고충실도 함수 $f_H(x) = (6x-2)^2 \sin(12x-4)$와, 이를 축소·이동한 저충실도 함수 $f_L = 0.5 f_H + 10(x - 0.5) - 5$로 실험한다. 고충실도 4점, 저충실도 11점이다. 여기서는 $\rho$와 $\delta$의 선형 추세를 최소제곱으로 정하고, 남은 차이를 GP로 보정하는 간단한 버전을 쓴다.

f_hi = lambda x: (6 * x - 2)**2 * np.sin(12 * x - 4)          # 고충실도: 비싸다
f_lo = lambda x: 0.5 * f_hi(x) + 10 * (x - 0.5) - 5           # 저충실도: 싸지만 편향
x_hi, x_lo = np.array([0.0, 0.4, 0.6, 1.0]), np.linspace(0, 1, 11)
x_te = np.linspace(0, 1, 101)
mu_hi_only, _ = gp(x_hi, f_hi(x_hi), x_te, ell=0.2, sf=5.0, sn=1e-3)
mu_lo = lambda x: gp(x_lo, f_lo(x_lo), x, ell=0.1, sf=5.0, sn=1e-3)[0]
A = np.column_stack([mu_lo(x_hi), np.ones(4), x_hi])          # (4, 3) f_hi ≈ rho·f_lo + (a + b x) + δ(x)
rho, a, b = np.linalg.lstsq(A, f_hi(x_hi), rcond=None)[0]
trend = lambda x: rho * mu_lo(x) + a + b * x
mu_mf = trend(x_te) + gp(x_hi, f_hi(x_hi) - trend(x_hi), x_te, ell=0.2, sf=1.0, sn=1e-3)[0]
for name, mu in [("고충실도 4점만      ", mu_hi_only), ("저충실도 그대로     ", f_lo(x_te)),
                 ("저 11점 + 고 4점 결합", mu_mf)]:
    print(f"{name}  RMSE {np.sqrt(np.mean((mu - f_hi(x_te))**2)):.3f}")
print(f"추정 rho = {rho:.2f} (구성상 참값 2)")

고충실도 4점만        RMSE 5.604
저충실도 그대로       RMSE 6.141
저 11점 + 고 4점 결합  RMSE 0.122
추정 rho = 2.00 (구성상 참값 2)

고충실도 4점만으로는 함수의 굴곡을 따라갈 수 없고(5.604), 저충실도 모델은 모양은 맞지만 크기와 기울기가 틀렸다(6.141). 둘을 결합하면 오차가 0.122로 줄어든다. 다만 이 예에서는 두 충실도의 관계가 선형이 되도록 만들었다는 점에 주의한다. 실제 문제에서는 $\rho$가 입력에 따라 변하거나, 저충실도 모델이 특정 영역(예: 충격파가 생기는 고속 영역)에서 고충실도와 전혀 다르게 행동할 수 있다. 그런 영역에서는 저충실도 데이터가 오히려 잘못된 모양을 강요한다. 결합 전에 두 충실도를 같은 입력에서 몇 번 실행해 상관을 확인한다.

능동 학습

왜 필요한가. 실행 예산이 정해져 있으면 다음 실행을 어디서 할지가 대리모델의 품질을 정한다. 능동 학습(active learning)은 지금까지의 대리모델을 써서 다음 실행할 입력을 고른다. 실험 설계 전반은 25.3장에서 다루고, 여기서는 대리모델을 만드는 관점만 본다.

정의. 후보 입력 $\mathbf{x}$마다 획득 함수(acquisition function) $a(\mathbf{x})$를 계산하고 최댓값을 다음 실행점으로 고른다. 대리모델 전체의 정확도를 높이려면 불확실성이 가장 큰 곳을 고른다.

최적값을 찾는 것이 목적이면 기대 개선량(expected improvement) 같은 획득 함수를 쓰는 베이지안 최적화가 된다(25.3장). 획득 함수가 무엇을 최대화하는지가 실험의 목표와 맞아야 한다.

코드. 초기 실행 4번에서 시작해 10번을 더 고른다. 무작위로 고르는 것과 GP 사후분산이 가장 큰 곳을 고르는 것을, 초기점을 바꿔 가며 10번씩 비교한다.

sim = lambda x: np.sin(6 * x) + 0.5 * np.sin(17 * x)          # 비싼 시뮬레이터라고 가정
x_cand = np.linspace(0, 1, 201)                               # (201,) 후보 입력
def run(pick, n_total=14, n0=4, seed=1):
    r = np.random.default_rng(seed)
    x_tr = r.uniform(0, 1, n0)                                # 초기 실행 4번
    while x_tr.size < n_total:
        x_tr = np.append(x_tr, pick(x_tr, r))
    mu, _ = gp(x_tr, sim(x_tr), x_cand, ell=0.08, sn=1e-3)
    return np.sqrt(np.mean((mu - sim(x_cand))**2))
random_pick = lambda x_tr, r: r.uniform(0, 1)
maxvar_pick = lambda x_tr, r: x_cand[np.argmax(gp(x_tr, sim(x_tr), x_cand, ell=0.08, sn=1e-3)[1])]
for name, pick in [("무작위   ", random_pick), ("최대 분산", maxvar_pick)]:
    errs = [run(pick, seed=s) for s in range(10)]
    print(f"{name}  14번 실행 후 RMSE  평균 {np.mean(errs):.3f}  최악 {np.max(errs):.3f}")

무작위     14번 실행 후 RMSE  평균 0.169  최악 0.428
최대 분산  14번 실행 후 RMSE  평균 0.012  최악 0.018

최대 분산 선택은 평균 오차를 열 배 이상 줄이고, 최악의 경우도 훨씬 안정적이다. 무작위 선택은 점이 한쪽에 몰려 빈 구간을 남기는 일이 잦다.

13.7장에서 본 것처럼 하이퍼파라미터를 고정하면 GP 사후분산은 관측값과 무관하다. 그래서 위의 최대 분산 선택은 사실상 공간을 고르게 채우는 설계이고, 길이척도를 잘못 정하면 그 간격이 틀린다. 실제로는 실행할 때마다 하이퍼파라미터를 다시 추정하므로 관측값이 간접적으로 선택에 영향을 준다. 초기 몇 번은 Latin hypercube 같은 공간 채움 설계로 시작하고, 하이퍼파라미터가 안정된 뒤 능동 학습으로 넘어가는 것이 흔한 절차다.

공정한 비교

하이브리드 모델이 "더 낫다"는 주장은 비교 조건이 같을 때만 의미가 있다. 다음을 고정한 뒤 비교한다.

마지막 두 행이 특히 자주 빠진다. 위의 셀 2도 시드 하나의 결과라 파라미터 학습의 내삽 오차(0.001)가 실제보다 좋게 보였다. 또 물리 전용 baseline을 문헌값 그대로 두는 것은 공정하지 않다. 파라미터만 다시 맞춘 물리 모델이 가장 강한 baseline인 경우가 많고, 이 baseline을 이기지 못하는 하이브리드 모델은 복잡도만 더한 것이다.

참고 문헌

• Kennedy, O'Hagan (2000). Predicting the output from a complex computer code when fast approximations are available. Biometrika, 87(1), 1–13. — 다중 충실도 GP의 자기회귀 구조 $f_H = \rho f_L + \delta$
• Kennedy, O'Hagan (2001). Bayesian calibration of computer models. Journal of the Royal Statistical Society: Series B, 63(3), 425–464. — 파라미터 보정과 모델 불일치 함수를 함께 다루는 베이지안 틀
• Rasmussen, Williams (2006). Gaussian Processes for Machine Learning. MIT Press. — GP 대리모델의 표준 교재