読者です 読者をやめる 読者になる 読者になる

東京に棲む日々

データ分析、統計、ITを勉強中。未だ世に出ず。

NN基礎からCNNと画像解析入門 - Chainer

大分前にちょこっと画像解析と畳込みニューラルネット(CNN)を勉強してて、周りに話そうかなと思ってスライド作ったけど、忙しくてそれどころじゃなくなったので、とりあえず保存。

Chainerのコードとかそのうちまとめて保存する予定...。

 

 

p.2 -11までがCNNに進むために理解が必要なニューラルネットワーク(NN)基礎。

p.12-23が画像データとCNNとChainerでの実装例。

 

 

 

参考書は、おなじみの以下書籍とマニュアル。

パターン認識と機械学習 上

パターン認識と機械学習 上

 

 

深層学習 Deep Learning (監修:人工知能学会)

深層学習 Deep Learning (監修:人工知能学会)

 

 

Chainerマニュアル

 

 

A Small Study on Segmentation by Clustering Analysis

業務で、Webログや契約者データを用いて、クラスター分析を使ったセグメンテーション(顧客分類)に関して考える必要があったので、それに関するメモ。
 

1)セグメンテーションに関して

以下の文献で、参考になることが多かったので、引用とメモ。

中村 2008 “マーケット・セグメンテーション - 購買履歴データを用いた販売機会の発見”

マーケット・セグメンテーション―購買履歴データを用いた販売機会の発見 (専修大学商学研究所叢書)

マーケット・セグメンテーション―購買履歴データを用いた販売機会の発見 (専修大学商学研究所叢書)

 

 

セグメンテーションアプローチ(顧客分類し、グループ別に異なったコミュニケーション)の際の注意点。使用する変数に関して。

 

 
スライドp.2
セグメンテーションアプローチにおいて、これらの評価指標があいまいな場合、実行不可能、もしくは結果に繋がらないと考えられる。

スライドp.3
実際に分析する場合、データ取得が容易だからと言って、性別や年齢といった変数を用いても効果は低い。
コンバージョンに直接繋がる影響の大きい変数は、一般的に抽象的かつ取得困難。


2)クラスター分析に関しての実験
「何でこのように分けたのだ」、「どのような変数がビジネス上重要なのか」といった質問に、クラスター分析は非常に弱い。クラスター分析は教師なし学習。データの一つの見方を提供する分析であるり、使用する変数によって安定性を失ってしまう分析であり、使い方を気をつけないと、間違った結論を導いてしまう可能性があるので注意が必要だと思う。

スライドp.3-7では、あるビジネス目標を達成するためのセグメンテーションにクラスター分析を用いる場合、なんでも変数を投入してクラスター分析してしまっていいのか、に関して実験してみた。
  

 

結論としては、セグメンテーションアプローチをするといった何らかのビジネス目標がある場合、データの取得が容易な変数だけを集めてクラスター分析にかけると言ったアプローチは危険で、クラスター分析に投入する変数は吟味すべきであると考える。

スライドp.8では、簡単であるが、クラスター分析によるセグメンテーションアプローチの提案。

 

Stackingに関して(2)

前回の続き。
highschoolstudent.hatenablog.com


サンプルデータの作成。あまり細かいことは考えず適当なパラメータ設定で適当に作成。

library(MASS)		# 乱数の作成(多変量正規分布)
library(glmnet)		# Ridge回帰
library(lme4)		# 混合モデル
library(lattice)	# lme4による分析結果の可視化
n_obs <- 2000	# サンプル数
n_r <- 20		# rの水準数
n_x <- 30		# Xの列数

n_r_each <- n_obs/n_r	# 各rの水準のサンプル数


###  r(random effect)の作成  ###

s <- 5	# random effectの標準偏差

# rのmeanの生成(20水準)
set.seed(123)
mean_r <- floor(rnorm(n=n_r, mean=0, sd=15))
mean_r

# mean_rのそれぞれの要素を20回繰り返す
mean_r_vec <- rep(mean_r, times=1, each=n_r_each)	# n_obs行のベクトル

# 乱数(sd=sの正規分布)を加える
set.seed(123)
r_vec <- mean_r_vec + rnorm(n=n_obs, mean=0, sd=s)

# 各rのID列の作成
r_id=rep(1:n_r, times=1, each=n_r_each)


###  X(2000×30)の作成  ###

mu <- rep(0, n_x)		# 各変数の平均は0
sig <- diag(n_x)		# 相関係数行列(30×30)
sig[sig==0] <- 0.1		# 30の変数間に0.1の相関を与える
sig_var <- sig*4		# 分散共分散行列。分散を4とする
set.seed(123)
x_mtrx <- mvrnorm(n_obs, mu=mu, Sigma=sig_var)	# 1000*30の行列

X <- data.frame(x_mtrx)


###  目的変数yの作成  ###

# すべての係数(β)は1とする。重み付け無しで30列を足し合わせる
logit_p <- apply(X, MARGIN=1, sum)

# randon effectを加える
logit_r_p <- logit_p + r_vec

# scaleしておく
p <- 1 / (1+exp(-1*scale(logit_r_p)))		

# yの生成
set.seed(123)
y <- rbinom(n=n_obs, size=1, prob=p)


### data.frameにまとめる ###

data1 <- data.frame(r_id, mean_r_vec, X, p, y)
# 列IDを付けておく
data1 <- data.frame(OBS_ID=1:nrow(data1), data1)
# データ型の変更
data1[, "r_id"] <- as.factor(data1[, "r_id"])
data1[, "y"] <- as.factor(data1[, "y"])
data1[, "OBS_ID"] <- as.character(data1[, "OBS_ID"])
> str(data1)
'data.frame':   2000 obs. of  35 variables:
 $ OBS_ID    : chr  "1" "2" "3" "4" ...
 $ r_id      : Factor w/ 20 levels "1","2","3","4",..: 1 1 1 1 1 1 1 1 1 1 ...
 $ mean_r_vec: num  -9 -9 -9 -9 -9 -9 -9 -9 -9 -9 ...
 $ X1        : num  -0.267 2.431 -0.354 -1.632 -2.469 ...
 $ X2        : num  0.0797 -0.6252 -3.987 -0.219 0.9425 ...
 $ X3        : num  -2.678 1.573 0.518 -3.111 0.12 ...
 $ X4        : num  3.29 -0.87 -7.05 1.1 -1.74 ...
 $ X5        : num  3.235 -0.832 2.129 -1.876 2.369 ...
 $ X6        : num  -1.846 -0.7481 0.0862 0.16 -2.193 ...
 $ X7        : num  2.383 -2.76 2.649 0.34 -0.901 ...
 $ X8        : num  -1.725 -1.1587 -1.2911 0.325 0.0011 ...
 $ X9        : num  0.0625 1.6856 0.11 -0.6957 3.046 ...
 $ X10       : num  -0.0845 2.1115 -2.3013 3.282 1.1797 ...
 $ X11       : num  0.6461 1.1761 -3.7507 1.4679 0.0739 ...
 $ X12       : num  1.024 2.317 -1.685 -0.367 -2.538 ...
 $ X13       : num  -1.27 -3.483 -0.253 0.399 2.68 ...
 $ X14       : num  -2.569 -0.559 -3.942 -3.517 -0.387 ...
 $ X15       : num  2.5588 -2.0188 0.0141 0.8221 0.2959 ...
 $ X16       : num  2.227 2.08 -0.565 0.555 0.261 ...
 $ X17       : num  -3.413 1.206 1.578 0.226 -1.649 ...
 $ X18       : num  3.189 -0.286 -1.177 1.409 -1.41 ...
 $ X19       : num  4.053 0.355 0.173 -2.395 0.607 ...
 $ X20       : num  -0.356 -0.632 -1.924 0.464 1.393 ...
 $ X21       : num  2.0457 -0.4395 0.7572 -2.1578 -0.0979 ...
 $ X22       : num  2.877 1.097 -1.977 1.85 0.151 ...
 $ X23       : num  -0.447 -0.177 -1.276 -1.001 -0.904 ...
 $ X24       : num  -0.506 0.924 -0.718 0.292 1.549 ...
 $ X25       : num  2.343 0.305 -3.97 -3.443 -0.222 ...
 $ X26       : num  -1.731 -0.769 2.971 2.385 -0.441 ...
 $ X27       : num  1.6675 -0.0147 -2.1689 1.4975 -2.7424 ...
 $ X28       : num  0.0978 3.6201 -2.9536 1.036 1.4886 ...
 $ X29       : num  -2.901 0.138 -3.443 -1.862 -1.06 ...
 $ X30       : num  0.1366 -0.6662 0.0819 3.1395 -0.1976 ...
 $ p         : num  0.491 0.429 0.162 0.374 0.364 ...
 $ y         : Factor w/ 2 levels "0","1": 1 2 1 2 2 1 1 1 1 1 ...
plot(data1[,"r_id"],data1[,"p"], xlab="r_id", ylab="p")

f:id:High_School_Student:20151230135232j:plain
rによって、適当に散らばったデータ。


Train(Stacking用に分割), Testセットの作成。

Train:Test=1600:400。Trainは8分割し、各Stacking用データは200とする。

ID列(dataSet_ID)を作成し、管理。
Stacking用に8分割したTrainセットは1,2,3,4,5,6,7,8とIDを付け、Testセットを99とする。

# Train(Stacking), Test IDの作成
# Train -> 1,2,3,4,5,6,7,8
# Test -> 99
each_ID <- c(rep(1,10),rep(2,10),rep(3,10),rep(4,10),rep(5,10),rep(6,10),rep(7,10),rep(8,10),rep(99,20))
dataSet_ID <- rep(each_ID, n_r)

data2 <- data.frame(data1, dataSet_ID)

# Train, Test dataの分割
trainData <- data2[data2[,"dataSet_ID"]!=99, ]
testData <- data2[data2[,"dataSet_ID"]==99, ]


では、Ridge回帰で、1stステージモデルを作成する。
glmnetパッケージを使用。
highschoolstudent.hatenablog.com

# Trainデータの説明変数
exp_vars <- as.matrix(trainData[,c("X1","X2","X3","X4","X5","X6","X7","X8","X9","X10",
					"X11","X12","X13","X14","X15","X16","X17","X18","X19","X20",
					"X21","X22","X23","X24","X25","X26","X27","X28","X29","X30")]
)

# Trainデータの目的変数
#target_var <- trainData[,"y"]
target_var <- as.vector(trainData[,"y"])

# Testデータの説明変数
exp_vars_test <- as.matrix(testData[,c("X1","X2","X3","X4","X5","X6","X7","X8","X9","X10",
					"X11","X12","X13","X14","X15","X16","X17","X18","X19","X20",
					"X21","X22","X23","X24","X25","X26","X27","X28","X29","X30")]
)

まず、Cross Validationで最適なλを見つける。

cv.glmnetは、k+1回数glmnetを走らせる。
k回でλを決定し、残り1回でパラメータを推定する。

前回述べた通り、この分割はStackingの分割と同とする必要はない。
しかし、foldid引数にIDが指定できる、8分割で十分だろうという事で、Cross ValidatioとStackingの分割は同じとして処理する。

set.seed(123)
fit1st_ridgeCV <- cv.glmnet(
	x=exp_vars,
	y=target_var,
	family="binomial",
	#nfold=8,
	foldid=trainData[,"dataSet_ID"],	# CVのIDを指定
	keep=TRUE,	# CVのIDを保存
	alpha=0
)

plot(fit1st_ridgeCV)

# CV(k=8)による、最適なλ
lambdaMin <- fit1st_ridgeCV$lambda.min

f:id:High_School_Student:20151230140019j:plain

> lambdaMin
[1] 0.1540261

見つけたλ(lambdaMin)によって、2ndモデル用の変数(Ridge回帰による予測値)を作成。

Trainデータの処理。
8個に分けた各Stackingセットを一つ抜いて、残りのデータで固定したλ(lambdaMin)でglmnetを実行し、パラメータを推定、抜いておいたセットにこの推定したモデルをあてはめて、予測値を計算する。
これをStackingに分けた数(8回)実行する。

# Stackingセット(= CVセット)のID
cv_id <- as.integer(levels(as.factor(trainData[,"dataSet_ID"])))
# 1 2 3 4 5 6 7 8

# 8回のあてはめ結果オブジェクトの格納用リスト
list_res <- list()

# 予測値格納用のdf(col1:OBS_ID、col2:予測obsのID、col3:予測値)
pred_y <- NULL

for(id in cv_id){
	#print(paste("cv_id:",id))

	# パラメータ推定用データ
	train <- trainData[trainData[,"dataSet_ID"]!=id, ]
	#print("train - cv_id")
	#print(table(train[,"dataSet_ID"]))

	x <- as.matrix(train[,c("X1","X2","X3","X4","X5","X6","X7","X8","X9","X10",
				"X11","X12","X13","X14","X15","X16","X17","X18","X19","X20",
				"X21","X22","X23","X24","X25","X26","X27","X28","X29","X30")])
	y <- as.vector(train[,"y"])
	
	# あてはめ用データ
	target <- trainData[trainData[,"dataSet_ID"]==id, ]
	#print("target - cv_id")
	#print(table(target[,"dataSet_ID"]))
	
	x_target <- as.matrix(target[,c("X1","X2","X3","X4","X5","X6","X7","X8","X9","X10",
					"X11","X12","X13","X14","X15","X16","X17","X18","X19","X20",
					"X21","X22","X23","X24","X25","X26","X27","X28","X29","X30")])

	fitRidge <- glmnet(	# 指定したλ(lambdaMin)でリッジ回帰
		x=x, 
		y=y, 
		family="binomial", 
		lambda=lambdaMin, 
		alpha=0 )
		
	list_res[[id]] <- fitRidge	# あてはめ結果オブジェクトの保存
	
	# あてはめように抜いてあったデータに、パラメータ推定用データで推定したモデルをあてはめて、予測値を推定。あてはめデータも保存。
	pred_y <- rbind(pred_y, data.frame(target[,"OBS_ID"], mean(target[,"dataSet_ID"]), predict(fitRidge, newx=x_target)))
	
	#print("---------------------")
}
names(pred_y) <- c("OBS_ID", "target_id", "pred_input")


# リッジ回帰による予測で作成した列を含めたdata.frame(Trainデータ)の作成
trainData2 <- merge(x=trainData, y=pred_y, by.x="OBS_ID", by.y="OBS_ID", sort=FALSE, all=TRUE)

# 列の絞り込みと順序の入れ替え
trainData2 <- trainData2[,
c("OBS_ID","r_id","mean_r_vec",
"X1","X2","X3","X4","X5","X6","X7","X8","X9","X10",
"X11","X12","X13","X14","X15","X16","X17","X18","X19","X20",
"X21","X22","X23","X24","X25","X26","X27","X28","X29","X30",
"p","y","dataSet_ID","pred_input")]

Trainデータの処理。
Trainデータ全体で、固定したλ(lambdaMin)でglmnetを実行し、パラメータを推定、それにTestデータをあてはめ予測値を作成。

fitRidge_test <- glmnet(	# 指定したλ(lambdaMin)でリッジ回帰
	x=exp_vars,
	y=target_var,
	family="binomial", 
	lambda=lambdaMin, 
	alpha=0
)

# Testデータでの予測値
ypred <- predict(fitRidge_test, newx=exp_vars_test)

# リッジ回帰による予測で作成した列を含めたdata.frame(Testデータ)の作成
testData2 <- data.frame( testData, pred_input=as.vector(ypred) ) 

TrainとTestデータを一つにまとめる。

data3 <- rbind(trainData2, testData2)

# pred_input(線形結合)をpに直す
ypredp <- 1 / (1 + exp(-1*data3[,"pred_input"]))

data3 <- data.frame(data3, pred_input_p=ypredp)

以上、1stステージモデル。

では、2ndステージモデルに入る。

# Train, Test dataの分割
Sec_trainData <- data3[data3[,"dataSet_ID"]!=99, ]
Sec_testData <- data3[data3[,"dataSet_ID"]==99, ]

以前、混合モデルを使用した際は、glmmMRを使ったが、実際はlme4パッケージの方がメジャーなようなので、こちらを使う。
highschoolstudent.hatenablog.com

predict関数もサポートしており、変量効果の各水準に対する予測値(BLUP(Best Linior Unbiased Predictor))を考慮した予測が可能。

詳細に関しては、lme4: Mixed-effects modeling with R 参照。

ランダム切片モデルとする。

# Intercept model
fit_mixed <- glmer(
	y ~
	pred_input +
	(1|r_id),	# 1は切片の意味
	
	data=Sec_trainData,
	family=binomial, 
	verbose=1
)
> summary(fit_mixed)
Generalized linear mixed model fit by maximum likelihood (Laplace Approximation) [
glmerMod]
 Family: binomial  ( logit )
Formula: y ~ pred_input + (1 | r_id)
   Data: Sec_trainData

     AIC      BIC   logLik deviance df.resid 
  1978.0   1994.2   -986.0   1972.0     1597 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-2.8686 -0.8164 -0.3042  0.8195  3.2091 

Random effects:
 Groups Name        Variance Std.Dev.
 r_id   (Intercept) 0.3299   0.5743  
Number of obs: 1600, groups:  r_id, 20

Fixed effects:
             Estimate Std. Error z value Pr(>|z|)    
(Intercept) -0.003755   0.139643  -0.027    0.979    
pred_input   1.168448   0.094078  12.420   <2e-16 ***
---
Signif. codes:  0***0.001**0.01*0.05 ‘.’ 0.1 ‘ ’ 1

Correlation of Fixed Effects:
           (Intr)
pred_input 0.017 


ranef関数に、結果オブジェクトを引数で渡すと、変量効果の各水準に対する予測値(BLUP(Best Linior Unbiased Predictor))を出力できる。

condVar=TRUEとすると、各変量効果の水準の予測値の95%予測区間も出力する。

dotplot(ranef(fit_mixed, condVar=TRUE))

f:id:High_School_Student:20151230140647j:plain

では、2ndステージモデルの予測値を出力する。

# Trainデータの予測値
predTrain_fit_mixed <- 1 / ( 1 + exp(-1*predict(fit_mixed, newdata=Sec_trainData)) )

# Testデータの予測値
predTest_fit_mixed <- 1 / ( 1 + exp(-1*predict(fit_mixed, newdata=Sec_testData)) )

# data.frameにまとめる
Sec_trainData2 <- data.frame(Sec_trainData, pred_y=predTrain_fit_mixed )
Sec_testData2 <- data.frame(Sec_testData, pred_y=predTest_fit_mixed)

data4 <- rbind(Sec_trainData2, Sec_testData2)


上で使用したprediction関数であるが、引数re.form=NULL(デフォルト)としておくと、変量効果の予測値(BLUP)を含めた予測となる。

allow.new.levelsという引数があり、allow.new.levels=TRUE(デフォルトはFALSE)とすると、新しいデータに対する予測の際、モデル推定のデータになかった変量効果の水準が含まれた場合、変量効果の予測を考慮しない値が返される。
デフォルト(FALSE)だと、モデル推定のデータになかった変量効果の水準が含まれた場合、Errorが返るとのこと。

> head(data4)
  OBS_ID r_id mean_r_vec         X1          X2         X3         X4         X5
1      1    1         -9 -0.2670289  0.07965359 -2.6781329  3.2920218  3.2347046
2      2    1         -9  2.4309568 -0.62516429  1.5726213 -0.8702877 -0.8316664
3      3    1         -9 -0.3538252 -3.98699188  0.5181232 -7.0503276  2.1287953
4      4    1         -9 -1.6324691 -0.21901088 -3.1107837  1.1008922 -1.8763152
5      5    1         -9 -2.4687447  0.94246812  0.1196235 -1.7423328  2.3693425
6      6    1         -9  1.2053775 -2.38738177  0.7406409 -2.9720772 -1.2437253

### 省略 ###

         X27         X28        X29         X30         p y dataSet_ID  pred_input
1  1.6675496  0.09780506 -2.9014770  0.13657163 0.4909926 0          1  0.49382763
2 -0.0146797  3.62010783  0.1383445 -0.66618073 0.4288412 1          1  0.05697817
3 -2.1689243 -2.95355951 -3.4428897  0.08190659 0.1619845 0          1 -0.89058592
4  1.4974864  1.03595553 -1.8618701  3.13950267 0.3741005 1          1  0.02697507
5 -2.7424436  1.48861771 -1.0603866 -0.19760381 0.3637205 1          1 -0.26863659
6  0.1702958  1.57345489 -2.8409584 -1.65512597 0.1465225 0          1 -1.05824906
  pred_input_p    pred_y
1    0.6210077 0.5308346
2    0.5142407 0.4044539
3    0.2909889 0.1833040
4    0.5067434 0.3960387
5    0.4332418 0.3170407
6    0.2576442 0.1557724

pred_yが最終予測値である。


なお、パフォーマンスに関してはAUCで確認している。
結果だけだが、Trainデータで0.7404149、Testデータで0.7289979となる。

Stackingを行わなかった場合のモデルとの比較等に関しては省略。


以上。

Stackingに関して (1)

予測モデルのコンペで使われているStackingという手法に関して。

モデルをEnsembleしまくって、マルチステージモデルにして予測精度を高めるときに使われる手法とのこと。

実務家な私も、業務でマルチステージモデル(2ステージ)を組む必要があったので、そのためのメモ。

kaggleなどで発明され、まだ教科書にはなっていない方法論のようで、以下を参照。

KAGGLE ENSEMBLING GUIDE

 

具体的には以下のような課題(データ)があった。

f:id:High_School_Student:20151228140804j:plain

目的変数y(0,1のBinomial)、説明変数rとXでモデルを作りたい。
rは個々のセールスマンを示すようなカテゴリカルデータで、変量効果として扱いたい。
Xは営業日誌のようなテキストデータから単語抽出して、TF-IDF変換した横長なデータ。

 

混合モデル(Logistic GLMM)を作成したいが、横長なXをそのままロジスティックモデルの変数としてそのまま加えるのは不適切、といった理由で、まずXをRidge回帰によって処理をしてから、その予測値を2段階モデルの説明変数として、rが変量効果となる混合モデルを作成することで対処する。

f:id:High_School_Student:20151228140836j:plain

Ridge回帰に関しては、以下を参照。

highschoolstudent.hatenablog.com

混合モデルに関しては、以下を参照。

highschoolstudent.hatenablog.com

 

マルチステージモデルを組む際、ちゃんとしたStackingの手順を踏まないと予測値が上がらない、と、こういったのに詳しい同僚に言われたので、上の参考文献(KAGGLE ENSEMBLING GUIDE)を参考にしつつ、具体的手順を記してみる。

コアとなる考え方だが、最終モデル(例では2nd Model)以前では、ペアとなっている目的変数をモデルパラメータ推定の際使ってはいけない、と言うことだと思う。

注:文献を読んでの私なりの理解なので、間違ってたり、他にもやり方がある可能性有り

 

f:id:High_School_Student:20151228140911j:plain

データをTrain(学習)データとTest(テスト)データに分ける。

Trainデータをさらに、パラメータ推定のためにいくつかに分ける。図では説明を簡単にするために2つに分けたが、もう少し細かく分割しても良いと思う。
なお、この分割は、モデルの複雑さを決定するためのk分割交差検証のためではない。

 

f:id:High_School_Student:20151228140924j:plain

まず、Trainデータ全体でk分割交差検証法(kは適当な値)などを用い、モデルの複雑さを決定する。今回の例では、Ridge回帰のλ(Complexity Parameter)を決定する。

決定したモデルの複雑さ(λ)を用い、X1とy1でモデル(Mod1)を作成(パラメータを推定)。そのモデルをX2にあてはめ、予測値x2を作成。

決定したモデルの複雑さ(λ)を用い、X2とy2でモデル(Mod2)を作成(パラメータを推定)。そのモデルをX1にあてはめ、予測値x1を作成。

決定したモデルの複雑さ(λ)を用い、X1,X2とy1,y2でモデル(Mod')を作成(パラメータを推定)。そのモデルをTestデータのX'にあてはめ、予測値x'を作成。

予測により作成する各変数(x1、x2、x')は、ペアとなる目的変数(x1はy1、x2はy2、x'はy')を使用してない、と言うことになる。

 

f:id:High_School_Student:20151228141021j:plain

Trainデータにおいて、作成した予測値(x1,x2)と、取っておいたその他の変数(変量効果とするr1,r2)を説明変数、Trainデータの目的変数(y1,y2)を使って、最終モデル(Mod'')を作成する。

これをTestデータの説明変数(r',x')にあてはめ、最終的な予測値(y*)を計算する。

 


では次回、シミュレーションによって作成したサンプルデータでやってみる。

 

単なるSQLメモ

 

このようなSampleデータがあったとする。

f:id:High_School_Student:20151020132347j:plain

 

Group別に、時間(Min)が最大となる行を取り出したい。

 

クエリー実行後。

f:id:High_School_Student:20151020132403j:plain

 

SELECT T1.*
FROM Sample T1
INNER JOIN(
     SELECT Group, MAX(Min) AS MAX_Min
     FROM Sample
     GROUP BY Group
) T2 ON T1.Group = T2.Group
        AND T1.Min = T2.MAX_Min

この場合、時間(Min)に重複がないことを仮定している。

 

またこういった、時間が分かれている(Min、Sec)Sampleデータの場合。

f:id:High_School_Student:20151020132541j:plain

 

クエリー実行後。

f:id:High_School_Student:20151020132551j:plain


SELECT T1.*
FROM Sample T1
INNER JOIN(
     SELECT Group, MAX( Min||Sec ) AS MAX_MinSec
     FROM Sample
     GROUP BY Group
) T2 ON T1.Group = T2.Group
        AND T1.Min || T1.Sec = T2.MAX_MinSec

 

以上、単なるSQLメモ。

 

 

モデル評価基準 追加 - リフトチャートに関して2 - R{ROCR}

highschoolstudent.hatenablog.com

highschoolstudent.hatenablog.com

 

上2つの記事に関連した追加。


そもそもリフトチャートと呼ばれるものは定義があいまいなようで、人によって使っているものが細かく違ったりする。

基本的に似たようなものであるのだが、”モデル評価基準 追加 - リフトチャートに関して2 - R{ROCR}”で作成したのとまたちょっと違ったリフトチャートを描いてみる。

 

業務上、このチャートを標準とすることとなった。確かに、シンプルで理解もしやすいのでお勧め。


以下の手順で作成できる。
1.予測値を小さいものから順に並べる。
2.いくつかのグループに分ける。(今回の19行のデータ数ではグループ数を大きくできないが、10とか20グループに分けるやり方が良く使われるらしい。)各グループ内でのデータ数(行数)は等しい。(必ずしも割り切れる訳ではないので、適当な調整を加える。)
3.各グループでの観測値の平均値を計算する。(観測値は0/1なので、各グループでの1の割合が計算される。)
4.横軸をグループ番号、縦軸を3.で計算した観測値の平均値をプロットする。すべての観測値での平均値もプロットしておくとグラフが理解しやすい。


データはこれまでと同じものを使用。
d_roc1 <- read.delim("clipboard", header=TRUE)
d_roc1

 

   pred1 pred2 observed
1   0.95  0.90        1
2   0.90  0.80        1
3   0.80  0.75        1
4   0.60  0.70        1
5   0.50  0.60        1
6   0.85  0.95        0
7   0.75  0.85        0
8   0.70  0.65        0
9   0.65  0.55        0
10  0.55  0.50        0
11  0.45  0.45        0
12  0.40  0.40        0
13  0.35  0.35        0
14  0.30  0.30        0
15  0.25  0.25        0
16  0.20  0.20        0
17  0.15  0.15        0
18  0.10  0.10        0
19  0.05  0.05        0

 

関数を定義する。 

fn_plot_Lift20bins <- function(pred, obs, aveline=TRUE, nbins=20, color="black"){
    df0 <- data.frame(pred=pred, obs=obs)
    df1 <- df0[order(df0$pred),]
    n <- nrow(df1)

    if(nbins*2>n) {stop(">>> Set smaller nbin. nbins*2 <= nrow")}

    meanobs <- mean(df1$obs)
    nrep <- floor(n/nbins)    # 各binに入る行数
    ndup <- n-nbins*nrep    # nbinでnが割り切れない場合の対処

    bin <- c()
    for(i in 1:nbins){
        if(i<=ndup){ bin <- c(bin, rep(i, nrep+1)) }    # nbinでnが割り切れない場合,最初のbinに1つづつ加えていく
        else{ bin <- c(bin, rep(i, nrep)) }
    }

    df2 <- data.frame(df1, bin=bin)
    #print( df2 )

    ratio <- tapply(df2$obs, df2$bin, mean)

    df3 <- data.frame(bin=1:nbins, ratio)
    print( df3 )

    plot( df3$bin, df3$ratio, type="b" , xlab="Bin", ylab="Ratio", col=color)
    abline(h=meanobs, col="yellow")
    abline( v=1:nbins, lty=3 )
}

 

pred: モデル予測値のベクトル

obs: 観測値のベクトル

aveline: 全観測値の平均線をプロットするかどうか

nbins: グループ数(データ数に対し大きいとエラーになる)

color: 線の色

 

pred1に対するリフトチャート。

fn_plot_Lift20bins(pred=d_roc1$pred1, obs=d_roc1$observed, aveline=TRUE, nbins=5)

  bin     ratio
1   1 0.0000000
2   2 0.0000000
3   3 0.5000000
4   4 0.2500000
5   5 0.6666667

f:id:High_School_Student:20150919145100j:plain


5グループに分けたので、各グループは3か4のデータ数となっている。(4*4+3*1=19)

各グループでの観測値の平均のデータフレームを表示し、それのプロットも出力される。
黄色の横線は、全観測値での平均値。

低い方から並べていっているので、基本的に右上がりのグラフとなる。

良いモデルでは、このリフトチャートは単調増加(monotone increasing)でカーブが急なものとなる。

今回の5グループのプロットだと、途中で平均値が落ちてしまっているが、3グループくらいにして、グループ内でのデータ数を確保すると単調増加のカーブになる。

 

定義した関数は重ね描きには対応してないので、今後修正の必要あり。

 

以上。

 

 

Gradient Boosting Treeを使ってみる - R{gbm} (Part.2)

前回に続き、いくつかの変数変換を試してみる。

highschoolstudent.hatenablog.com

 

Partial Dependence Plotで気づいたことを元に、以下の説明変数に対して変換を実行。

"age" - 数値変数だとノイズに敏感に見えるため、カテゴリー化してみる。
"fnlwgt" - 外れ値にキャップをかけてみる。
"capital.gain" - 外れ値にキャップをかけてみる。
"capital.loss" - 元の変数の意味が分からなくて解釈に困るが、カテゴリー化して様子を見てみる。
"hours.per.week" - 外れ値にキャップをかけてみる。
"native.country" - カテゴリカル変数の多水準問題に2つの方法(両方とも数値化する方法)で取り組んでみる。

 

キャップと言うのは、連続変数に上限下限を設けること。
例えば、"身長"という変数があって180cmで上限のキャップをかけた場合、それ以上高い観測値はすべて180cmに変換する。

 

Tree系のモデルは説明変数の外れ値にロバストであり、キャップをかけることによってパフォーマンスが大幅に改善したりするとは考えられない気がするが、それでも念のために外れ値の処理はしておいたほうが良いとのこと(同僚より)なので、実行しておく。

 

カテゴリー化の理由だが、Plotを見た感じで学習データにオーバーフィッティングしている感じがしたので、カテゴリー化することによりスムーズにしようというのが一つで、もう一つが解釈しやすくなる(Plotで見やすくなる)のではないかという期待から。

 

いずれにせよそこまで科学的な根拠があると言うわけではないが、良くやられているので、とりあえずやってみる、といった要素が強い。

 

キャップ関数を定義。変数(データフレームの列)を引数とし、ベクトルを返す。

# Capをかける関数、カラムと最小/最大値を指定
fn_cap <- function(num_col, low=min(num_col), high=max(num_col)){
     num_col[num_col<=low] <- low
     num_col[num_col>=high] <- high

     return(num_col)

}


カテゴリー化の関数を定義。変数(データフレームの列)を引数とし、ベクトルを返す。
分け方だが、上記理由を考慮した完全な主観。

# age(numeric)をカテゴリー化する関数
fn_age_category <- function(age_col){
     age_col[age_col<20] <- "-19"
     age_col[age_col>=20 & age_col<25] <- "20-24"
     age_col[age_col>=25 & age_col<30] <- "25-29"
     age_col[age_col>=30 & age_col<35] <- "30-34"
     age_col[age_col>=35 & age_col<40] <- "35-39"
     age_col[age_col>=40 & age_col<45] <- "40-44"
     age_col[age_col>=45 & age_col<50] <- "45-49"
     age_col[age_col>=50 & age_col<55] <- "50-54"
     age_col[age_col>=55 & age_col<60] <- "55-59"
     age_col[age_col>=60 & age_col<65] <- "60-65"
     age_col[age_col>=65] <- "65-"

     return(age_col)
}

 

# capital.lossをカテゴリー化する関数
fn_capLoss_category <- function(closs_col){
     closs_col[closs_col==0] <- "0"
     closs_col[closs_col>=0 & closs_col<1000] <- "<1000"
     closs_col[closs_col>=1000 & closs_col<1250] <- "<1250"
     closs_col[closs_col>=1250 & closs_col<1500] <- "<1500"
     closs_col[closs_col>=1500 & closs_col<1750] <- "<1750"
     closs_col[closs_col>=1750 & closs_col<2000] <- "<2000"
     closs_col[closs_col>=2000] <- ">=2000"

     return(closs_col)
}

 

データ加工は、adult2_2から出発。

後で実行する多水準カテゴリカル変数の変換のときに必要であるからだが、学習データ(学習+検証)とテストデータを示すフラグ変数を作成しておく。


前回同様、24420までが学習+検証(0とする)、24421から最後の32561までがテスト(1とする)。

test_col <- c(rep(0,24420),rep(1,(nrow(adult2_2)-24420)))
adult2_3 <- data.frame(adult2_2, Test=test_col)

head(adult2_3)

  age workclass fnlwgt    education education.num marital.status
1  34   Private 191856 Some-college            10        Married
2  40   Private 205047      HS-grad             9        Married
3  35   Private 230054      7th-8th             4        Married
4  27   Private 279608      5th-6th             3        Married
5  34   Private 226872    Bachelors            13       Divorced
6  61 State-gov 151459         10th             6  Never-married
         occupation  relationship  race    sex capital.gain capital.loss
1   Exec-managerial          Wife White Female         7298            0
2      Craft-repair Not-in-family White   Male            0            0
3 Machine-op-inspct       Husband White   Male            0            0
4   Farming-fishing       Husband White   Male            0            0
5             Sales Not-in-family White Female            0            0
6     Other-service Not-in-family Black Female            0            0
  hours.per.week native.country income_class income Test
1             40  United-States         >50K      1    0
2             40  United-States         >50K      1    0
3             40         Mexico        <=50K      0    0
4             40         Mexico        <=50K      0    0
5             40  United-States        <=50K      0    0
6             38  United-States        <=50K      0    0

最終列に0/1データの"Test"カラムが作成される。

 

ageのカテゴリー化。
age2 <- fn_age_category(adult2_3[,"age"])


fnlwgtに500000でキャップ。
fnlwgt2 <- fn_cap(adult2_3[,"fnlwgt"], high=500000)


capital.gainに10000でキャップ。
capital.gain2 <- fn_cap(adult2_3[,"capital.gain"], high=10000)


capital.lossのカテゴリー化。
capital.loss2 <- fn_capLoss_category(adult2_3[,"capital.loss"])


hours.per.weekに80でキャップ。
hours.par.week2 <- fn_cap(adult2_3$hours.per.week, high=80)


データフレームにまとめる。

adult2_4 <- data.frame(adult2_3, age2, fnlwgt2, capital.gain2, capital.loss2, hours.par.week2)

 


では、カテゴリカル変数の多水準問題(High Cardinality Problem)への対処を考えてみる。
今回は逆に、数値化することによって、オーバーフィッティングへの対処を試みる。

一つ目の方法が、単純にランキング化する方法。Simple Rankingと呼ぶことにする。

学習データにおいて(これがTest列を作成しておいた理由)水準別に目的変数を集計し、値が高いものから順にランキングを付ける。

"native.country"列は42水準。この42水準別に目的変数である"income"が1となる割合を計算する。高いものから順に並べ、1から42のランクを付ける。

関数群を定義しておく。

 

############################
###### Simple Ranking ######

# 水準別ランキングをdfとして返す
fn_ranking_table <- function(d_f, y_name, x_name, x_new_name){

     # 水準別平均値のテーブル
     tbl <- tapply(d_f[,y_name], d_f[,x_name], mean)

     # ランキング
     rnk <- rank(-tbl, na.last="keep", ties.method="average")

     # dfとして返す
     tbl_df <- data.frame(names(rnk), rnk)
     names(tbl_df) <- c(x_name, x_new_name)
     return( tbl_df )
}

 

# (df, "目的変数(numeric)名", "説明変数(カテゴリー)名", "Rankingによる新しい変数名", "テスト変数名")
# テスト変数 - 0:train&validation 1:test
fn_SimpleRanking2 <- function(d_f, y_name, x_name, x_new_name, test_name=NULL){

     #テストデータ(水準1)を除く
     if(is.null(test_name)){     # テストデータを指定しない場合
          d_f_train <- d_f
     }
     else if(sum(d_f[,test_name]==1)>0) {      # テストデータ用の水準(1)を指定した場合
          d_f_train <- d_f[d_f[,test_name]!=1,]
     }
     else { stop("--- Test data set??? ---") }     # 不明な変数を指定した場合

     #Rankingテーブルの作成
     tbl_df <- fn_ranking_table(d_f_train, y_name, x_name, x_new_name)

     return( merge(d_f, tbl_df, by=x_name, all.x=TRUE) )
}

############################
############################

fn_SimpleRanking2(データフレーム, "目的変数(numeric)名", "説明変数(カテゴリー)名", "Rankingによる新しい変数名", "テスト変数名(0:学習, 1:テスト)")として実行する。

temp_df <- fn_SimpleRanking2(adult2_4, y_name="income", x_name="native.country", x_new_name="native.country2_SR", test_name="Test")

 


では、もう一つの方法を実行。
学習データにおいて、説明変数のある水準を除いたときの目的変数の平均値をその水準の値とする。

例えば、"native.country"変数でJapan水準を除いたときの"Income"変数の平均値をJapanの値とする。

さらに、オーバーフィッティングを押さえるため、適当なノイズを加える。

Mean of Response Variable as Explanetory Variable by Leave One Out、略してLOOと呼ぶことにする。

この方法は、同僚に教えてもらった。

 

###########################
########### LOO ###########

# 説明変数の水準を指定すると、その水準を除いた目的変数の平均値を返す関数
# (df, "目的変数(numeric)名", "説明変数(カテゴリー)名", 説明変数の水準)
fn_mean_oneLevel <- function(d_f, y_name, x_name, x_level){

     return( mean(d_f[d_f[,x_name]!=x_level,y_name]) )
}

 

# fn_mean_oneLevel関数を、指定した説明変数の水準分イテレーション実施して、dfを返す関数
# (df, "目的変数(numeric)名", "説明変数(カテゴリー)名", "LOOによる新しい変数名")
fn_mean_oneLevel_table <- function(d_f, y_name, x_name, x_new_name){

     Num <- c()
     for(lvl in levels(d_f[,x_name])){

          Num <- c(Num, fn_mean_oneLevel(d_f, y_name, x_name, lvl))
     }

     loo_df <- data.frame(levels(d_f[,x_name]), Num)
     names(loo_df) <- c(x_name, x_new_name)
     return( loo_df )
}

 

# 元のdfにLOOを付け、返す関数
# (df, "目的変数(numeric)名", "説明変数(カテゴリー)名", "LOOによる新しい変数名", ノイズの大きさ, "テスト変数名")
# テスト変数 - 0:train&validation 1:test
# 学習データのみにある水準が含まれ、テストデータにその水準が含まれない場合の処理はしていない
fn_LOO2 <- function(d_f, y_name, x_name, x_new_name, noise_ratio=0.1, test_name=NULL){

     #テストデータ(水準1)を除く
     if(is.null(test_name)){      # テストデータを指定しない場合
          d_f_train <- d_f
     }
     else if(sum(d_f[,test_name]==1)>0) {      # テストデータ用の水準(1)を指定した場合
          d_f_train <- d_f[d_f[,test_name]!=1,]
     }
     else { stop("--- Test data set??? ---") }      # 不明な変数を指定した場合

     #LOOテーブルの作成
     loo_table <- fn_mean_oneLevel_table(d_f_train, y_name, x_name, x_new_name)

     # merge
     d_f_merge <- merge(d_f, loo_table, by=x_name, all.x=TRUE, sort=FALSE)

 

     # ノイズのSD
     sdiv <- sd(d_f_merge[,x_new_name])
     # sdivをベースに正規分布からノイズを作成、noise_ratioで大きさを調整
     noise <- rnorm(nrow(d_f), mean=0, sd=sdiv*noise_ratio)
     # ノイズを加える
     d_f_merge[,x_new_name] <- d_f_merge[,x_new_name] + noise

     return( d_f_merge )
}

###########################
###########################

fn_LOO2(データフレーム, "目的変数(numeric)名", "説明変数(カテゴリー)名", "LOOによる新しい変数名", ノイズの大きさ, "テスト変数名(0:学習, 1:テスト)") として実行。


'ノイズの大きさ'とは、0から1の間の値を取り、0だと全くノイズを与えないという意味。

LOOによって各水準の値を計算した後、それらのばらつき(SD)を基準として、そのSDの何割かをノイズとして加えておく。

rnorm(n, mean=0, SD*'ノイズの大きさ')といったノイズを発生させ、LOOで計算した変数に加えておく。


adult2_5 <- fn_LOO2(temp_df, y_name="income", x_name="native.country", x_new_name="native.country2_LOO", noise_ratio=0.1, test_name="Test")

ノイズは0.1にしてみた。


変数が増えてややこしくなってきたので、順序を入れ替えて整理しておく。

adult3 <- adult2_5[,c("Test", "age", "age2", "workclass", "fnlwgt", "fnlwgt2", "education", "education.num", "marital.status", "occupation", "relationship", "race", "sex", "capital.gain", "capital.gain2", "capital.loss",  "capital.loss2", "hours.per.week", "hours.par.week2", "native.country", "native.country2_SR",

"native.country2_LOO", "income_class","income")]


str(adult3)

'data.frame':   32561 obs. of  24 variables:
 $ Test               : int  0 1 0 0 0 1 0 0 0 0 ...
 $ age                : int  26 41 35 20 27 35 38 37 37 36 ...
 $ age2               : Factor w/ 11 levels "-19","20-24",..: 3 6 5 2 3 5 5 5 5 5 ...
 $ workclass          : Factor w/ 7 levels "?","Federal-gov",..: 4 4 2 1 4 4 4 3 5 4 ...
 $ fnlwgt             : int  56929 207685 35309 200061 269354 187415 220237 244803 188774 31438 ...
 $ fnlwgt2            : num  56929 207685 35309 200061 269354 ...
 $ education          : Factor w/ 16 levels "10th","11th",..: 10 12 10 16 10 5 9 12 2 12 ...
 $ education.num      : int  13 9 13 10 13 3 11 9 7 9 ...
 $ marital.status     : Factor w/ 5 levels "Divorced","Married",..: 3 2 3 3 2 2 2 2 2 1 ...
 $ occupation         : Factor w/ 13 levels "?","Adm-clerical",..: 4 13 12 1 11 7 3 3 11 13 ...
 $ relationship       : Factor w/ 6 levels "Husband","Not-in-family",..: 2 1 2 4 1 1 1 1 2 5 ...
 $ race               : Factor w/ 4 levels "Asian-Pac-Islander",..: 2 4 1 4 4 1 4 4 4 4 ...
 $ sex                : Factor w/ 2 levels "Female","Male": 2 2 2 1 2 2 2 2 2 2 ...
 $ capital.gain       : int  0 0 0 0 0 0 0 0 0 0 ...
 $ capital.gain2      : int  0 0 0 0 0 0 0 0 0 0 ...
 $ capital.loss       : int  0 0 0 0 0 0 0 0 0 0 ...
 $ capital.loss2      : Factor w/ 6 levels "<1000","<1250",..: 1 1 1 1 1 1 1 1 1 1 ...
 $ hours.per.week     : int  50 45 28 40 25 50 40 40 60 43 ...
 $ hours.par.week2    : int  50 45 28 40 25 50 40 40 60 43 ...
 $ native.country     : Factor w/ 42 levels "?","Cambodia",..: 1 1 1 1 1 1 1 1 1 1 ...
 $ native.country2_SR : num  19 19 19 19 19 19 19 19 19 19 ...
 $ native.country2_LOO: num  0.241 0.24 0.241 0.242 0.241 ...
 $ income_class       : Factor w/ 2 levels "<=50K",">50K": 1 1 1 1 1 1 2 1 1 1 ...
 $ income             : int  0 0 0 0 0 0 1 0 0 0 ...

 

では、再びモデルを作成してみる。

学習(+検証)、テストデータの分割。
train3 <- adult3[adult3[,"Test"]==0,]
test3 <- adult3[adult3[,"Test"]==1,]


"native.country2_SR"と"native.country2_LOO"で、別々にモデルを作ってみる。
(元々の変数重要度が、最初に作ったモデルでも4.893%しかないので、大きく違ってくるとは考えにくい)

まず、Simple Ranking,"native.country2_SR"。

キャップをかけたりカテゴリー化した変数ももちろん用いる。
x_col_names <- c("age2", "workclass", "fnlwgt2", "education", "education.num", "marital.status", "occupation", "relationship", "race", "sex", "capital.gain2", "capital.loss2", "hours.par.week2", "native.country2_SR")

x_cols <- train3[, x_col_names]
y_col <- train3[,"income"]

 

train3データセットのうち、最初の16264行を学習、残りを検証。これは前回と同じ。
for_train <- 16264


前回のPatial Dependence Plotの結果を見て、単調増加(monotone increasing)に見えた変数には、その制約を課してみる。
"education.num"と"capital.gain"をmonotone increasingとする。ベクトルに1を入れる。
単調減少(monotone decreasing)の場合は-1を入れる。
"native.country2_SR"は、学習(+検証)データにおいて目的変数の値が高いものから順に並べれれているので、monotone decreasingの制約をかける。
単調増加/減少にすることにより、オーバーフィッティングの回避が行われるのだと思われる。

monotone <- c(0,0,0,0,1,0,0,0,0,0,1,0,0,-1)

0はこの制約をかけないという意味。

引数var.monotoneに、このベクトルを渡す。


mod2_1 <- gbm.fit(
     x = x_cols,
     y = y_col,
     distribution = "bernoulli",
     var.monotone = monotone,
     bag.fraction = 0.5,
     n.trees = 5000,
     interaction.depth = 7,
     n.minobsinnode = 100,
     shrinkage = 0.005,
     nTrain = for_train
)

 


変数重要度。
summary(mod2_1) 

                                  var    rel.inf
relationship             relationship 30.1044827
capital.gain2           capital.gain2 15.2426028
occupation                 occupation 14.1707533
education                   education 11.5529362
age2                             age2  8.6704758
education.num           education.num  4.4322248
capital.loss2           capital.loss2  3.9274653
hours.par.week2       hours.par.week2  3.6641404
fnlwgt2                       fnlwgt2  3.1119212
workclass                   workclass  2.3148916
marital.status         marital.status  1.4242388
native.country2_SR native.country2_SR  0.8392975
sex                               sex  0.3368926
race                             race  0.2076771

native.countryをランク化したことによって、重要度が落ちた。他はあまり変化がない。


最適なTreeの数。
gbm.perf(mod2_1, plot.it=TRUE)
3132

 

では、予測を行う。

学習データ。
x_cols_train <- x_cols[1:for_train,]
y_col_train <- y_col[1:for_train]
検証データ。
x_cols_varidation <- x_cols[(for_train+1):nrow(x_cols),]
y_col_varidation <- y_col[(for_train+1):nrow(x_cols)]


学習データに対する予測値。
pred_mod2_1_train <- predict(mod2_1, newdata=x_cols_train, n.trees=gbm.perf(mod2_1, plot.it=FALSE), type="response")

検証データに対する予測値。
pred_mod2_1_varidation <- predict(mod2_1, newdata=x_cols_varidation, n.trees=gbm.perf(mod2_1, plot.it=FALSE), type="response")

テストデータに対する予測値。
test_x <- test3[,x_col_names]
test_y <- test3[,"income"]

pred_mod2_1_test <- predict(mod2_1, newdata=test_x, n.trees=gbm.perf(mod2_1, plot.it=FALSE), type="response")

 

学習データでのAUC。
fn_AUC(pred=pred_mod2_1_train, obs=y_col_train)
0.9332073

検証データでのAUC。
fn_AUC(pred=pred_mod2_1_varidation, obs=y_col_varidation)
0.9198591

テストデータでのAUC。
fn_AUC(pred=pred_mod2_1_test, obs=test_y)
0.9192887

fn_AUCに関しては、前回参照。

 

個人的な経験だが、今回示したようなオーバーフィッティングの回避を試みると、学習とテストデータのパフォーマンス(AUCなど)の幅が狭くなり、検証とテストデータのパフォーマンスもかなり近くなるといった傾向がある。
前回の結果と見比べてみてほしい。

 

 

一応、LOOも試しておく。

x_col_names <- c("age2", "workclass", "fnlwgt2", "education", "education.num", "marital.status", "occupation", "relationship", "race", "sex", "capital.gain2", "capital.loss2", "hours.par.week2", "native.country2_LOO")

x_cols <- train3[, x_col_names]
y_col <- train3[,"income"]

for_train <- 16264

 

monotone <- c(0,0,0,0,1,0,0,0,0,0,1,0,0,1)
"native.country2_LOO"はmonotone increasingとしておく。

同様にgbmを実行し、結果をmod2_2に格納。


変数重要度。
summary(mod2_2)

                                    var    rel.inf
relationship               relationship 30.0720371
capital.gain2             capital.gain2 15.2900749
occupation                   occupation 14.2201974
education                     education 11.1203147
age2                               age2  8.7793286
education.num             education.num  5.1393013
capital.loss2             capital.loss2  3.9263051
hours.par.week2         hours.par.week2  3.6922352
fnlwgt2                         fnlwgt2  3.0737037
workclass                     workclass  2.3722413
marital.status           marital.status  1.3842674
sex                                 sex  0.3466526
native.country2_LOO native.country2_LOO  0.3382367
race                               race  0.2451039

 

最適なTreeの数。

gbm.perf(mod2_2, plot.it=TRUE)
3328

 

予測を行い、AUCを計算する。

 

学習データでのAUC。
0.9331306

検証データでのAUC。
0.9196043

テストデータでのAUC。
0.9192822

 

以上、今回はSimple RankingやLOO関数のメモを個人的理由とする。