回帰で扱う代理損失

入力を$x \in \mathcal{X}$，出力を$y \in \mathbb{R}$とする．分類問題と違って，$m$の定義の仕方に工夫が必要．例えば以下のように$m$を定義する．

$\begin{aligned} m = \epsilon - |y - f(x)| \end{aligned}$

ただし，$\epsilon \geq 0$．そうすると，

$\begin{aligned} \mathbb{I}\left[m < 0\right] \end{aligned}$

で不正解となる．分類問題と同じ代理損失が使える．

$\begin{aligned} \max\{0, 1 - \epsilon + |y - f(x)|\} \end{aligned}$

これは，$\epsilon$-insensitive lossと呼ばれる．

Written with StackEdit.

分類問題で扱う代理損失

入力を$x \in \mathcal{X}$，出力を$y \in \{-1, 1\}$とする．$f(x) \geq 0$のときは正例と，$f(x) < 0$のときは負例と予測する分類器を考える．$m = yf(x)$とすると，$m \geq 0$のとき正解で，$m < 0$のとき，不正解となる．不正解のときは$1$を，正解の時は$0$になるような損失関数を以下のように定義する．

$\begin{aligned} \mathbb{I}\left[m < 0\right] \end{aligned}$

ただし，$\mathbb{I}$は指示関数である．有名な損失関数として，hinge lossや，exponential lossがある．

$\begin{aligned} \ell_{hinge} &= \max\{0, 1 - m\}\\ \ell_{exp} &= \exp \left(- m\right) \end{aligned}$

流れ的には，目的関数を指示関数を使って上のように定義して，代理損失で置き換えるって感じかな．

Written with StackEdit.

ヒンジ損失の劣勾配

ヒンジ損失の劣勾配を求めようと思う．と言っても，[1]のPointwise maximumのテクニックを使えば一発なのだが．まずはこれを示そう．以下の関数を考える．

$\begin{aligned} f(x) = \max_{i=1,\ldots,m} f_i(x) \end{aligned}$

ある$x$に対して，$f(x) = f_k(x)$だったとする．$g \in \partial f_k(x)$とする．

$\begin{aligned} f(y) \geq f_k(y) \geq f_k(x) + g^T (y - x) = f(x) + g^T (y - x) \end{aligned}$

より，$g \in \partial f(x)$である．
以下で定義されるヒンジ損失を考える．

$\begin{aligned} \ell_{hinge}(\bm{w}) = \max\{0, 1 - y \bm{w}^T \bm{x}\} \end{aligned}$

上記の事実を踏まえると，

$\begin{aligned} \bm{g} = \left\{\begin{array}{ll} \bm{0}_d & (1 - y \bm{w}^T \bm{x} \geq 0) \\ -y \bm{x} & (1 - y \bm{w}^T \bm{x} < 0) \end{array}\right. \end{aligned}$

となる．

1. https://see.stanford.edu/materials/lsocoee364b/01-subgradients_notes.pdf

Written with StackEdit.

L1-Normの劣微分

前にL1-Normの劣微分
について書いたが，間違っていたので修正．

[1]のPointwise maximumの項に書いてあるテクニックを使う．$\bm{x} \in \mathbb{R}^d$として，$\|\bm{x}\|_1$を以下のように表す（結構感動的）．

$\begin{aligned} f(\bm{x}) = \|\bm{x}\|_1 = \max\{\bm{s}^T \bm{x} \vert \bm{s} \in \{-1, 1\}^d\} \end{aligned}$

$\bm{x}$に対して，$\|\bm{x}\|_1 = \bm{s}^T \bm{x}$となったとする．この時，$\bm{s}$は以下のようになる．

$\begin{aligned} s_j = \left\{\begin{array}{ll} 1 & (x_j > 0) \\ -1 & (x_j < 0) \\ -1 \ \text{or} \ 1 & (x_j = 0) \end{array}\right. \end{aligned}$

$\bm{s}^T \bm{x}$のsubgradientを$\bm{g}$とすると，$\bm{g} \in \partial f(\bm{x})$である．最後に，$\bm{g} = \bm{s}$である．

参考:

1. https://see.stanford.edu/materials/lsocoee364b/01-subgradients_notes.pdf

Written with StackEdit.

劣勾配法（Lasso）

今回は劣勾配法を使って，L1正則化最小二乗学習（Lasso）を実装してみたいと思う．最適化問題を，以下のように定義する．

$\begin{aligned} \min_{\bm{w}} \ \frac{1}{2} \sum_{i=1}^n \left(y_i - \bm{w}^T \bm{x}_i\right)^2 + \lambda \|\bm{w}\|_1 \end{aligned}$

第一項は$w_j$で微分可能なので，その勾配を使い，第二項は列勾配を使う．これで更新式を出すと，

$\begin{aligned} w_j \leftarrow w_j - \eta \left(- \sum_{i=1}^n (y_i - \bm{w}^T \bm{x}_i) x_{ij} + \lambda sgn(w_j)\right) \end{aligned}$

となる．まとめると，

$\begin{aligned} \bm{w} \leftarrow \bm{w} - \eta \left(- \bm{X}^T \left(\bm{y} - \bm{X} \bm{w}\right) + \lambda sgn(\bm{w})\right) \end{aligned}$

となる．これを実装したクラスが以下．

class Lasso(object):

    def __init__(self, eta=1e-4, reg=1, max_iter=1000):
        self.eta = eta
        self.reg = reg
        self.max_iter = max_iter

    def fit(self, X, y):
        w = np.zeros(X.shape[1])
        for _ in range(self.max_iter):
            g = - X.T.dot(y - X.dot(w)) + self.reg * np.sign(w)
            w = w - self.eta * g
        self.w = w
        return self

    def predict(self, X):
        return X.dot(self.w)

1-d example

1次元で実験してみる．

import numpy as np
import matplotlib.pyplot as plt

random_state = 0
rnd = np.random.RandomState(random_state)

n, d = 100, 2
sigma = 0.5
x = np.linspace(-2, 2, n)
y = x + rnd.normal(0, sigma, size=n)

model = Lasso()
model.fit(np.c_[x], y)

y_pred = model.predict(np.c_[x])
plt.scatter(x, y)
plt.plot(x, y_pred, "r-", linewidth=3)

結果は以下．きちんと学習されているっぽい．

100-d example

次にL1正則化を施したことにより，解がスパースになることを確認するために，100次元のデータで実験する．コードは例えば以下．

from sklearn.datasets import make_regression
from sklearn.model_selection import train_test_split
from sklearn.metrics import mean_squared_error

X, y = make_regression(n_samples=500, n_features=100, n_informative=10, random_state=random_state)

X_train, X_test, y_train, y_test = train_test_split(X, y, random_state=random_state)

model = Lasso()
model.fit(X_train, y_train)

y_pred = model.predict(X_test)

print(np.sqrt(mean_squared_error(y_test, y_pred))) # 0.008980784567524594

きちんと解がスパースになっていることを確認してみる．例えば，

np.around(model.w, 3)

とすると，

array([ 0.   , -0.   , -0.   ,  0.   , -0.   ,  0.   ,  0.   , -0.   ,
       -0.   ,  0.   ,  0.   , -0.   ,  0.   ,  0.   , -0.   , -0.   ,
        0.   ,  0.   , -0.   ,  0.   ,  0.   , -0.   , 99.673,  0.   ,
        0.   , -0.   ,  0.   ,  0.   , -0.   ,  0.   , 41.788,  0.   ,
       -0.   , -0.   ,  0.   , -0.   ,  0.   , -0.   ,  0.   , -0.   ,
       -0.   , 74.23 , -0.   ,  0.   ,  0.   ,  0.   ,  0.   ,  0.   ,
       -0.   , -0.   ,  0.   , -0.   ,  0.   ,  0.   ,  0.   ,  0.   ,
       -0.   , -0.   ,  0.   , -0.   ,  0.   , 37.469, -0.   ,  0.   ,
        0.   ,  0.   ,  0.   ,  0.   ,  0.   ,  0.   , -0.   , -0.   ,
       30.374, 44.89 ,  0.   , 50.938,  0.   ,  0.   , -0.   , -0.   ,
       -0.   ,  0.   ,  0.   , 98.669, -0.   ,  0.   , 27.71 ,  0.   ,
        0.   , -0.   , -0.   , -0.   , -0.   , 86.562, -0.   ,  0.   ,
       -0.   , -0.   ,  0.   ,  0.   ])

となって，n_informativeの数だけ非ゼロになっている．

Written with StackEdit.