分布的ロバスト最適化をロジスティク回帰に適用しcvxpyで実装

Last updated at 2023-12-14Posted at 2022-10-17

はじめに

分布的ロバスト最適化をロジスティク回帰に適用した論文がNIPS2015に提案されました。
cvxpyで実装してみたので備忘録として残します。

分布的ロバスト最適化については、みるかさんの記事が参考になります。
最適輸送理論から導入して、分布的ロバスト最適を紹介されているのは、とても流れがわかりやすく素晴らしいなぁと思います。

ロジスティク回帰のモデル

標本の集合$\lbrace(x_1,y_1)\cdots(x_n,y_n)\rbrace$が与えられているとする。標本$i$の属性を$x\in\mathbb{R}^d$、標本$i$のラベルを$y_i\in\lbrace\pm1\rbrace$とする。このとき、$x$が与えられたときの$y$の確率は、

$$
\mathrm{Prob}(y|x)=[1+\exp(-y\langle\beta,x\rangle)]^{-1}.
$$

ここで、$\beta\in\mathbb{R}^d$は重みベクトルである。$\mathrm{Prob}(+1|x)+\mathrm{Prob}(-1|x)=1$となっているので確かに確率になっている。

分布的ロバスト最適化による定式化

ロジスティク回帰の分布的ロバスト最適化による最適化問題は下式で表されれる。

\begin{eqnarray}
\inf_{\lambda,s_i} &  &\lambda\epsilon+\dfrac{1}{N}\sum_{i=1}^N s_i\\
s.t. & &
l_{\beta}(\hat{x}_i,\hat{y}_i)\leq s_i \ \ (\forall i \leq N)\\
     & &
l_{\beta}(\hat{x}_i,-\hat{y}_i)-\lambda\kappa\leq s_i \ \ (\forall i \leq N)\\
     & &
\|\beta\|_{*}\leq\lambda     
\end{eqnarray}

ただし、$l_{\beta}(x,y)=\log(1+\exp(-y,\langle\beta ,x\rangle))$である。

cvxpyによる実装

import numpy as np
import cvxpy as cp
from sklearn.metrics import confusion_matrix

class dr_logistic_regression():
  def __init__(self,epsilon,kappa,pnorm):
    self.epsilon = epsilon
    self.kappa   = kappa
    self.pnorm   = pnorm

  def fit(self,x,y):
    self.N = x.shape[0]
    self.M = x.shape[1]
    self.x = x
    self.y = y
    beta_   = cp.Variable((self.M,1))
    lambda_ = cp.Variable((1,1))
    s       = cp.Variable((self.N,1))

    constraints = [
                 cp.logistic( cp.multiply(-self.y,self.x@beta_) ) <= s 
                ,cp.logistic( cp.multiply(+self.y,self.x@beta_) ) - lambda_*self.kappa <= s
                ,cp.norm(beta_, self.pnorm)<=lambda_
                  ]

    objective = lambda_*self.epsilon + 1/self.N*cp.sum(s)

    problem = cp.Problem(cp.Minimize(objective),constraints)
    problem.solve(verbose=False)
    self.beta_ = beta_.value
    return problem.value,beta_.value,lambda_.value,s.value
    
  def infer(self,x):
      y_est = []
      for i in 1/(1+np.exp(-x@self.beta_)):
        if i>=0.5:
          y_est.append(1)
        if i<0.5:
          y_est.append(-1)
      self.y_est = y_est
      return y_est
    
  @staticmethod
  def confusion_matrix(y,y_est):
    C = confusion_matrix(y, y_est)
    return C

コードは下記に置きました。よかったら使ってみてください。

参照

You get articles that match your needs
You can efficiently read back useful information
You can use dark theme

What you can do with signing up