冷めたコーヒー

東京へのアクセスは容易なのでお誘いはお受けします．↩
2 ヶ月前に退去の旨を伝えなければならない物件でした．厳しいですね…．↩
今回は相見積もりサイトに登録して決定しました．これはあまりおすすめできる方法ではありませんが，結果的にかなり低価格で抑えることができたので結果オーライです．↩
開封したときの利便性を考慮して詰めると遅々として進まないので，とにかく詰めることだけに集中すると良いと思いました．↩
立ち合いが必要だと思ったのですが，気付いたらガスが使えなくなっていました．オートロックの物件だったので，少なくとも部屋にいる必要があったとは思うのですが，何の連絡もなかったです．↩
伝達のミスで，冷蔵庫の大きさと自転車の有無に齟齬がありましたが，特に問題はありませんでした．↩
元々引越し当日に行う予定ではありませんでした．引越しの時間が早まったので当日に時間を変更しました．Web での手続きによって立ち合いの時間を変更できることになっているのですが，直近での変更だったので電話で行いました．↩

GW 期間中ということで時間が取れたので，最近勉強していた選択公理についての記事を Mathlog に書いてみました．記事は以下になります：

mathlog.info

今後書くときのために，今回ハマったポイントについてメモしておきたいと思います．

\labelおよび\refコマンドは使用可能．ただし，自動プレビューにおいて\labelコマンドを使用すると数式の反映がされないようなので注意．
公理に対応する囲みボックスは存在しない．今回は定義用のボックスで代替した．
参考文献は「参考文献」のタブから追加することが可能で，ナンバリングが自動で行われる．
よく使用するコマンドは「マクロ」のタブからマクロを追加することが可能．例えば，\left\{ - \right\} なんかは \newcommand{\lr}[1]{\left\{#1\right\}} のようにして追加することができる．使うときは，\lr{A_n} のように使うことができる．[] 内の数字がコマンドに対する引数の個数を表し，引数は #1 のように指定する．

また，以下はLaTeXの問題だが，気になったので書いておきます．

集合の内包的記法の書き方をセミコロンで記述したが，セミコロン前後の空白をどのように制御すべきか不明だった（今回は\;をセミコロンの前後に挟んだ）．
$X_\lambda \in \Lambda \quad (\forall \lambda \in \Lambda)$ などの数式において，カッコの位置をどのように制御すべきか不明だった（今回は\quadコマンドを使用した）．

はじめに
問題設定
勾配降下法
近接勾配法
近接作用素
おわりに
Refference
更新ログ

この記事は「数理最適化 Advent Calendar 2020」の24日目の記事です．

23 日目は @Atsushi_twi さんによる数理最適化初心者のための(線形)割当問題の概要とscipy.optimize.linear_sum_assignmentによる解法です．
25 日目は...空いております！余力のある方は是非に！

Advent Calendar 初参戦です．今回は「数理最適化」がテーマということで，凸計画問題に対する近接勾配法について書こうと思います*1．連続最適化や機械学習に馴染みのある方にとっては新鮮味のない話題かもしれませんが，どうぞよろしくお願いします．

はじめに

数理最適化において，凸集合*2上で凸関数*3を最小化する*4問題は「凸計画問題」あるいは「凸最適化問題」などと呼ばれる．線型計画問題，二次計画問題，半正定値計画問題，二次錐計画問題などは凸計画問題の代表的な例であり，馴染み深いと感じる方も多いと思う．凸計画問題は非線形最適化問題に分類されるが，局所的最適解が大域的最適解に一致するなどの性質上，理論・応用の両面から古くから研究されてきた最も重要な最適化問題のクラスの一つである．

凸計画問題において，制約条件が凸集合に限定されず空間全体*5であるとき，問題は次のように記述される： \begin{align} \min_{x\in\mathbb{R}^n} f(x)\tag{$\diamond$} \end{align} ここに，$f\colon\mathbb{R}^n\to\mathbb{R}$ は凸関数であり，一般に目的関数と呼ばれる．もし関数 $f$ が微分可能であるならば，後述する勾配降下法などのアルゴリズムによって最適解に収束する点列を構成することができる．一方，目的関数が微分できない点を含む場合，関数の勾配の情報が使えないため，そのようなアルゴリズムを適用することができないことは想像に容易いだろう．例えば，機械学習などで頻出するLASSO（$\ell_1$ 正則化）における正則化項や，行列の補完問題（Matrix completion problem）における核ノルム（nuclear / trace norm）は微分できない点を含むことでよく知られている．

なお，本稿は数理最適化がテーマなので，$n$ 次元ベクトルが頻出する．多くの文献ではベクトルを $\mathbf{x}$ あるいは $\boldsymbol{x}$ と表記されることが多い*6が，本稿では太字にせず $x$ と表すこととする．また，本稿を通して，$\|\cdot\|_2$ はユークリッドノルム，$\|\cdot\|_1$ は $\ell_1$ ノルムを表す*7．行列 $X\in\mathbb{R}^{m\times n}$ に対して，フロベニウスノルム $\|X\|_{F}$ および核ノルム $\|X\|_{\mathrm{tr}}$ はそれぞれ \begin{align} \|X\|_{F} := \sqrt{\sum_{i=1}^{\min \{m,n\}} \left\{ \theta_i (X) \right\}^2 },\quad \|X\|_{\mathrm{tr}} := \sum_{i=1}^{\min \{m,n\}} \theta_i (X) \end{align} で定義される．ここで，$\theta_i(X)$ は $X$ の特異値を表す*8．

最後に，指示関数（indicator function）を定義しておく．閉凸集合 $\mathcal{C}$ に対して， \begin{align} 1_\mathcal{C}(x) := \left\{\begin{array}{ll} 0, &\quad \mathrm{if} \; x \in \mathcal{C}, \\ \infty, &\quad \mathrm{otherwise} \end{array}\right. \end{align} と定義される関数 $1_\mathcal{C} \colon \mathbb{R}^n \to \mathbb{R} \cup \{\infty\}$ を集合 $\mathcal{C}$ の指示関数という*9．すなわち，入力ベクトル $x$ が集合 $\mathcal{C}$ に含まれていたら $0$ を，含まれていなかったら $\infty$ を返す関数である．この指示関数を用いると，制約条件を含んだ最適化問題 $\min_{x\in\mathcal{C}} f(x)$ が次のように等価に書き換えられる*10． \begin{align} \min_{x\in\mathbb{R}^n} f(x) + 1_{\mathcal{C}}(x) \end{align}

問題設定

本稿で取り上げる最適化問題は以下のように定式化される問題である： \begin{align} \min_{x \in \mathbb{R}^n} f(x) := g(x) + h(x) \tag{$\ast$} \end{align} ここに，$g\colon\mathbb{R}^n\to\mathbb{R}$ は微分可能な凸関数，$h\colon\mathbb{R}^n\to\mathbb{R}$ は必ずしも微分可能とは限らない凸関数とする．

先に例に挙げたLASSOを例にとると，次のような対応関係になる． \begin{align} \min_{\beta\in\mathbb{R}^{n}}f(\beta) = \underbrace{\frac{1}{2} \| y - X \beta \|_2^2}_{g(\beta)} + \underbrace{\lambda \| \beta \|_1}_{h(\beta)} \end{align}

また，行列の補完問題を例にとると，次のような対応関係になる． \begin{align} \min_{B\in\mathbb{R}^{m\times n}}f(B) = \underbrace{\frac{1}{2} \sum_{(i,j)\in\Omega} \left(Y_{ij} - B_{ij}\right)^2 }_{g(B)} + \underbrace{\lambda \| B \|_{\mathrm{tr}}}_{h(B)} \end{align}

収束解析をする上では，関数 $g$ や $h$ にもう少し条件を付け加える必要があったりするのだが，今回は概説ということなので詳しく触れないこととする．

勾配降下法

問題($\diamond$)の目的関数が微分可能な凸関数である場合，次のような手続きによって，最適解に収束するような点列を生成することができる． \begin{align} x^{k+1} = x^{k} + t_k \cdot d_k, \quad k=0,1,2,\dots \end{align} ここで，$t_k$ は「ステップ幅」と呼ばれる正の実数であり，$d_k$ は降下方向を表すベクトルである．ベクトル $d_k$ として勾配ベクトル $-\nabla f(x^k)$ を採用したアルゴリズムが「最急降下法」に相当する．

少し違った見方をしてみよう．現在の点を $\overline{x}$ として，関数 $f$ の点 $\overline{x}$ の周りで Taylor 展開して得られる関数： \begin{align} f(x) \simeq f(\overline{x}) + \nabla(\overline{x})^\top (x-\overline{x}) + \frac{1}{2} (x-\overline{x})^\top \nabla^2 f(\overline{x}) (x-\overline{x}) \end{align} は点 $\overline{x}$ の周りにおける関数 $f$ の近似した関数になっている．上式においてヘッセ行列 $\nabla^2 f(\overline{x})$ を $tI_n$ とおいた式： \begin{align} \tilde{f}(x) = f(\overline{x}) + \nabla f(\overline{x})^\top (x-\overline{x}) + \frac{1}{2t} \|x-\overline{x}\|_2^2 \end{align} は $x$ に関する二次式であり，最小となる点が存在する．そこで，その最小となる点を求めて次の点として採用することを考えてみる．具体的には，以下のように式変形をすればよい． \begin{align} x^{\star} &= \mathop{\mathrm{arg~min}}_{x} \left\{ \tilde{f}(x) \right\} \\ &= \mathop{\mathrm{arg~min}}_{x} \left\{ f(\overline{x}) + \nabla f(\overline{x})^\top (x-\overline{x}) + \frac{1}{2t} \|x-\overline{x}\|_2^2 \right\} \\ &= \mathop{\mathrm{arg~min}}_{x} \left\{ 2t \cdot \nabla f(\overline{x})^\top x + \|x \|^2_2 - 2\overline{x}^\top x \right\} \\ &= \mathop{\mathrm{arg~min}}_{x} \left\{ \|x \|^2_2 - 2 (\overline{x} - t \cdot \nabla f(\overline{x}))^\top x \right\} \\ &= \mathop{\mathrm{arg~min}}_{x} \left\{ \|x - (\overline{x} - t \cdot \nabla f(\overline{x})) \|^2_2 \right\} \end{align} 以上より，関数 $\tilde{f}$ を最小とする $x$ は以下のようになる： \begin{align} x^{\star} = \overline{x} - t \cdot \nabla f (\overline{x}) \end{align} これは勾配法において勾配ベクトル $d$ を $-\nabla f(x)$ とした最急降下法に相当する．

近接勾配法

問題($\ast$)に対する近接勾配法は以下の更新式によって点列を更新するアルゴリズムである． \begin{align} x_{k+1} &= \textrm{prox}_{t_k,h} ( x_k - t_k \nabla g(x_k) ),\\ \textrm{prox}_{t,h}(y) &= \mathop{\mathrm{arg~min}}_{x} \left\{ h(x) + \frac{1}{2t} \|x-y\|_2^2 \right\}. \end{align} 上式における $\textrm{prox}$ は近接作用素（proximal operator / mapping）と呼ばれる作用素 *11である．意味合いとしては，現在地 $y$ の近くで関数を最小化する問題である．

本稿では上記の更新式がどのような経緯を辿って導出されたものなのか見ていきたい．解くべき問題について，関数 $h$ は微分不可能な点を含むが，関数 $g$ については微分可能なので $g$ に関する勾配の情報のみを用いることを考える．先のアナロジーにより，現在の点 $x^k$ における関数 $f$ の二次式による近似を考えると，次のような手続きを考えたらよさそうである： \begin{align} x^{k+1} &= \mathop{\mathrm{arg~min}}_{x} \left\{ g(x^{k}) + \nabla g(x){^\top} (x^{k}-x) + \frac{1}{2t} \|x^{k}-x\|_2^2 + h(x) \right\} \end{align} この式を変形していくと次のようになる： \begin{align} x^{k+1} &= \mathop{\mathrm{arg~min}}_{x} \left\{ g(x^{k}) + \nabla g(x){^\top} (x-x^{k}) + \frac{1}{2t_k} \|x-x^{k}\|_2^2 + h(x) \right\} \\ &= \mathop{\mathrm{arg~min}}_{x} \left\{ \frac{1}{2t_k} \left( \| x - x^k \|_2^2 + 2t_k \nabla g (x^k)^\top x \right)+ h(x) \right\} \\ &= \mathop{\mathrm{arg~min}}_{x} \left\{ \frac{1}{2t_k} \left\| x - (x^k - t_k \nabla g (x^k)) \right\|_2^2 + h(x) \right\} \end{align} ここで，近接作用素の定義を用いると， \begin{align} x^{k+1} = \mathrm{prox}_{t_k,h} (x^k - t_k \nabla g (x^k)) \end{align} となり，所望の更新式が得られた．

近接作用素

近接勾配法の更新式において，近接作用素 $\mathrm{prox}$ を逐次計算する必要がある．この計算が重たいと効率的に解の更新を行うことができない．本稿では応用上重要な以下の 3 つの関数に対する近接作用素が具体的にどのように構成されるのかについて見ていく．

$\ell_1$ ノルム
核ノルム
指示関数

まず，$\ell_1$ ノルムについて見ていく．関数 $h$ が $h(x) = \lambda\|x\|_1$ で与えられるとき，近接作用素は次のように表される． \begin{align} \mathrm{prox}_{t,h} (x) = \mathop{\mathrm{arg~min}}_{y\in\mathbb{R}^n} \left\{ \frac{1}{2t} \left\| x - y \right\|_2^2 + \lambda \|y\|_1 \right\} \end{align} 簡単な変形により，上記の右辺を最小にするようなベクトル $y^\star$ は解析的に求めることができ，一般にソフト閾値作用素（soft-thresholding operator）と呼ばれる以下の関数と等価となる． \begin{align} \left[ \mathrm{prox}_{t,h} (x) \right]_i = \mathrm{sgn} (x_i) \max \{ |x_i| - \lambda , 0 \} \end{align} ここに，$\mathrm{sgn}(\cdot)$ は符号関数である．

次に，核ノルムについてである．関数 $h$ が $h(X) = \|X\|_{\mathrm{tr}}$ で与えられるとき，近接作用素は次のように表される． \begin{align} \mathrm{prox}_{t,h} (X) = \mathop{\mathrm{arg~min}}_{Y\in\mathbb{R}^{m\times n}} \left\{ \frac{1}{2t} \left\| X - Y \right\|_{F}^2 + \lambda \|Y\|_{\mathrm{tr}} \right\} \end{align} 上式は行列 $X$ の特異値に対するソフト閾値作用素に帰着される．すなわち，行列 $X$ を $X = U\Theta V^\top$ と特異値分解したとき，近接作用素は \begin{align} \mathrm{prox}_{t,h} (X) = U \tilde{\Theta}_{\lambda} V^\top \end{align} で与えられる．ここで，$\tilde{\Theta}_{\lambda}$ は対角成分が $\left(\tilde{\Theta}_{\lambda} \right)_{ii} = \max \{ \theta_i (X) - \lambda , 0 \} , \; (i=1,2,\dots, \min \{m,n\} )$ で与えられる対角行列を表す．

最後に，指示関数について見ていこう．関数 $h$ が $h(x) = 1_{\mathcal{C}} (x)$ で与えられるとき，近接作用素は次のように表される． \begin{align} \mathrm{prox}_{t,h} (X) &= \mathop{\mathrm{arg~min}}_{y\in\mathcal{C}} \left\{ \frac{1}{2t} \left\| x - y \right\|_{2}^2 + 1_{\mathcal{C}}(y) \right\} \\ &= \mathop{\mathrm{arg~min}}_{y\in\mathcal{C}} \left\| x - y \right\|_{2}^2 =: P_{\mathcal{C}}(x) \end{align} 要するに，入力ベクトル $x$ とのユークリッド距離が最小となるような $\mathcal{C}$ 内の点を見つける問題に相当し，一般に凸射影と呼ばれる写像 $P_{\mathcal{C}}$ に他ならない*12．

おわりに

本稿では，凸計画問題に対する近接勾配法に関する初歩的な内容について解説した．途中で断ったように，今回はアルゴリズムの収束解析について触れていない．そのため，目的関数の条件について凸性と微分可能性についてしか言及していない．様々な文献をみていくと，収束解析をするために，関数 $g$ にリプシッツ連続性や平滑性を仮定として置いている場合がある．その上で，リプシッツ定数などの情報を用いてステップ幅などを決定するアルゴリズムを構成し，収束解析がなされているようである．

Refference

本稿を記述する上で参考にした資料を以下に挙げる．より詳しい内容を知りたい場合は以下の文献を参照されたい．

更新ログ

（12/24 8:18）指示関数の記述を整理しました。集合論などで出てくる指示関数と最適化の文脈で出てくる指示関数を混同してしまっていました。今回の場合は、入力ベクトルが集合に含まれていたら $0$ を返し、含まれていなかったら $\infty$ を返す関数の方が正しいです。指摘していただいた @weary_samurai さん、 @paruma184 さん、@jirohhana2 さん、ありがとうございました。

*1:12月に入ってから勉強した内容になります．間違い等ございましたら Twitter:@mirucaaura まで知らせて頂けると幸いです．

*2:集合 $\mathcal{C}$ が任意の $x,y\in\mathcal{C}$ と $\lambda \in [0,1]$ に対して $\lambda x + (1-\lambda)y \in \mathcal{C}$ が成り立つとき，集合 $\mathcal{C}$ は凸集合と呼ばれる．

*3:関数 $f\colon\mathbb{R}^n\to\mathbb{R}\cup\{\infty\}$ が任意の $x,y\in\mathbb{R}^n$ と $\lambda \in [0,1]$ に対して $f(\lambda x + (1-\lambda)y) \leqq \lambda f(x) + (1-\lambda)f(y)$ を満たすとき，関数 $f$ は凸関数と呼ばれる．

*4:もちろん最大化する問題を考えても差し支えない．

*5:本稿の議論はそのままヒルベルト空間に拡張できる点が多いが，応用上分かりやすいようにユークリッド空間で考えることとする．

*6:ちなみに，$\boldsymbol{x}$ ははてなブログにおいて「\boldsymbol{x}」と記述すると出力される．

*7:過剰な配慮かもしれないが，$n$ 次元ベクトル $x = (x_1,\dots, x_n)^\top$ に対して，ユークリッドノルムは $\|x\|_2 = \sqrt{x_1^2 + \dots + x_n^2}$ で定義され，$\ell_1$ ノルムは $\|x\|_1 = |x_1| + \dots + |x_n|$ で定義される．

*8:はてなブログにおいて sigma が表示されないので代わりに $\theta$ を用いている．

*9:閉凸集合に限定するような書き方をしたが，一般の集合に対して定義される概念である．また，空でない閉凸集合に対する支持関数は閉真凸関数になることが示せる．

*10:一見すると制約条件がないような問題に見えるかも知れないが，$x$ が集合 $\mathcal{C}$ に制限されているので，制約を含んだ最適化問題である．

*11:作用素は写像と思っていただいて構わない．

*12:凸射影が効率的に計算できるかは凸集合の形状に依存する．

冷めたコーヒー

Weniger, aber besser

関数の凸性とヘッセ行列の半正定値性の関係（定理 2.30）

参考文献

引越しログ

Summary

時系列

所感

久しぶりに Mathlog に記事を書きました

応用情報技術者試験の受験記

試験勉強

所感

近接勾配法概説