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

一流の大人（ビジネスマン、政治家、リーダー…）として知っておきたい、教養・社会動向を意外なところから取り上げ学ぶことで“気付く力”を伸ばすブログです。データ分析・語学に力点を置いています。　→現在、コンサルタントの雛になるべく、少しずつ勉強中です(※2024年1月21日改訂)。

MENU

Rによるデータサイエンス（16/21）

R

　Rについて

Rによるデータサイエンス(第2版):データ解析の基礎から最新手法まで

Rによるデータサイエンス(第2版):データ解析の基礎から最新手法まで

Amazon

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

前回

power-of-awareness.com

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)) # 誤判別率

次回

power-of-awareness.com

補足　スペック情報

エディション 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
プライバシーポリシー お問い合わせ