「大人の教養・知識・気付き」を伸ばすブログ

※今月(8月)は一部コンテンツを隔週更新にします(夏休みです…)。 一流の大人(ビジネスマン、政治家、リーダー…)として知っておきたい、教養・社会動向を意外なところから取り上げ学ぶことで“気付く力”を伸ばすブログです。

MENU

Rによるデータサイエンス(16/XX)

 Rについて

をベースに学んでいく。
 今回はカーネル法サポートベクターマシン(PP.206-215)を扱う。

17. カーネル法サポートベクターマシン

 非線形的なデータ構造を線形構造に変換できれば、線型的なデータ解析手法で非線形データを扱うことが可能である。データを異なる空間に射影することで、非線形構造を線形構造に変換することが可能な場合がある。
 カーネル法とはカーネル関数を用いてデータを表現し直す方法で、\Phi(\boldsymbol{x})によりデータ\boldsymbol{x}を異なる空間へ写像し、データの中に含まれている情報を表現する。

17.1 カーネル主成分分析

 カーネル主成分分析は非線形主成分分析とも呼ばれ、大まかには以下の手順で計算する:

   (1) カーネル関数K(\boldsymbol{x},\boldsymbol{z})を決める。
   (2) データからの写像行列K_{n\times n}を求める。
   (3) K_{n\times n}固有値固有ベクトルを求める。
   (4) 固有値固有ベクトルを正規化する。
##########################
### カーネル主成分分析 ###
##########################

library("kernlab")

op <- par(mfrow = c(2,2), mai = c(rep(0,3,4)))

x <- as.matrix(iris[,1:4])

#
iris_kpc1 <- kpca(x, kernel = "rbfdot", kpar = list(sigma = 0.1), features = 2)
plot(rotated(iris_kpc1), col = as.integer(iris[,5]))

#
iris_kpc3 <- kpca(x, kernel = "polydot", kpar = list(degree = 1), features = 2)
plot(rotated(iris_kpc3), col = as.integer(iris[,5]))

#
iris_kpc5 <- kpca(x, kernel = "polydot", kpar = list(degree = 5), features = 2)
plot(rotated(iris_kpc5), col = as.integer(iris[,5]))

#
iris_kpc6 <- kpca(x, kernel = "polydot", kpar = list(degree = 6), features = 2)
plot(rotated(iris_kpc6), col = as.integer(iris[,5]))

par(op)

17.2 サポートベクターマシン

 サポートベクターマシンは分類と回帰問題を主としたデータ解析手法であり、線形分離が可能な高次元の仮説空間において線形的なアプローチで学習を行う手法である。
 学習データ(\boldsymbol{x}_1,y_1),\cdots,(\boldsymbol{x}_n,y_n)があるとすると、\boldsymbol{x}_k={}^{t}(x_{k1},x_{k2},\cdots,x_{kp})は個体の特徴ベクトルで、y_kは目的変数であり、回帰問題では数値、分類問題ではクラスのラベルを与える。線形回帰と線形判別の問題では、



\begin{aligned}
y=f(\boldsymbol{x})=\displaystyle{\sum_{i=1}^{p}\omega_ix_i }+b={}^{t}\boldsymbol{\omega}\boldsymbol{x}+b
\end{aligned}


で定式化できる。ここで\boldsymbol{\omega}={}^{t}(\omega_1,\cdots,\omega_p)は各説明変数に与えるウェイトである。



 判別方法として


\begin{aligned}
y=\begin{cases}
1&,\ {}^{t}\boldsymbol{\omega}\boldsymbol{x}+b\geq1,\\-1&,\ {}^{t}\boldsymbol{\omega}\boldsymbol{x}+b\leq-1,\\
\end{cases}
\end{aligned}


を考える。
 平面{}^{t}\boldsymbol{\omega}\boldsymbol{x}+b=1,{}^{t}\boldsymbol{\omega}\boldsymbol{x}+b=-1の間隔であるマージンを最大にすることとする。これはすなわち\|\boldsymbol{\omega}\|を最小にすることに等しい。そのような解を求める第1ステップは、\mathrm{Lagrange}乗数\boldsymbol{\alpha}=(\alpha_1,\cdots,\alpha_n)の関数を導入し、\boldsymbol{\omega},bを推測する。すなわち



\begin{aligned}
L(\boldsymbol{\omega},\boldsymbol{\alpha},b)=\displaystyle{\frac{1}{2}\|\boldsymbol{\omega}\|}+\displaystyle{\sum_{i=1}^{n}\alpha_i\{1-y_i({}^{t}\boldsymbol{\omega}\boldsymbol{x}+b)\} }
\end{aligned}


を考える。
 まず本式の両辺を\boldsymbol{\omega},b偏微分した方程式



\begin{aligned}
\displaystyle{\frac{\partial L}{\partial \boldsymbol{\omega}}}&=\boldsymbol{\omega}-\displaystyle{\sum_{i=1}^{n}\alpha_iy_i\boldsymbol{x}_i}=0,\\
\displaystyle{\frac{\partial L}{\partial b}}&=-\displaystyle{\sum_{i=1}^{n}\alpha_iy_i}=0
\end{aligned}


を解くことで



\begin{aligned}
\boldsymbol{\omega}&=\displaystyle{\sum_{i=1}^{n}\alpha_iy_i\boldsymbol{x}_i},\\
\displaystyle{\sum_{i=1}^{n}\alpha_iy_i}&=0
\end{aligned}


が得られる。これらを代入することで



\begin{aligned}
L(\boldsymbol{\alpha})=\displaystyle{\sum_{i=1}^{n}\alpha_i}-\displaystyle{\frac{1}{2}\sum_{i=1}^{n}\sum_{j=1}^{n}\alpha_i\alpha_jy_iy_j{}^{t}\boldsymbol{x}_i\boldsymbol{x}_j }
\end{aligned}


を得る。
 次にL(\boldsymbol{\alpha})を最大にするような\boldsymbol{\alpha}=\hat{\boldsymbol{\alpha}}を求める、すなわち


\begin{aligned}
\max\left\{L(\boldsymbol{\alpha})\right\}=\max\left\{\displaystyle{\sum_{i=1}^{n}\alpha_i}-\displaystyle{\frac{1}{2}\sum_{i=1}^{n}\sum_{j=1}^{n}\alpha_i\alpha_jy_iy_j{}^{t}\boldsymbol{x}_i\boldsymbol{x}_j }\right\}
\end{aligned}

を、制約条件\alpha_i\geq0, \displaystyle{\sum_{i=1}^{n}\alpha_iy_i}=0の下で解く。これにより



\begin{aligned}
\hat{\boldsymbol{\omega}}&=\displaystyle{\sum_{i=1}^{n}\hat{\alpha}_iy_i\boldsymbol{x}_i},\\
\hat{b}&=\displaystyle{\frac{1}{2}\left(\max_{y=-1}\left\{{}^{t}\boldsymbol{\omega}\boldsymbol{x}\right\}+\min_{y=1}\left\{{}^{t}\boldsymbol{\omega}\boldsymbol{x}\right\}\right)},\\
\hat{f}(\boldsymbol{x})&=\displaystyle{\sum_{i=1}^{n}\hat{\alpha}_iy_i\boldsymbol{x}_i}+\hat{b}
\end{aligned}


を得る。多クラス分類でも同様の考え方で導出できる。

###########
### SVM ###
###########

library("kernlab")

data(spam)
dim(spam)
table(spam[,58])

### ランダムに学習データとテストデータを生成する
set.seed[50]

tr_num <- sample(4601, 2500)

spam_train <- spam[tr_num,]
spam_test <- spam[-tr_num,]

# SVMの構築
(spm_svm <- ksvm(type~., data = spam_train, cross = 3))

# 
spm_pre <- predict(spm_svm, spam_test[,-58])
spm_tbl <- table(spam_test[,58],spm_pre)

(1-sum(diag(spm_tbl))/sum(spm_tbl)) # 誤判別率

補足 スペック情報

エディション Windows 10 Home
バージョン 20H2
プロセッサ Intel(R) Core(TM) i5-1035G4 CPU @ 1.10GHz 1.50 GHz
実装 RAM 8.00 GB
システムの種類 64 ビット オペレーティング システム、x64 ベース プロセッサ
R バージョン 4.1.3 (March, 2022)
RStudio バージョン 2022.02.2 Build 485
プライバシーポリシー お問い合わせ