rkm0959

http://www.secmem.org/blog/2021/06/15/SCLI_Framework_and_Applications/

http://www.secmem.org/blog/2021/04/15/Cataylst_Framework/   

http://www.secmem.org/blog/2021/05/04/VarianceReductionCatalyst/ 

github.com/rkm0959/rkm0959_presents/blob/main/acceleration_convex_optimization.pdf

빠르게 Nesterov's Estimate Sequence를 정리한다. 출처는 "Introductory Lectures on Convex Optimization"

정의 : $\{\phi_k(x)\}$와 $\{\lambda_k\}$, $\lambda_k \ge 0$이 함수 $f(x)$의 estimate sequence라는 것은, $\lambda_k \rightarrow 0$이고 임의의 $x \in \mathbb{R}^n$과 $k \ge 0$에 대하여 $$\phi_k(x) \le (1-\lambda_k) f(x) + \lambda_k \phi_0(x)$$라는 것이다. 이 경우, 만약 적당한 수열 $\{x_k\}$가 있어 $$f(x_k) \le \phi_k^{*} \equiv \min_{x \in \mathbb{R}^n} \phi_k(x)$$가 성립한다면, $$f(x_k) \le \phi_k(x_{*}) \le (1-\lambda_k) f_{*} + \lambda_k \phi_0(x_{*})$$가 되어 이를 정리하면 결과적으로 $$f(x_k) - f_{*} \le \lambda_k (\phi_0(x_{*}) - f_{*})$$를 얻는다. 그러니 $x_k$를 잘 잡고, $\lambda_k$가 빠르게 $0$으로 수렴하기만 하면 된다.

이제 estimate sequence를 하나 만들어보자. $f \in \mathcal{F}_{\mu, L}$이라 가정한다. 즉, $f$는 $\mu$-strongly convex, $L$-smooth function.

$\phi_0$를 임의의 함수, $\{y_k\}$를 임의의 $\mathbb{R}^n$ 위 수열, $\{\alpha_k\}$를 $\alpha_k \in (0, 1)$과 $\sum_{k} \alpha_k \rightarrow \infty$를 만족하는 수열이라 하자. $\lambda_0 = 1$이라 하자.

이제 $\{\phi_k\}$와 $\{\lambda_k\}$를 다음과 같이 재귀적으로 정의한다. $$\begin{equation*} \begin{aligned} \lambda_{k+1} &= (1 - \alpha_k) \lambda_k \\ \phi_{k+1}(x) & = (1-\alpha_k)\phi_k(x) + \alpha_k(f(y_k) + \langle \nabla f(y_k), x - y_k \rangle + \frac{\mu}{2} ||x-y_k||^2) \end{aligned} \end{equation*}$$

이는 estimate sequence가 된다. 증명은 $\phi_{k+1}(x) \le (1-\alpha_k) \phi_k(x) + \alpha_k f(x)$임을 활용하여 귀납법. 특히, $\sum_k \alpha_k \rightarrow \infty$ 조건에서 $\lambda_k \rightarrow 0$이 나온다. 이렇게 되면 estimate sequence의 큰 family를 하나 얻었으니, 여기서 입맛에 맞게 $\lambda_k$가 빠르게 감소하도록 $\{\alpha_k\}$와 $\{y_k\}$를 잘 잡아주면 되겠다.

$\phi_0$도 잡아줘야 하는데, 뒤에 계속해서 quadratic function과 섞어먹고 있으니, $\phi_0$도 quadratic하게 잡자.

$\phi_0 = \phi_0^{*} + \frac{\gamma_0}{2} ||x-v_0||^2$이라고 정의하면, $\phi_k$가 전부 quadratic function이 된다. 귀찮은 계산을 하면, $$ \begin{equation*} \begin{aligned} \gamma_{k+1} &= (1-\alpha_k) \gamma_k + \alpha_k \mu \\ v_{k+1} &= \frac{1}{\gamma_{k+1}} \left( (1-\alpha_k) \gamma_k v_k + \alpha_k \mu y_k - \alpha_k \nabla f(y_k) \right) \\ \phi_{k+1}^{*} &= (1-\alpha_k) \phi_k^{*} + \alpha_k f(y_k) - \frac{\alpha^2_k}{2 \gamma_{k+1}} ||\nabla f(y_k)||^2 \\ & + \frac{\alpha_k(1-\alpha_k) \gamma_k}{\gamma_{k+1}} \left( \frac{\mu}{2} ||y_k - v_k||^2 + \langle \nabla f(y_k), v_k - y_k \rangle \right) \end{aligned} \end{equation*}$$라 하면 $\phi_k(x) = \phi_k^{*} + \frac{\gamma_k}{2} ||x - v_k||^2$이 된다.

이 framework에 속하는 알고리즘에서 계산하기 가장 어렵게 생긴 것은 $x_k$를 설정하는 단계이다. 첫 $x_0$야 $\lambda_0 = 1$ 조건에서 아무렇게나 잡아도 충분하니, $x_k$에서 $x_{k+1}$로 쉽게 넘어가는 방법이 있으면 좋을 것 같다. $\phi_k^{*} \ge f(x_k)$를 가정하자. 위 식과 $f$의 convexity를 활용하면, $$\begin{equation*} \begin{aligned} \phi_{k+1}^{*} & \ge (1-\alpha_k) f(x_k) + \alpha_k f(y_k) - \frac{\alpha_k^2}{2\gamma_{k+1}} ||\nabla f(y_k)||^2 + \frac{\alpha_k (1-\alpha_k) \gamma_k}{\gamma_{k+1}} \langle \nabla f(y_k), v_k - y_k \rangle \\ & \ge f(y_k) - \frac{\alpha_k^2}{2\gamma_{k+1}} ||\nabla f(y_k)||^2 + (1-\alpha_k) \langle \nabla f(y_k), x_k - y_k + \frac{\alpha_k \gamma_k}{\gamma_{k+1}} (v_k - y_k) \rangle \end{aligned} \end{equation*}$$ 그러니 우리는 $$f(y_k) - \frac{\alpha_k^2}{2\gamma_{k+1}} ||\nabla f(y_k)||^2 + (1-\alpha_k) \langle \nabla f(y_k), x_k - y_k + \frac{\alpha_k \gamma_k}{\gamma_{k+1}} (v_k - y_k) \ge f(x_{k+1})$$인 $x_{k+1}$이 잡기 쉬운 값이었으면 좋겠다. 그런데 우리는 이미 강력한 Descent Lemma, 즉 $$f \left(x - \frac{1}{L} \nabla f(x) \right) \le f(x) - \frac{1}{2L} || \nabla f(x)||^2$$이라는 식을 알고 있으니, 이를 위 부등식에 끼워맞추면 좋겠다는 생각을 할 수 있다. 이러면 필요한 식이 $$\frac{\alpha_k^2}{2 \gamma_{k+1}} = \frac{1}{2L}, \quad x_k - y_k + \frac{\alpha_k \gamma_k}{\gamma_{k+1}} (v_k - y_k) = 0, \quad x_{k+1} = y_k - \frac{1}{L} \nabla f(y_k)$$

이를 모두 순서에 맞게 잘 합치면, 다음과 같은 iteration을 얻는다. $\gamma_0 > 0$과 $x_0$를 임의로 잡고, $v_0 = x_0$라 하자.

$$ \begin{equation*} \begin{aligned} \alpha_k \in (0, 1) &: L\alpha_k^2 = (1-\alpha_k) \gamma_k + \alpha_k \mu \\ \gamma_{k+1} &= (1-\alpha_k)\gamma_k + \alpha_k \mu \\ y_k &= \frac{\alpha_k \gamma_k v_k + \gamma_{k+1} x_k}{\gamma_k + \alpha_k \mu} \\  x_{k+1} &= y_k - \frac{1}{L} \nabla f(y_k) \\  v_{k+1} &= \frac{1}{\gamma_{k+1}} \left( (1-\alpha_k)y_k v_k + \alpha_k \mu y_k - \alpha_k \nabla f(y_k) \right) \end{aligned} \end{equation*}$$를 반복한다. 여기서 $\gamma_0 = L$을 주고 $\lambda_k$의 값에 대한 bound를 열심히 계산하면, 이 방법이 acceleration을 보임을 확인할 수 있다. 

이제 복잡한 식을 간단하게 하기 위해서 열심히 식 조작을 하면 

$$ \begin{equation*} \begin{aligned} v_{k+1} &= x_k + \frac{1}{\alpha_k} (x_{k+1} - x_k) \\ \alpha_{k+1}^2 &= (1-\alpha_{k+1}) \alpha_k^2 + \alpha_{k+1} \cdot \mu /L \\ \beta_k &= \frac{\alpha_k (1-\alpha_k)}{\alpha_k^2 + \alpha_{k+1}} \\ y_{k+1} &= x_{k+1} + \beta_k (x_{k+1} - x_k) \end{aligned} \end{equation*} $$을 얻는다. 결국 모든 식이 다 사라지고 남는 것은 $\alpha_k, x_k, y_k$ 뿐이다. 결과는 $$\begin{equation*} \begin{aligned} x_{k+1} &= y_k - \frac{1}{L} \nabla f(y_k) \\ \alpha_{k+1}^2 &=  (1-\alpha_{k+1})\alpha_k^2 + \alpha_{k+1} \cdot \mu / L \\ \beta_k &= \frac{\alpha_k (1-\alpha_k)}{\alpha_k^2 + \alpha_{k+1}} \\ y_{k+1} &= x_{k+1} + \beta_k (x_{k+1} - x_k) \end{aligned} \end{equation*} $$ 특히 $\alpha_0 \ge \sqrt{\mu/L}$이면 acceleration이 일어난다. 또한, $\alpha_0 = \sqrt{\mu/L}$이면 $\alpha, \beta$의 값이 고정되어 fixed momentum인 경우가 나타난다.

물론 아예 $\alpha_0 = \sqrt{\mu / L}$로 잡아버리면 $\mu = 0$인 경우에 위험하므로 주의하자. 이렇게 해서 유도 과정이 끝난다.

다음에 기회가 되면 여기서 아이디어를 일부 얻어간 Catalyst에 대해 알아보고자 한다.

이어서 간다. 이번 글에서 다룰 주제는

• OSPEP를 이용한 DRS의 tight contraction factor
• OSPEP를 이용한 optimal parameter selection
• CVXPY를 이용한 위 주제들과 관련된 구현 몇 가지

Part 3. Tight Analytic Contraction Factors for DRS

결과부터 소개한다. 우리가 앞에서 준비했던, 논문의 정리 4.1.만 소개한다. 그래도 논문 링크를 열고 정리 4.1.~4.4.를 한 번 확인하는 것을 추천.

결과도 매우 중요하지만, 아무래도 일단 관심이 가는 것은 대체 이렇게 복잡한 것을 어떻게 증명했냐는 것이다. 논문의 말을 옮기면,

• 기본적으로 Mathematica 같은 Computer Algebra System은 필요하다
• OSPEP에서 얻은 Primal/Dual Problem을 활용한다. 이들을 풀기 위해서, 다음과 같은 접근을 사용한다.
• Primal : 핵심 아이디어는 worst-case operator가 $\mathbb{R}^2$ 위에서 존재할 것이라는 추측이다. 이는 Primal Problem의 해 $G^{*}$의 rank가 2 이하라는 것과 동치인데, 이 추측의 근거는 complementary slackness라고 한다. 이러한 2차원 위에서의 worst-case analysis는 SDP 대신 non-convex quadratic programming 문제로 분석할 수 있고, 이 문제에 대응되는 Karush–Kuhn–Tucker condition을 이용하여 문제를 해결했다고 한다.
• Dual : 역시 여기서도 Dual Problem의 해로 얻어지는 행렬 $S^{*}$가 rank 2 이하라는 가정을 한다. 근거는 Primal과 동일하다. 이 가정으로 $\rho^2$과 $\lambda_\mu^A$, $\lambda_\beta^B$ 사이의 관계식을 추가적으로 얻고, 이를 constraint로 하여 $\rho^2$을 최소화했다고 한다. 말이 쉽지만 정말 어려워 보인다..
• 이렇게 Primal/Dual 문제를 (analytical) 푼 다음, 그 결과가 같음을 확인하면 자동으로 증명이 끝난다.

Appendix를 구경하면 이 결과의 엄밀한 증명을 직접 볼 수 있다. ㄹㅇ 가슴이 웅장해진다... 

Part 4 : Optimal Parameter Selection with OSPEP

여기서는 operator class $\mathcal{Q}_1$, $\mathcal{Q}_2$를 알고 있을 때 parameter $\alpha, \theta$를 optimal 하게 설정하는 방법을 다룬다.

이는 결국 $\rho^2_{*} (\alpha, \theta)$를 해당 $\alpha, \theta$에 대한 optimal contraction factor라고 할 때, $$\rho^2_{*} = \inf_{\alpha, \theta > 0} \rho^2_{*}(\alpha, \theta)$$와 대응되는 $\alpha_{*}$, $\theta_{*}$를 찾는 문제와 같다. 원래 $\theta \in (0, 2)$라는 조건이 있으나, 이 조건은 relax 해도 무방하다. 

$\rho^2_{*}(\alpha, \theta)$의 값 자체는 SDP를 풀어서 계산할 수 있으니, 이차원 위에서 미분 불가능한 함수를 최소화하는 일반적인 알고리즘이 있다면 이를 사용하여 문제를 해결할 수도 있을 것이다. 그러나, 조금 더 작업을 거치면 더욱 편하게 문제를 풀 수 있다. $\theta$를 최적화 문제 안에 집어넣자.

Dual Problem $$\begin{equation*} \begin{aligned} \text{minimize} \quad & \rho^2 \\ \text{subject to} \quad & \lambda_\mu^A, \lambda_\beta^B \ge 0 \\ \quad & S = -M_O - \lambda_\mu^A M_\mu^A - \lambda_\beta^B M_{\beta}^B + \rho^2 M_I \succeq 0 \end{aligned} \end{equation*}$$을 다시 생각해보자. 우리가 $\alpha, \theta$를 고정된 값이 아니라 최적화 문제의 변수로 두지 못하는 이유는, 행렬 $S$의 entry가 $\alpha, \theta$의 nonlinear term을 포함하고 있기 때문이다. $\alpha$는 아예 분모에 있고, $\theta$는 $\theta^2$ 항이 $S$에 있다. 여기서 $\theta$에 대한 quadratic term을 제거하기 위해, Schur Complement를 사용한다.

$$\begin{equation*} \begin{aligned} S &= -M_O - \lambda_\mu^A M_\mu^A - \lambda_\beta^B M_{\beta}^B + \rho^2 M_I \\ S' &= - \lambda_\mu^A M_\mu^A - \lambda_\beta^B M_{\beta}^B + \rho^2 M_I \end{aligned} \end{equation*} $$라고 정의하면, $$M_O = (1, \theta, -\theta)^T (1, \theta, -\theta)$$이므로, Schur Complement의 성질에서 $S \succeq 0$는 $v = (1, \theta, -\theta)^T$라 했을 때 $$\overline{S} = \left[ \begin{matrix} S' & v \\ v^T & 1 \end{matrix} \right] \succeq 0$$과 동치이다. 또한, $S'$은 각 entry에 $\theta$에 대한 non-linear term이 없으므로, $\theta$까지 SDP에 넣을 수 있다.

즉, 새로운 문제를 $$\begin{equation*} \begin{aligned} \text{minimize} \quad & \rho^2 \\ \text{subject to} \quad & \lambda_\mu^A, \lambda_\beta^B, \theta \ge 0 \\ \quad & S' = - \lambda_\mu^A M_\mu^A - \lambda_\beta^B M_{\beta}^B + \rho^2 M_I \\ \quad & v = (1, \theta, -\theta)^T \\ \quad & \overline{S} = \left[ \begin{matrix} S' & v \\ v^T & 1 \end{matrix} \right] \succeq 0 \end{aligned} \end{equation*}$$로 바꾸면, $\theta$에 대한 최적화까지 자동으로 된다. 이제 남은 것은 $\alpha$에 대한 최적화다.

이렇게 해서 얻은 주어진 $\alpha$에 대한 optimal contraction factor를 $\rho^2_{*}(\alpha)$라고 하자. 이를 최소화하면 된다.

변수 하나에 대한 (미분 불가능한) 함수의 minimization은 다양한 방법이 있다. 그런데, 실험 결과에 따르면

•  함수 $\rho^2_{*}(\alpha)$는 연속함수이다. (별로 놀랍지는 않은 사실)
• 함수 $\rho^2_{*}(\alpha)$는 unimodal 함수이다. (놀라운 사실)

그러므로, 이를 최소화하기 위해서 (PS를 했던 사람의 입장에서 쉽게 생각하면) 삼분탐색을 사용해도 좋다.

사실 scipy나 matlab 등에서 이미 derivative-free optimization이 구현이 되어있기 때문에, 이걸 가져다 쓰면 된다.

Part 5 : Programming OSPEP with CVXPY

이건 논문에는 없고 (정확하게는 CVXPY로 하지 않고 MATLAB 등으로 한 것 같다) 그냥 내가 심심해서 하는 것이다.

논문의 저자들이 사용한 코드를 보고 싶다면, github.com/AdrienTaylor/OperatorSplittingPerformanceEstimation로 가자. 

내가 여기서 짠 코드들은 github.com/rkm0959/rkm0959_implements에 가면 있다. 

Task 1 : Verifying Theorem 4.1. Numerically

CVXPY를 이용하여 Theorem 4.1.을 검증해보자. 방법은 다음과 같다.

• $\alpha, \theta, \mu, \beta$를 적당히 랜덤하게 선정한다.  
• 대응되는 Primal/Dual Problem을 CVXPY로 해결한다. SDP를 풀어야하니, blas+lapack을 설치해야 한다.
• 또한, Theorem 4.1.을 이용하여 optimal contraction factor를 계산한다. 
• Theorem 4.1.에는 Case가 총 5가지가 있다. 이때, 어떤 Case에서 계산이 되었는지도 기록한다.
• Primal Problem, Dual Problem, Theorem 4.1.로 얻은 답이 일치하는지 확인한다.
• 위 과정을 여러 번 반복한 후, 5가지의 Case를 각각 몇 번 다루었는지를 확인한다.

이를 구현한 코드는 다음과 같다. 

이렇게 하여 Theorem 4.1.의 결과를 검수할 수 있다. SDP solver가 문제를 잘 해결할 수 있도록 극단적인 값들은 넣지 않았다.

각 코드의 의미는 거의 자명하니 (식을 그대로 옮긴 것에 불과하다) 자세한 설명은 생략하도록 하겠다. Task 1은 이렇게 끝.

Task 2 : Optimal Parameter Selection with OSPEP

이번에는 Part 4에서 다루었던 내용을 구현해보자. 우리가 할 일은

• 주어진 $\mu, \beta$에 대하여 optimal parameter $\alpha, \theta$를 계산
• $\rho^2_{*} (\alpha)$의 그래프를 그리고, 실제로 unimodal function인지 눈으로 확인
• 물론, 이러한 내용 및 그래프는 논문에도 있으며, DRS가 아니라 DYS로 설명되어 있다.

Part 4의 내용을 그대로 CVXPY로 옮기면, 다음과 같은 코드를 얻을 수 있다.

여기서 얻는 결론은 $\mu = 0.53$, $\beta = 1.35$일 때, optimal parameter는 $$\alpha_{*} = 1.5949751, \quad \theta_{*} = 1.4244099$$이고, 이때 얻어지는 optimal contraction factor는 $$\rho_{*} = 0.5010598$$이라는 것이다. 또한, 이 경우에 대한 $\rho^2_{*} (\alpha)$의 그래프를 그려서 unimodal function임을 확인할 수 있다. 그래프는 아래와 같다.

당연하지만, $x$축은 $\alpha$의 값, $y$축은 $\rho^2_{*}(\alpha)$의 값이다. 이렇게 Task 2와 논문 정리도 끝 :)

이제 본격적으로 논문을 읽어보자. 논문에서는 예시 문제를 Davis-Yin Splitting (DYS)로 두고 문제를 해결한다.

하지만 이 글에서는 예시 문제를 Douglas-Rachford Splitting (DRS)로 두고 해결하고자 한다. 이유는

• 일단 나도 손계산 연습을 해보고 싶고, operator가 3개인 DYS보다는 2개인 DRS가 간단함
• 앞선 글에서 DRS를 어느 정도 설명했으니 독자가 이해하기도 편할 것
• 이 논문의 Theorem 4에서 DRS에 대한 정리를 하는데, 직접 코드를 짜서 verify 해보고 싶음
• 그리고 솔직히 논문 내용 그대로 옮기면 재미없음 ㄹㅇㅋㅋ 대신 흐름은 그대로 유지할 것

이번 글에서는 다음 주제를 다룬다.

• Operator Splitting Performance Estimation Problem을 SDP로 변환하는 방법
• OSPEP로 얻어진 SDP의 Dual Problem과 그 의미 및 활용법

Part 1 : Operator Splitting Performance Estimation Problem to SDP

먼저 DRS를 다시 볼 필요가 있다. 우리는 $R_{\alpha A}R_{\alpha B}$가 nonexpansive operator이므로, 이를 $I$와 함께 섞어서 $$\frac{1}{2} I + \frac{1}{2} R_{\alpha A} R_{\alpha B}$$로 fixed point iteration을 돌렸다. 하지만 굳이 이렇게 섞을 필요는 당연히 없고, $\theta \in (0, 2)$에 대하여 $$\frac{2 - \theta}{2} I + \frac{\theta}{2} R_{\alpha A} R_{\alpha B}$$로 섞어도 된다. 여기에 $R_{\alpha A} = 2J_{\alpha A} - I$, $R_{\alpha B} = 2J_{\alpha B} - I$를 대입하면 $$I - \theta J_{\alpha B} + \theta J_{\alpha A} (2J_{\alpha B} - I)$$를 우리의 operator로 얻는다. 이제 목표는 이 operator가 "얼마나 contractive"한지 파악하는 것이다.

우리는 Performance Estimation Problem의 작은 문제에서 $\mathcal{F}_{\mu ,L}$에서 문제를 생각했다.

이번에 다루어 볼 operator class를 알아보자. 앞에서 정의한 네 가지 operator class들과 이들의 교집합을 생각한다.

• maximal monotone operator의 집합 $\mathcal{M}$
• $\mu$-strongly monotone이면서 maximal monotone인 operator의 집합 $\mathcal{M}_\mu$
• $L$-Lipschitz operator의 집합 $\mathcal{L}_L$
• $\beta$-cocoercive operator의 집합 $\mathcal{C}_{\beta}$

이제 본격적으로 문제를 서술할 준비가 되었다. $\mathcal{Q}_1, \mathcal{Q}_2$가 각각 operator class고, (단, 이들은 $\mathcal{M}$의 subset) $$T(z ; A, B, \alpha, \theta) = I - \theta J_{\alpha B} + \theta J_{\alpha A} (2J_{\alpha B} - I)$$이라 하자. 우리가 풀고자 하는 최적화 문제는, 고정된 $\alpha, \theta > 0$에 대하여, $$\begin{equation*} \begin{aligned} \text{maximize} \quad & \frac{ ||T(z;A,B,\alpha, \theta) - T(z';A,B,\alpha,\theta)||^2}{||z-z'||^2} \\ \text{subject to} \quad & A \in \mathcal{Q}_1, B \in \mathcal{Q}_2, z \neq z' \end{aligned} \end{equation*}$$이다. 즉, 최악의 경우에서도 $T$를 $L$-Lipschitz라고 부를 수 있는 가장 최적의 $L$을 구하는 것이다.

이제 문제를 조금씩 변형해가며 SDP로 만들어보자. Performance Estimation Problem과 크게 다르지 않다.

Step 1 : $T(z; A, B, \alpha, \theta)$를 표현하기 위해서, 다음과 같이 변수를 설정한다. $$z_B = J_{\alpha B} z, \quad z_A = J_{\alpha A}(2z_B - z), \quad T(z;A,B,\alpha, \theta) = z - \theta (z_B - z_A)$$

Step 2 : 분모에 있는 $||z-z'||^2$가 불편하다. 이를 일반성을 잃지 않고 $1$로 강제하고 싶은 마음이 자연스럽다.

다행히도 이는 가능하다. Operator class $\mathcal{Q}$가 임의의 $\gamma > 0$에 대하여 $$A \in \mathcal{Q} \iff (\gamma^{-1}I)A(\gamma I) \in \mathcal{Q}$$라면 $\mathcal{Q}$를 homogeneous 하다고 한다. 우리가 다루는 operator class는 전부 homogeneous.

이제 $\gamma = ||z-z'||$라고 하면, $z$들을 전부 $\gamma^{-1}z$로 교체하고 $A, B$를 $(\gamma^{-1} I)A(\gamma I)$, $(\gamma^{-1} I)B(\gamma I)$로 교체하자.

이래도 문제 없음은 homogeneous 성질에서 나온다. 또한, Step 1에서 얻은 식 정리를 그대로 사용해도 괜찮다는 것을 열심히 계산하면 알 수 있다.

결국 $||z-z'||=1$인 경우에서만 문제를 해결해도 무방하다는 사실을 얻는다. 제약 조건을 하나 추가하고, 분모를 날리자.

중간 점검을 하자면, 결국 우리가 풀 문제는 $$\begin{equation*} \begin{aligned} \text{maximize} \quad & ||z - \theta (z_B - z_A) - z' + \theta(z'_B - z'_A)||^2 \\ \text{subject to} \quad & A \in \mathcal{Q}_1, B \in \mathcal{Q}_2, ||z-z'||^2 = 1 \\ \quad & z_B = J_{\alpha B} z, z_A = J_{\alpha A}(2z_B - z) \\ \quad & z'_B = J_{\alpha B} z', z'_A = J_{\alpha A}(2z'_B- z') \end{aligned} \end{equation*}$$

Step 3 : 이제 슬슬 operator class $\mathcal{Q}$가 눈에 거슬린다. Interpolation의 아이디어를 적용하면 좋겠다.

뒤에서 Interpolation을 적용할 때 알게 되겠지만, 필요한 것은 두 점에 대한 Interpolation 조건이다.

즉, Graph에 $(x_1, q_1), (x_2, q_2)$가 속하는 operator $T \in \mathcal{Q}$가 존재할 조건만이 필요하다. 결론을 제시하면,

• $\mathcal{M}$ : $\langle q_1 - q_2, x_1 - x_2 \rangle \ge 0$
• $\mathcal{M}_\mu$ : $\langle q_1 - q_2 , x_1 - x_2 \rangle \ge \mu ||x_1 - x_2||^2$
• $\mathcal{C}_{\beta}$ : $\langle q_1 - q_2, x_1 - x_2 \rangle \ge \beta || q_1 - q_2||^2$
• $\mathcal{L}_L$ : $||q_1- q_2||^2 \le L^2 ||x_1 - x_2||^2$
• 단, $L$-Lipschitz 조건을 넣고 싶다면 무조건 $\mathcal{M} \cap \mathcal{L}_L$ 형태로 operator class를 잡아야 한다
• $\mathcal{Q}$가 위 operator class의 교집합이라면, 단순히 대응되는 조건을 모두 가져오면 된다.
• 즉, 두 점이 operator class의 조건을 만족시키면 무조건 interpolation도 된다는 것이 결론. 점의 개수가 많아지면 성립하지 않는다.
• TMI : monotone operator는 일반적인 경우에서도 maximal monotone operator로 interpolate 할 수 있다 (Zorn's Lemma)

논문의 Theorem 4.1을 직접 검증해보기 위해서, 대응되는 세팅인 $\mathcal{Q}_1 = \mathcal{M}_\mu$, $\mathcal{Q}_2 = \mathcal{M} \cap \mathcal{C}_\beta$를 생각해보자.

• $B$가 지나야 하는 점은 $(z_B, (z-z_B)/\alpha)$와 $(z'_B, (z'-z'_B)/\alpha)$
• 그러니 얻는 식은 $\langle ((z-z')-(z_B - z'_B))/\alpha, z_B - z'_B \rangle \ge \beta ||((z-z')-(z_B - z'_B))/\alpha||^2$
• $A$가 지나야 하는 점은 $(z_A, (2z_B - z-z_A)/\alpha)$와 $(z'_A, (2z'_B-z'-z'_A)/\alpha)$
• 그러니 얻는 식은 $\langle (2(z_B - z'_B) - (z-z') - (z_A - z'_A))/\alpha, z_A- z'_A \rangle \ge \mu ||z_A - z'_A||^2$

추가적으로, $||z-z'||^2 = 1$이라는 조건도 가지고 있다. 눈치챘겠지만, $z-z'$, $z_A-z_A'$, $z_B - z'_B$과 같이 두 변수의 차로 식이 표현된다.

이를 치환하고 정리하면, 우리가 원하는 것은 결국 $$\begin{equation*} \begin{aligned} \text{maximize} \quad & ||z - \theta (z_B - z_A)||^2 \\ \text{subject to} \quad & ||z||^2=1 \\ \quad & \langle z-z_B, z_B \rangle \ge (\beta/\alpha) ||z-z_B||^2 \\ \quad & \langle 2z_B - z , z_A \rangle \ge (1+\alpha \mu)||z_A||^2 \end{aligned} \end{equation*}$$를 얻는다. 이제 SDP를 꺼낼 준비가 완료되었다. 마찬가지로 Gram Matrix $$G = \left[ \begin{matrix} \langle z, z \rangle & \langle z_A, z \rangle & \langle z_B, z \rangle \\ \langle z, z_A \rangle & \langle z_A, z_A \rangle & \langle z_B, z_A \rangle \\ \langle z, z_B \rangle & \langle z_A, z_B \rangle & \langle z_B, z_B \rangle \end{matrix} \right]$$를 꺼내자. 작업하고 있는 공간의 차원이 3 이상이라는 가정하에서, 우리가 아는 것은 $G$가 Positive Semidefinite이라는 사실 뿐이다.

이제 다음과 같이 행렬을 정의하자. $$ M_O = \left[ \begin{matrix} 1 & \theta & -\theta \\ \theta & \theta^2 & -\theta^2 \\ -\theta & -\theta^2 & \theta^2 \end{matrix} \right], \quad M_I = \left[ \begin{matrix} 1 & 0 & 0 \\ 0 & 0 & 0 \\ 0 & 0 & 0 \end{matrix} \right]$$ $$ M_\beta^B = \left[ \begin{matrix} -\beta/\alpha & 0 & \beta/\alpha + 1/2 \\ 0 & 0 & 0 \\ \beta/\alpha + 1/2 & 0 & -\beta/\alpha - 1 \end{matrix} \right], \quad M_\mu^A = \left[ \begin{matrix} 0 & -1/2 & 0 \\ -1/2 & -(1+\alpha \mu) & 1 \\ 0 & 1 & 0 \end{matrix} \right]$$ 그러면 최종적으로 우리가 원하는 SDP 형태를 얻는다. $$\begin{equation*} \begin{aligned} \text{maximize} \quad & Tr(M_OG) \\ \text{subject to} \quad & G \succeq 0, Tr(M_IG) = 1 \\ \quad & Tr(M_\mu^A G), Tr(M_\beta^BG) \ge 0 \end{aligned} \end{equation*}$$

특히, 이 SDP formulation의 해를 가지고 "worst case"에 대응되는 operator $A, B$를 복원할 수 있다.

$Tr(M_O G)$가 최대가 되는 경우의 $G$를 찾았다면, Gram Matrix의 성질을 이용하여 $z, z_A, z_B$를 복원한 다음 interpolate 하면 끝이다.

Part 2 : The Dual Problem and its Interpretations/Usage

이제 얻은 SDP의 dual problem을 생각해보자. 전형적인 계산을 하면 결국 $$\begin{equation*} \begin{aligned} \text{minimize} \quad & \rho^2 \\ \text{subject to} \quad & \lambda_\mu^A, \lambda_\beta^B \ge 0 \\ \quad & S = -M_O - \lambda_\mu^A M_\mu^A - \lambda_\beta^B M_{\beta}^B + \rho^2 M_I \succeq 0 \end{aligned} \end{equation*}$$를 얻는다. 또한, "웬만하면" (regularity condition) Strong Duality가 성립함은 Slater's Condition으로 보일 수 있다.

이는 논문의 정리 3.3.으로, 나름 전형적인 증명이고 테크니컬하므로 생략. 중요한 것은 이 Dual Problem을 어떻게 바라보느냐, 그리고 어떻게 써먹느냐다.

눈치를 챘겠지만, notation을 $M_{\mu}^A$, $M_{\beta}^B$로 둔 이유는 사실 이 행렬들, 정확하게는 이 행렬들과 관련된 trace에 대한 조건이 $A$의 $\mu$-strongly monotone, $B$의 $\beta$-cocoercivity를 encode 한 것이기 때문이다. 이와 같이 operator class에 $A, B$가 속하게 되면 이로부터 얻어지는 부등식들이 있고, 우리가 증명할 수 있는 것은 기본적으로 이러한 부등식들의 nonnegative combination이다. 이 생각을 이어가면 결국 얻을 수 있는 결론은

• dual variable $\lambda_\mu^A$, $\lambda_\beta^B$은 $A, B$에 대응되는 부등식들에 대한 가중치이다. (그러니 $0$ 이상이다)
• feasible point에서의 $\rho$ 값은 해당 가중치를 사용하였을 때 얻어지는 contraction factor이며, 이는 optimal 값에 대한 upper bound다.

즉, "가중치를 이렇게 넣었더니 contraction factor가 이렇게 나왔다"는 사실 "dual problem의 feasible point를 찾았다"와 같다.

강조하자면, Primal/Dual Problem의 직관적인 의미는 각각

• Primal : 변수로 들어가는 것은 "최악의 경우에 해당하는 operator", 얻어지는 것은 "이런 contraction factor는 기대할 수 없다"는 결론
• Dual : 변수로 들어가는 것은 "증명에 사용될 가중치", 얻어지는 것은 "이런 contraction factor는 얻어질 수 있다"는 결론
• 특히, Duality Gap이 "웬만하면" (Slater's Condition) 0이므로, 정말 아무런 걱정도 할 필요가 없음

이제 Dual Problem의 의미를 알았으니, 활용하는 법도 파악할 수가 있다.

• Contraction Factor에 대한 Upper Bound를 확보할 수 있다. 애초에 dual problem에서 항상 나오는 결론이다.
• Dual Problem을 풀면 "증명에 사용될 가중치"를 얻을 수 있다! 그러니 여러 parameter에 대해서 문제를 수치적으로 풀어보고, 이를 통해서 어떤 가중치가 사용되는지 확인한 후, 이에 대한 analytical 공식을 때려맞춘뒤 엄밀한 증명을 시작할 수 (정확하게는, 시도할 수) 있다.
• Convex Optimization에서 수렴성 증명 중 상당 부분은 정말 단순 계산인데, 그 계산각을 재는 게 어렵다. (후에 Acceleration을 다루면 Lyapunov Function에서 또 이 이야기가 나올 것이다) 부등식들을 섞어먹을 가중치를 선정하는 것도 난이도가 상당하다. Dual Problem은 여기서 도움을 준다.

글 3개에 걸쳐서 arxiv.org/pdf/1812.00146.pdf 논문을 간단하게 정리한다.

저번 글과 비슷한 제목과 비슷한 아이디어를 갖고 있으나, 그 context가 약간 다르다.

저번 글에서는 Gradient Descent라는 거의 모두가 아는 context가 있어서 설명이 쉬웠으나, 이번에는 약간 다르다. 

일단 Operator Splitting이 무엇인지를 알아야 context를 이해할 수 있다. 그러니 이번 글에서는 이를 빠르게 설명한다.

• Disclaimer : 사실 아래에서 $\mathbb{R}^n$으로 말하는 것은 Real Hilbert Space로 보는 것이 맞는 것으로 알고 있다. 굳이 $\mathbb{R}^n$으로 쓰는 이유는 내가 실해석도 모르고 함수해석도 모르고 여러모로 모르는 게 많기 때문. 올해 공부를 하면 어떻게 해결할 수 있지 않을까 기대중이다.
• 아래 내용의 대부분은 www.math.snu.ac.kr/~ernestryu/courses/optimization.html의 6, 7강에서 공부할 수 있다.

set-valued operator는 $\mathbb{R}^n$의 한 점을 $\mathbb{R}^n$의 subset으로 보내는 operator이다.

만약 set-valued operator의 값이 empty set이거나 크기가 1인 집합이라면, 이를 single-valued operator라고 한다.

Example (Subgradient) : $f : \mathbb{R}^n \rightarrow \mathbb{R} \cup \{\pm \infty\}$가 convex 함수이고 모든 $y \in \mathbb{R}^n$에 대하여 $$f(y) \ge f(x) + \langle g, y-x \rangle$$이 성립하면 $g$를 $f$의 subgradient at $x$라고 부르며, subgradient의 모음을 $\partial f(x)$라고 쓴다. 이는 set-valued operator.

• $0 \in \partial f(x)$라면 $x$가 $f$의 값을 최소화하며, 그 역도 성립한다. 
• $\partial (f+g)(x) \supseteq \partial f(x) + \partial g(x)$임은 자명. 단, 집합끼리의 덧셈은 Minkowski Sum.
• $\text{dom} f \cap \text{int } \text{dom} g \neq \emptyset$이면 $\partial (f+g)(x) = \partial f(x) + \partial g(x)$. 증명은 여기. 단 domain은 $f(x) < \infty$인 점들의 집합.

함수에 대하여 Convex, $L$-smooth, $\mu$-strongly convex 등을 정의하는 것을 operator에도 비슷하게 할 수 있다.

operator $A$가 있다. 각 $x, y \in \mathbb{R}^n$와 $u \in Ax$, $v \in Ay$에 대하여 ($A(x), A(y)$를 $Ax, Ay$라고 간단하게 표기한다)

• $\langle u-v, x-y \rangle \ge 0$이면, $A$를 monotone operator라 한다.
• $\langle u-v, x-y \rangle \ge \mu ||x-y||^2$이면, $A$를 $\mu$-strongly monotone operator라 한다.
• $\langle u-v, x-y \rangle \ge \beta||u-v||^2$이면, $A$를 $\beta$-cocoercive라 한다.
• $||u-v|| \le L||x-y||$라면 $A$를 $L$-Lipschitz라 한다. 물론, $\mu, \beta, L > 0$.
• 특히, $1$-Lipschitz인 operator를 nonexpansive operator라고 한다.

operator $A$에 대응되는 graph를 $\{(x, u) : u \in Ax\}$라고 정의하자.

또한, operator $A$가 다음 조건을 만족하면, $A$를 maximal monotone 하다고 정의한다.

• $A$는 monotone operator
• operator $B$가 다음 조건을 만족하면, $A = B$
• $B$는 monotone operator이며, $B$의 graph는 $A$의 graph를 포함한다.

operator $A$의 inverse를 $\{(u, x) : u \in Ax\}$를 graph로 갖는 operator라고 정의하자.

$A$의 resolvent operator를 $J_A = (I+A)^{-1}$로 정의한다. 물론, $I$는 identity operator이다.

또한, reflected resolvent operator를 $R_A = 2J_A - I$로 정의한다. resolvent는 뒤에서 중요한 역할을 한다.

최적화 문제와의 관련성을 본격적으로 논의해보자. 기본적으로, 다음과 같은 세팅을 거친다.

• 최적화 문제를 operator에 대한 문제로 바꾼다. 예를 들어, $f$의 최솟값을 구하는 문제는 $\partial f$의 zero를 찾는 문제가 된다.
• 이를 operator $T$의 fixed point를 찾는 문제로, 즉 $z \in Tz$인 $z$를 찾는 문제로 바꾼다. 이를 fixed point encoding이라고 부른다.
• operator $T$의 fixed point를 찾기 위해서, $z^{k+1} = Tz^k$를 반복한다. 이것이 수렴하면 성공이다.

operator $T$가 identity operator $I$와 nonexpansive operator $N$의 convex combination이라면 (단, $I$에 대응되는 비율이 nonzero) $T$를 averaged operator라고 부른다. 이 경우, $T$가 fixed point를 가진다면 $z^{k+1} = Tz^k$가 잘 수렴함을 증명할 수 있다. 증명 방식은 나중에 다루도록 하겠다.

특히, $T$가 $L$-Lipschitz고 $L<1$이라면, 수렴 속도가 기하급수적으로 빠르다는 것은 쉽게 증명할 수 있는 사실이다.

이제 fixed point encoding의 예시를 몇 개 보이는 것으로 이론적인 준비를 마치도록 하겠다.

Proximal Point Method : $A$가 maximal monotone이면 다음이 성립한다.

• $J_A$가 averaged operator : $R_A = 2J_A - I$가 nonexpansive 함을 쉽게 보일 수 있다. 특히, $J_A$는 single-valued.
• Minty Surjectivity Theorem : $J_A$의 domain은 $\mathbb{R}^n$이다. 즉, 모든 $x \in \mathbb{R}^n$에 대하여 $J_Ax$는 nonempty.

또한, $f$가 convex, closed, proper (CCP) function이라면, $\partial f$가 maximal monotone operator임을 보일 수 있다.

• Convex : $f(\theta x + (1-\theta)y) \le \theta f(x) + (1-\theta)f(y)$ for $\theta \in [0, 1]$
• Closed : $\{(x, \alpha) \in \mathbb{R}^n \times \mathbb{R} : f(x) \le \alpha\}$ is closed
• Proper : $f(x) \neq -\infty$ always, and $f(x) < \infty$ somewhere

이제 준비를 할 수 있다. $f$가 CCP function이라고 하면, 

• $f$의 최솟값을 구하는 문제를 $\partial f$의 zero를 찾는 문제로 바꿀 수 있다.
• 이를 $J_{\alpha \partial f} = (I + \alpha \partial f)^{-1}$의 fixed point를 찾는 문제로 바꿀 수 있다. 단, $\alpha > 0$.
• 그러니 $x^{k+1} = J_{\alpha \partial f}(x^k)$를 반복하여 문제를 해결할 수 있다. 이것이 Proximal Point Method.

Douglas-Rachford Splitting (DRS) : 이번에는 maximal monotone operator $A, B$가 있다고 가정하고, $$ 0 \in (A+B)x$$인 $x$를 찾는 문제를 생각하자. 열심히 계산하면, 이는 결국 $$x = J_{\alpha B}z, \quad R_{\alpha A}R_{\alpha B} z = z$$인 $(x, z)$를 찾는 문제와 같음을 알 수 있다. ($\alpha > 0$) 그러니 fixed point encoding을 얻는다.

그런데 $R_{\alpha A} R_{\alpha B}$는 nonexpansive 하지만 averaged operator는 아닐 수 있다. 그러므로, 이 문제를 해소하기 위하여 $$z^{k+1} = \left( \frac{1}{2} + \frac{1}{2} R_{\alpha A} R_{\alpha B} \right) z^k$$라고 하면 이는 averaged operator를 사용하였으니 수렴이 보장되는 fixed point iteration이 된다. 이를 통해서 문제를 해결할 수 있다.

이처럼 operator 여러 개가 섞인 문제를 쪼개서 푸는 것을 operator splitting이라고 한다.

예고 : 목표는 Operator Splitting Method에 대해서 Performance Estimation을 하는 것이다. 예를 들어, DRS에서는 iteration을 averaged operator $$ \frac{1}{2} I + \frac{1}{2} R_{\alpha A} R_{\alpha B}$$를 이용하여 했었다. $A, B$가 단순히 maximal monotone이 아니라 추가적인 성질을 가지고 있다면, 이 operator가 $L$-Lipschitz가 ($L<1$) 될 수도 있을 것이다. 이때, 우리가 보장할 수 있는 최선의 $L$ 값은 얼마일까? 그리고 이러한 문제를 어떻게 응용할 수 있을까? 다음 글에서 알아보도록 하자.

이 글에서는 논문 www.di.ens.fr/~ataylor/share/PESTO_CDC_2017.pdf를 정리한다.

(지금 보니까 저 논문은 toolkit 소개고, 아이디어 자체는 arxiv.org/pdf/1502.05666.pdf를 보는 게 더 좋을 것 같다)

보통 우리는 문제가 주어지면

• 이를 해결하는 알고리즘을 설계하고
• 최악의 경우에 어떻게 작동하는지 분석

한다. 만약 PS를 하는 사람이라면 답을 정확하게 구하는 알고리즘을 설계한 후, 최악의 경우의 시간복잡도를 분석할 것이다.

비슷하게, 최적화에서는 보통 답을 iterative 하게 근사하는 알고리즘을 설계한 후, 최악의 경우의 수렴성을 분석한다.

이 논문은 "최악의 경우의 수렴성"을 분석하는 자동화된 toolkit을 제시하는 논문이다. toolkit 링크는 여기다. Matlab으로 작성되어 있다.

toolkit의 기본 작동 원리를 알기 위해서, 작은 문제를 하나 풀어보자.

문제 (Informal) : $\mu$-strongly convex, $L$-smooth function $F$가 있다. $F$를 최소화하는 문제를 풀고자 한다. Initial point $x_0$에서 Gradient Descent를 적용하여 $x_1 = x_0 - \gamma \nabla F(x_0)$를 얻었을 때, $||\nabla F(x_1)||$의 크기를 얼마나 잘 줄일 수 있을까? 최악의 경우를 비교하여, 최선의 $\gamma$를 구할 수 있을까?

• $\mu$-strongly convex, $L$-smooth의 의미는 이 글과 이 글에서 매우 잘 정리되어 있다. 
• $\mu$-strongly convex, $L$-smooth function의 집합을 $\mathcal{F}_{\mu, L}$이라고 쓴다.

이 문제를 엄밀하게 서술하려면, 이를 최적화 문제 형태로 서술해야 한다. $$\begin{equation*} \begin{aligned} \text{maximize}  \quad  & ||\nabla F(x_1)||^2 \\  \text{subject to} \quad & F \in \mathcal{F}_{\mu, L} \\ \quad & x_1 = x_0 - \gamma \nabla F(x_0) \end{aligned} \end{equation*}$$ 그런데 이렇게만 쓰면 답이 $\infty$인 것이 뻔하므로, 추가적인 제약조건이 필요하다. 예를 들어, $|| \nabla F(x_0) ||^2$의 크기에 제약조건을 걸 수 있다. 그러면 $$\begin{equation*} \begin{aligned} \text{maximize}  \quad  & ||\nabla F(x_1)||^2 \\  \text{subject to} \quad & F \in \mathcal{F}_{\mu, L} \\ \quad & x_1 = x_0 - \gamma \nabla F(x_0) \\ \quad & ||\nabla F(x_0)||^2 \le R^2 \end{aligned} \end{equation*}$$ 라는 최적화 문제를 얻는다. 그런데 $F \in \mathcal{F}_{\mu, L}$가 변수라는 점이 문제 해결을 어렵게 한다. 

여기서 멋진 아이디어가 사용된다. 먼저, 우리가 $F$라는 함수에 접근하는 점은 $x_0, x_1$ 두 곳이다. 

그러니 사실 우리는 $f_i = F(x_i)$, $g_i = \nabla F(x_i)$의 값에만 관심을 가져도 충분하다. 이제 최적화 문제를 쓰면 $$\begin{equation*} \begin{aligned} \text{maximize}  \quad  & ||g_1||^2 \\ \text{subject to} \quad & \exists F \in \mathcal{F}_{\mu, L} \text{  s.t.  } F(x_i) = f_i, \nabla F(x_i) = g_i \\ \quad & x_1 = x_0 - \gamma g_0, \text{  } ||g_0||^2 \le R^2 \end{aligned} \end{equation*}$$의 형태가 된다. 그러니, 우리가 관심을 갖게 되는 문제는

• 언제 $F(x_i) = f_i$, $\nabla F(x_i) = g_i$를 만족하는 $F \in \mathcal{F}_{\mu, L}$이 존재할까?

가 된다. 이런 문제는 Lagrange Interpolation의 향기가 나고, 실제로 Interpolation이라고 불린다. 

그리고 $F \in \mathcal{F}_{\mu, L}$로 interpolation이 되는 조건은 이 논문에서 이미 증명이 되었다.

결론은 각 $i \neq j$에 대하여 부등식 $$f_i \ge f_j + \langle g_j, x_i - x_j \rangle + \frac{1}{2(1-\mu / L)} \left( \frac{1}{L} ||g_i - g_j||^2 + \mu ||x_i - x_j||^2  - 2 \frac{\mu}{L} \langle g_j - g_i, x_j - x_i \rangle \right)$$이 성립하는 것이, 원하는 interpolation 조건을 만족하는 $F \in \mathcal{F}_{\mu, L}$이 존재할 필요충분조건이라는 것이다.

이제 최적화 문제의 constraint 첫 줄을 interpolation 조건으로 교체할 수 있다. 또한, $x_1 = x_0 - \gamma g_0$를 대입하여 $x_1$ 변수를 제거할 수 있다.

남은 문제를 바라보면 다루고 있는 변수가

• $f_0$, $f_1$ : 이 값들은 constraint에서 그냥 혼자 있는 값들이다.
• $g_0, g_1, x_0$ : 이 값들은 objective나 constraint나 무조건 다른 값과 내적이 되어있다. 

이 성질 때문에, 문제를 Semidefinite Programming (SDP) 형태로 바꿀 수 있다. Gram Matrix $$G = \left[ \begin{matrix} \langle x_0, x_0 \rangle & \langle g_0, x_0 \rangle & \langle g_1, x_0 \rangle \\ \langle x_0, g_0 \rangle & \langle g_0, g_0 \rangle & \langle g_1, g_0 \rangle \\ \langle x_0, g_1 \rangle & \langle g_0, g_1 \rangle & \langle g_1, g_1 \rangle \end{matrix} \right]$$를 생각하자. Gram Matrix의 성질에 의하여, 다음이 성립한다.

• $G$는 Positive Semidefinite이며, 특히 풀고 있는 문제가 $\mathbb{R}^d$ 위의 문제면 rank가 최대 $d$.
• Positive Semidefinite이며 rank가 최대 $d$인 행렬이 주어지면, 이를 Gram Matrix로 보고 vector $x_0, g_0, g_1$을 복원할 수 있음.

그러므로 우리의 모든 문제를 $f_0$, $f_1$과 $G$에 대한 문제로 볼 수 있으며, $G$에 대해 명시할 조건은 $G$가 Positive Semidefinite이라는 점과 rank가 $d$ 이하라는 점이다. 하지만 우리의 관심사는 일반적인 $d$에 대한 것이므로, rank 조건은 제거해도 무방하다. (당장 $d \ge 3$이기만 해도 rank 조건이 무의미하다)

이제 남은 문제를 보면 진짜 단순한 Semidefinite Programming 문제고, 이는 standard solver로 해결할 수 있다.

논문에 의하면, 이렇게 SDP로 문제를 변형하여 해결할 수 있는 경우는 다음을 포함한다.

• Class of Functions : Smooth Convex, Strongly Convex, Bounded Subgradients, Bounded Domains, Non-Convex Smooth 등
• Methods : First-Order Methods : Subgradient, Proximal Operator 등
• Performance Measure (objective 설정) : $||x_n - x_{*}||^2$, $||\nabla F(x_n)||^2$, $F(x_n) - F(x_{*})$ 등
• Initial Condition : Performance Measure와 동일함. 물론, $x_n$ 대신 initial point $x_0$에 대한 조건.

그리고 이 toolkit은 이러한 경우 대부분을 지원한다고 한다. 갖고 놀기 좋은 것 같은데, 아직 쓰는 법은 모른다 ㅎㅎ;

궁금한 점은 이 toolkit으로 답을 $\mu$, $L$에 대한 수식으로 얻으려면 수식이 변수로 들어가있는 Symbolic한 SDP를 해결해야 할텐데, 이건 대체 어떻게 한 건지 궁금하다. toolkit 보면 수식으로 답을 얻어주는 것 같은데... 관련 논문 읽어봐도 아이디어는 알겠으나 이 점이 계속 헷갈린다. 더 공부하는 걸로.

UPD : 직접 돌려본 결과, 수치적으로 계산해준다. 그냥 예시에 답이 이렇게 나와야 한다고 주석으로 써있는 것..

arxiv.org/pdf/1502.05666.pdf 논문의 추가적인 내용은

• Interpolation 관련 엄밀한 증명
• SDP로 문제를 변환하는 것과, 해당 문제의 성질에 대한 엄밀한 증명
• SDP로 변환된 문제의 Lagrangian Dual을 생각하여, Upper Bound를 계산
• 특히, Slater Condition을 이용하여 "웬만하면" Duality Gap이 없음을 증명
• Performance Estimation Problem의 여러 응용 : 이걸로 문제를 많이 해결/제시하는 것 같다