제약 최적화와 쌍대성

Lagrange 승수

등식 제약만 있는 문제부터 본다.

기하적 해석. 제약이 하나뿐인 2차원 문제를 생각한다. 제약 $h(\mathbf{x}) = 0$은 평면 위의 곡선이고, 우리는 이 곡선 위를 따라가며 $f$가 가장 작은 점을 찾는다. 곡선 위의 어떤 점에서 $f$의 등고선이 곡선을 가로지르면, 곡선을 따라 한쪽으로 움직여 더 낮은 등고선으로 갈 수 있으므로 그 점은 최적이 아니다. 따라서 최적해에서는 등고선이 제약 곡선에 접한다. 접한다는 것은 두 곡선의 법선이 평행하다는 뜻이고, $f$의 등고선의 법선은 $\nabla f$, 제약 곡선의 법선은 $\nabla h$다.

이 비례 상수 $\nu$를 Lagrange 승수(Lagrange multiplier)라 한다. 제약이 여러 개이면 $\nabla f$가 제약 법선들의 선형결합이 된다. 이 조건을 한 함수로 묶은 것이 Lagrangian이다.

$\nabla_{\mathbf{x}} L = \mathbf{0}$은 위의 평행 조건이고, $\nabla_{\boldsymbol{\nu}} L = \mathbf{0}$은 제약 $h_j(\mathbf{x}) = 0$ 자체다. 제약이 있는 문제가 $L$의 stationary point를 찾는 문제로 바뀐다. 이 조건은 최적해에서 제약 gradient들이 선형독립일 때(regularity 조건) 필요조건이다.

예로 $\min x^2 + 2y^2$ s.t. $x + y = b$를 푼다. $\nabla_{\mathbf{x}} L = (2x + \nu,\ 4y + \nu) = \mathbf{0}$과 $x + y = b$는 미지수 세 개의 선형 연립방정식이다.

import numpy as np

# min x^2 + 2y^2  s.t.  x + y = b.   L = x^2 + 2y^2 + nu (x + y - b)
def solve(b):
    K = np.array([[2.0, 0.0, 1.0],      # dL/dx = 2x + nu = 0
                  [0.0, 4.0, 1.0],      # dL/dy = 4y + nu = 0
                  [1.0, 1.0, 0.0]])     # x + y = b
    x, y, nu = np.linalg.solve(K, [0.0, 0.0, b])
    return x, y, nu, x**2 + 2 * y**2

x, y, nu, p_star = solve(1.0)
print(f"x = {x:.4f}, y = {y:.4f}, nu = {nu:.4f}, f* = {p_star:.4f}")
print("grad f =", np.array([2 * x, 4 * y]).round(4), "  grad h = [1. 1.]")
d = 1e-3
print(f"df*/db (수치) = {(solve(1 + d)[3] - solve(1 - d)[3]) / (2 * d):.4f}   -nu = {-nu:.4f}")

x = 0.6667, y = 0.3333, nu = -1.3333, f* = 0.6667
grad f = [1.3333 1.3333]   grad h = [1. 1.]
df*/db (수치) = 1.3333   -nu = 1.3333

$b = 1$일 때 해는 $(2/3, 1/3)$이고, 그 점에서 $\nabla f = (4/3, 4/3)$은 0이 아니지만 $\nabla h = (1, 1)$과 평행하다.

승수의 의미. 마지막 줄이 승수의 해석이다. 제약의 우변 $b$를 조금 바꾸면 최적값 $f^\star(b)$가 바뀌는데, 그 변화율이 $-\nu$와 같다.

여기서는 $f^\star(b) = \frac{2}{3}b^2$이므로 $b = 1$에서 $4/3$이고, 승수 $\nu = -4/3$의 부호를 바꾼 값과 같다. 승수는 "제약을 한 단위 풀어 주면 목적함수를 얼마나 개선할 수 있는가", 즉 제약의 값이다. 경제학에서는 이를 shadow price라 부른다. 등식 제약의 승수는 부호에 제한이 없다. $h = 0$을 $-h = 0$으로 써도 같은 제약이므로 승수의 부호는 표기에 달려 있다.

KKT 조건

부등식 제약 $g_i(\mathbf{x}) \le 0$이 있으면 상황이 하나 더 늘어난다. 최적해에서 제약은 활성(active, $g_i = 0$, 경계 위)이거나 비활성(inactive, $g_i < 0$, 경계 안쪽)이다. 비활성 제약은 최적해 근처에서 아무 역할도 하지 않으므로 없는 것처럼 다뤄야 하고, 활성 제약은 등식 제약처럼 다뤄야 한다. 그런데 어느 제약이 활성인지는 풀기 전에는 모른다. 이 두 경우를 한 번에 표현한 것이 KKT 조건(Karush–Kuhn–Tucker conditions)이다. Lagrangian을

로 두면, 최적해 $\mathbf{x}^\star$와 승수 $\boldsymbol{\lambda}^\star, \boldsymbol{\nu}^\star$는 다음 네 조건을 만족한다.

1. stationarity: $\nabla f(\mathbf{x}^\star) + \sum_i \lambda_i^\star \nabla g_i(\mathbf{x}^\star) + \sum_j \nu_j^\star \nabla h_j(\mathbf{x}^\star) = \mathbf{0}$
2. primal feasibility: $g_i(\mathbf{x}^\star) \le 0$, $h_j(\mathbf{x}^\star) = 0$
3. dual feasibility: $\lambda_i^\star \ge 0$
4. 상보적 여유(complementary slackness): $\lambda_i^\star\, g_i(\mathbf{x}^\star) = 0$

상보적 여유가 활성·비활성의 구분을 담는다. 각 부등식 제약에 대해 승수가 0이거나($\lambda_i = 0$, 제약이 해를 밀어내지 않음) 제약이 등호로 성립한다($g_i = 0$, 경계 위). dual feasibility $\lambda_i \ge 0$은 방향을 정한다. 부등식 제약은 feasible set 안쪽으로만 밀 수 있으므로, $-\nabla f$가 feasible set 바깥을 향해야(즉 내려가는 방향이 막혀 있어야) 최적이다.

KKT 조건은 적절한 제약 자격 조건(constraint qualification) 아래에서 최적해의 필요조건이고, 문제가 볼록이면(목적함수와 $g_i$가 볼록, $h_j$가 아핀) 충분조건이기도 하다. 볼록 문제에서는 KKT 조건을 만족하는 점을 찾으면 그것이 전역 최적해다.

10.1장의 첫 문제를 일반화해, 점 $\mathbf{p} = (2, 1)$에서 반공간 $\{\mathbf{x} : \mathbf{a}^\top\mathbf{x} \le c\}$까지 가장 가까운 점을 찾는다. stationarity $2(\mathbf{x} - \mathbf{p}) + \mu\mathbf{a} = \mathbf{0}$에서 $\mathbf{x} = \mathbf{p} - \frac{\mu}{2}\mathbf{a}$이고, 제약이 활성이면 $\mathbf{a}^\top\mathbf{x} = c$에서 $\mu$가 정해진다. $c$를 바꿔 가며 승수와 여유를 본다.

# min ||x - p||^2  s.t.  a^T x <= c   (반공간으로의 투영)
p, a = np.array([2.0, 1.0]), np.array([1.0, 1.0])
for c in [0.0, 1.0, 2.0, 3.0, 4.0]:
    mu = max(0.0, 2 * (a @ p - c) / (a @ a))        # 승수: 제약이 활성일 때만 양수
    x = p - mu / 2 * a                              # stationarity: 2(x - p) + mu a = 0
    slack = c - a @ x                               # 여유: 0 이상
    print(f"c = {c:.0f}  x* = {x}  mu = {mu:.1f}  여유 = {slack:.1f}  "
          f"mu*여유 = {mu * slack:.1f}  f* = {(x - p) @ (x - p):.2f}")

c = 0  x* = [ 0.5 -0.5]  mu = 3.0  여유 = 0.0  mu*여유 = 0.0  f* = 4.50
c = 1  x* = [1. 0.]  mu = 2.0  여유 = 0.0  mu*여유 = 0.0  f* = 2.00
c = 2  x* = [1.5 0.5]  mu = 1.0  여유 = 0.0  mu*여유 = 0.0  f* = 0.50
c = 3  x* = [2. 1.]  mu = 0.0  여유 = 0.0  mu*여유 = 0.0  f* = 0.00
c = 4  x* = [2. 1.]  mu = 0.0  여유 = 1.0  mu*여유 = 0.0  f* = 0.00

표의 각 줄이 KKT 조건의 한 경우다.

• $c = 0, 1, 2$: $\mathbf{p}$가 반공간 밖에 있어 제약이 활성이다(여유 0). 승수 $\mu$는 양수이고, 제약이 느슨해질수록($c$가 클수록) 줄어든다. $c = 1$은 10.1장의 예와 같은 해 $(1, 0)$이다.
• $c = 4$: $\mathbf{p}$가 이미 반공간 안에 있다. 여유가 1이므로 상보적 여유 때문에 $\mu = 0$이다.
• $c = 3$: 경계 사례다. 제약이 등호로 성립하지만 $\mu = 0$이다. 제약이 해에 닿아 있지만 밀지는 않는다.

등식 제약의 경우처럼 승수는 민감도다. 이 문제의 최적값은 $f^\star(c) = (3 - c)^2/2$($c \le 3$)이므로 $df^\star/dc = -(3 - c) = -\mu$다. 표에서 $c$를 1 늘릴 때마다 $f^\star$가 줄어드는 정도가 승수의 크기를 따라간다.

쌍대성

Lagrangian을 $\mathbf{x}$에 대해 먼저 최소화하면 승수만의 함수가 남는다. 이를 쌍대 함수(dual function)라 한다.

쌍대 함수에는 원 문제가 볼록이든 아니든 성립하는 두 성질이 있다.

• 쌍대 함수는 항상 오목하다. $\boldsymbol{\lambda}, \boldsymbol{\nu}$에 대한 아핀 함수들의 점별 최솟값이기 때문이다.
• 약쌍대성(weak duality): $\boldsymbol{\lambda} \ge \mathbf{0}$이면 $g(\boldsymbol{\lambda}, \boldsymbol{\nu}) \le p^\star$. feasible point $\mathbf{x}$에서는 $\lambda_i g_i(\mathbf{x}) \le 0$, $\nu_j h_j(\mathbf{x}) = 0$이므로 $L(\mathbf{x}, \boldsymbol{\lambda}, \boldsymbol{\nu}) \le f(\mathbf{x})$이고, $g$는 그 최솟값이므로 더 작다.

그래서 쌍대 함수의 값은 언제나 최적값의 하한이다. 가장 좋은 하한을 찾는 문제

를 쌍대 문제(dual problem)라 하고, 그 최적값 $d^\star$와 $p^\star$의 차이 $p^\star - d^\star \ge 0$을 쌍대 간극(duality gap)이라 한다. 간극이 0이면 강쌍대성(strong duality)이 성립한다고 한다. 볼록 문제에서는 Slater 조건(모든 부등식 제약을 엄격하게 만족하는 feasible point가 존재) 같은 가벼운 조건 아래 강쌍대성이 성립한다. 볼록이 아닌 문제에서는 일반적으로 간극이 있다.

앞 절의 등식 제약 문제($\min x^2 + 2y^2$ s.t. $x + y = 1$)에서 쌍대 함수를 직접 구하면, $\nabla_{\mathbf{x}} L = \mathbf{0}$에서 $x = -\nu/2$, $y = -\nu/4$이고 이를 $L$에 넣어 $g(\nu) = -\frac{3}{8}\nu^2 - \nu$를 얻는다.

# cell 1의 문제: g(nu) = min_{x,y} L(x, y, nu) = -3 nu^2 / 8 - nu
nus = np.linspace(-4, 2, 601)
g = -3 * nus**2 / 8 - nus
i = np.argmax(g)
print(f"max g(nu) = {g[i]:.4f}  (nu = {nus[i]:.2f})   p* = {p_star:.4f}")
print("모든 nu에서 g(nu) <= p* :", bool(np.all(g <= p_star + 1e-12)))
for nu0 in [-3.0, -1.0, 1.0]:
    print(f"g({nu0:+.0f}) = {-3 * nu0**2 / 8 - nu0:+.4f}  <= p*")

max g(nu) = 0.6667  (nu = -1.33)   p* = 0.6667
모든 nu에서 g(nu) <= p* : True
g(-3) = -0.3750  <= p*
g(-1) = +0.6250  <= p*
g(+1) = -1.3750  <= p*

$g$는 오목한 이차함수이고, 모든 $\nu$에서 $p^\star = 2/3$ 아래에 있다(약쌍대성). 최댓값은 $\nu = -4/3$에서 정확히 $2/3$이다(강쌍대성). 이 최대점의 $\nu$는 첫 절에서 구한 Lagrange 승수와 같다.

쌍대성이 연구에서 쓰이는 방식은 셋이다. 첫째, 하한 인증서다. 어떤 $\boldsymbol{\lambda} \ge \mathbf{0}$에서든 $g$를 계산하면 $p^\star$의 하한을 얻으므로, 현재 해 $f(\mathbf{x})$와 $g$의 차이로 "최적값에서 최대 얼마나 떨어져 있는가"를 보증할 수 있다. 볼록 최적화 solver가 멈출 때 보고하는 duality gap이 이것이다. 둘째, 다른 형태의 문제다. 서포트 벡터 머신(SVM)의 쌍대 문제는 데이터 사이의 내적만으로 쓰이므로 커널 방법(16.3장)으로 이어진다. 셋째, 승수의 해석이다. 강쌍대성이 성립하면 쌍대 문제의 해가 곧 각 제약의 민감도다.

제약과 penalty

제약을 직접 다루는 알고리즘은 복잡하다. 그래서 흔히 제약을 목적함수의 항으로 옮겨 무제약 문제로 바꾼다. 방법은 세 가지가 대표적이다.

• penalty 방법: 제약 위반에 벌점을 준다. 등식 제약은 $f(\mathbf{x}) + \frac{\rho}{2}\sum_j h_j(\mathbf{x})^2$, 부등식 제약은 $\frac{\rho}{2}\sum_i \max(0, g_i(\mathbf{x}))^2$를 더한다. feasible set 밖에서 출발해도 되지만, 유한한 $\rho$의 해는 제약을 정확히 만족하지 않는다.
• barrier 방법: 부등식 제약의 경계에 다가갈수록 무한대로 커지는 항, 예를 들어 $-\frac{1}{t}\sum_i \log(-g_i(\mathbf{x}))$를 더한다. 해가 항상 feasible set 안쪽에 머물고, $t$를 키우며 경계에 다가간다. 내부점 방법(interior-point method)의 기초다.
• augmented Lagrangian: Lagrangian에 이차 penalty를 더한 $L_\rho(\mathbf{x}, \boldsymbol{\nu}) = f + \sum_j \nu_j h_j + \frac{\rho}{2}\sum_j h_j^2$를 $\mathbf{x}$에 대해 최소화하고, 승수를 $\nu_j \leftarrow \nu_j + \rho\,h_j(\mathbf{x})$로 갱신하기를 반복한다. 승수가 제약을 맞추는 역할을 맡으므로 $\rho$를 무한히 키우지 않아도 된다.

첫 절의 문제를 penalty 방법과 augmented Lagrangian으로 풀어 제약 위반과 Hessian 조건수를 비교한다.

# cell 1의 문제를 penalty로: min x^2 + 2y^2 + (rho/2)(x + y - 1)^2 + nu (x + y - 1)
def penalized(rho, nu=0.0):
    H = np.array([[2 + rho, rho], [rho, 4 + rho]])          # (2, 2) Hessian
    x, y = np.linalg.solve(H, [rho - nu, rho - nu])
    return x + y - 1, np.linalg.cond(H)                     # 제약 위반, 조건수

for rho in [1, 10, 100, 1e4, 1e6]:
    viol, cond = penalized(rho)
    print(f"penalty  rho = {rho:>9.0f}  위반 = {viol:+.2e}  조건수 = {cond:.1e}")

nu = 0.0
for k in range(1, 6):
    viol, cond = penalized(10.0, nu)
    nu += 10.0 * viol                                       # 승수 갱신
    print(f"AL 반복 {k}  rho = 10  위반 = {viol:+.2e}  조건수 = {cond:.1e}  nu = {nu:.4f}")

penalty  rho =         1  위반 = -5.71e-01  조건수 = 2.1e+00
penalty  rho =        10  위반 = -1.18e-01  조건수 = 7.8e+00
penalty  rho =       100  위반 = -1.32e-02  조건수 = 6.8e+01
penalty  rho =     10000  위반 = -1.33e-04  조건수 = 6.7e+03
penalty  rho =   1000000  위반 = -1.33e-06  조건수 = 6.7e+05
AL 반복 1  rho = 10  위반 = -1.18e-01  조건수 = 7.8e+00  nu = -1.1765
AL 반복 2  rho = 10  위반 = -1.38e-02  조건수 = 7.8e+00  nu = -1.3149
AL 반복 3  rho = 10  위반 = -1.63e-03  조건수 = 7.8e+00  nu = -1.3312
AL 반복 4  rho = 10  위반 = -1.92e-04  조건수 = 7.8e+00  nu = -1.3331
AL 반복 5  rho = 10  위반 = -2.25e-05  조건수 = 7.8e+00  nu = -1.3333

penalty 방법에서 위반은 정확히 $-4/(4 + 3\rho)$로, $\rho$를 100배 키울 때마다 약 100배 줄어든다. 그러나 $\rho = 10^6$에서도 0이 아니다. 같은 비율로 Hessian의 조건수가 커진다. penalty 항이 제약의 법선 방향에만 $\rho$ 크기의 곡률을 더하기 때문이다. 10.3장에서 보듯 경사하강의 반복 수는 조건수에 비례하므로, 큰 $\rho$는 수렴을 그만큼 느리게 만든다. augmented Lagrangian은 $\rho = 10$에 고정한 채로 반복마다 위반을 약 8분의 1로 줄이고, 승수 $\nu$는 첫 절의 Lagrange 승수 $-4/3$으로 수렴한다.

정규화 항과 노름 제약의 동치. 머신러닝의 정규화는 penalty 방법과 모양이 같다. Ridge 회귀

와 노름 제약 문제

를 비교한다. 제약 문제의 Lagrangian은 $L(\mathbf{w}) + \frac{\mu}{2}(\lVert \mathbf{w} \rVert^2 - r^2)$이고(승수를 $\mu/2$로 잡았다), stationarity 조건은 $\nabla L(\mathbf{w}) + \mu\mathbf{w} = \mathbf{0}$이다. 이는 Ridge의 1계 조건 $\nabla L(\mathbf{w}) + \lambda\mathbf{w} = \mathbf{0}$과 $\mu = \lambda$일 때 똑같다. 따라서 $L$이 볼록일 때, Ridge 해 $\mathbf{w}_\lambda$는 $r = \lVert \mathbf{w}_\lambda \rVert$인 제약 문제의 해이고, 그 KKT 승수는 $\lambda$다. 제약이 활성이므로 상보적 여유도 만족한다. 거꾸로 제약이 활성인 제약 문제의 해는 승수 $\mu$를 $\lambda$로 둔 Ridge의 해다. 큰 $\lambda$는 작은 반지름 $r$에 대응한다.

다만 두 형태는 같은 해 집합을 주지만 같은 하이퍼파라미터는 아니다. $\lambda$와 $r$의 대응 관계는 데이터에 따라 달라진다. 데이터를 바꾸면 같은 $\lambda$가 다른 $r$에 대응한다. 물리적으로 의미 있는 상한(예: 파라미터 크기의 사전 지식)이 있다면 제약 형태로 쓰는 편이 해석하기 쉽고, 하이퍼파라미터를 교차검증으로 고를 거라면 penalty 형태가 계산하기 쉽다.

투영과 근위 연산

제약을 penalty로 바꾸지 않고 직접 다루는 1차 방법도 있다. 두 가지가 기본이다.

Projected gradient. feasible set $\mathcal{C}$가 볼록이고 그 위로의 투영 $\Pi_{\mathcal{C}}(\mathbf{v}) = \arg\min_{\mathbf{x} \in \mathcal{C}} \lVert \mathbf{x} - \mathbf{v} \rVert$를 쉽게 계산할 수 있으면, gradient step 뒤에 투영한다.

상자 $[l, u]$로의 투영은 좌표별로 잘라 내기(clip), 반지름 $r$인 공으로의 투영은 노름이 $r$보다 크면 $r/\lVert \mathbf{v} \rVert$배로 줄이기, 반공간으로의 투영은 앞 절에서 KKT로 구한 식이다. 반복마다 해가 feasible하다는 것이 장점이다. 실습에서 Ridge의 노름 제약 형태를 이 방법으로 푼다.

근위 연산. 목적함수가 매끄러운 부분 $f$와 매끄럽지 않은 부분 $R$의 합 $f(\mathbf{x}) + R(\mathbf{x})$일 때, $R$의 근위 연산(proximal operator)

를 쓴다. "$\mathbf{v}$에서 너무 멀리 가지 않으면서 $R$을 줄이는 점"이다. $R$이 집합 $\mathcal{C}$의 지시함수(안에서 0, 밖에서 $+\infty$)이면 근위 연산은 투영과 같다. 즉 근위 연산은 투영의 일반화다. proximal gradient 방법은 $f$에 대해 gradient step, $R$에 대해 근위 연산을 번갈아 적용한다.

$R(\mathbf{x}) = \lambda\lVert \mathbf{x} \rVert_1$이면 근위 연산은 좌표별 soft thresholding이다.

절댓값이 $t\lambda$ 이하인 성분은 정확히 0이 되고, 나머지는 $t\lambda$만큼 0 쪽으로 줄어든다. 이 반복을 Lasso에 적용한 것이 ISTA(iterative shrinkage-thresholding algorithm)다. 20개 계수 중 3개만 0이 아닌 참 모델에서 생성한 데이터로 확인한다.

rng = np.random.default_rng(0)
A = rng.normal(size=(50, 20))                               # (N, D)
w_true = np.zeros(20)
w_true[[2, 7, 11]] = [1.5, -2.0, 1.0]
y = A @ w_true + 0.1 * rng.normal(size=50)                  # (N,)

def soft(v, t):                                             # prox of t * ||.||_1
    return np.sign(v) * np.maximum(np.abs(v) - t, 0.0)

lam = 5.0
step = 1 / np.linalg.norm(A, 2) ** 2                        # 1 / L, L = sigma_max(A)^2
w = np.zeros(20)
for _ in range(500):                                        # ISTA: 0.5||Aw - y||^2 + lam ||w||_1
    w = soft(w - step * A.T @ (A @ w - y), step * lam)
w_ridge = np.linalg.solve(A.T @ A + lam * np.eye(20), A.T @ y)
print("L1 해의 0이 아닌 좌표 :", np.flatnonzero(w), " 값", w[w != 0].round(3))
print("L1 해의 정확한 0 개수 :", int(np.sum(w == 0)), "/ 20")
print("Ridge 해의 정확한 0 개수:", int(np.sum(w_ridge == 0)), "/ 20")

L1 해의 0이 아닌 좌표 : [ 2  7 11]  값 [ 1.355 -1.823  0.873]
L1 해의 정확한 0 개수 : 17 / 20
Ridge 해의 정확한 0 개수: 0 / 20

L1 해는 참 모델의 0이 아닌 좌표 2, 7, 11을 정확히 찾고 나머지 17개를 정확히 0으로 만든다. 같은 $\lambda$의 Ridge 해에는 정확한 0이 하나도 없다. L2의 근위 연산 $\mathbf{v}/(1 + t\lambda)$은 모든 성분을 같은 비율로 줄일 뿐이기 때문이다. L1 해의 0이 아닌 값(1.355, −1.823, 0.873)은 참값(1.5, −2.0, 1.0)보다 크기가 작다. soft thresholding이 살아남은 계수도 0 쪽으로 줄이기 때문이고, 이를 L1 정규화의 수축 편향(shrinkage bias)이라 한다.

참고 문헌

• Boyd, Vandenberghe (2004). Convex Optimization. Cambridge University Press. — Lagrangian, 쌍대 함수, 약·강 쌍대성, Slater 조건, KKT 조건
• Nocedal, Wright (2006). Numerical Optimization (2nd ed.). Springer. — 제약 최적화의 최적성 조건, penalty·augmented Lagrangian 방법
• Parikh, Boyd (2014). Proximal algorithms. Foundations and Trends in Optimization, 1(3), 127–239. — 근위 연산의 정의와 성질, proximal gradient 방법, Lasso 예