BUSINESS DATA SCIENCE 3章 Regularization③ 情報量基準を用いた評価(CI)
BUSINESS DATA SCIENCEの続き。
データなどは作者のgitにある。
最近「ビジネスデータサイエンスの教科書」という邦題で日本語版も出たみたいです
- 作者:マット・タディ
- 発売日: 2020/07/22
- メディア: 単行本
前回の記事
この章はなにか
高次元のモデルを作成するときに、オーバーフィットを避けた良い予測をできるモデル選択(変数選択)において注意すべき点が色々とある。
この節はなにか
OOS評価をする際にCVでは時間がかかったり、CV内でのデータ分割のランダム性による評価の変化が気になる場合はCVの代わりに、 Information Criteria(IC:情報量判断？)を使う方法があるので紹介。
Information Criteria(IC) provide an alternative to cross validation for model selection
Information Criteria(IC)はAIC, AICc, BICなどのように、予実の差を情報量基準を使って表現した評価方法となる。
AIC, AICc
AICはIS-devianceとモデルの自由度(df:model degree of freedom)を使った以下で表される。
解釈の仕方として、IS-devianceは小さくなっても複雑なモデルだった場合は、パラメータ数dfが増えることでAICは悪化する（大きくなる)。一般的には、devianceは-2log(L)で表される場合の方が多い(L:最大尤度)。
例えば、glm()のsummary出力のうち以下の部分に注目する(注：先ほどとは別のテキトーなデータを使ってる）。
## Null deviance: 18754 on 53 degrees of freedom ## Residual deviance: 15598 on 52 degrees of freedom ## AIC: 465.21
このとき、作成したモデルのAICは465.21となる。
なお、Null devianceは(Null model:期待値のモデル:最も当てはまり(deviance)の悪いモデルの)deviance - full model deviance(1対1対応のモデル:最も当てはまりのよいモデル)、Residual devianceは(作成したモデルの)deviance - full model deviance(1対1対応のモデル)となる。
また、 degrees of freedom は正確には degrees of freedom left after fit で、n - df となる(本書ではdf:自由度=パラメータ数として書いているが、一般的にはdf = n-パラメータ なような気がするがよくわからん。ややこしいが以下では本書に合わせてdf=パラメータ数とする)。
lassoではパラメータ数は係数が0ではない変数の数となる。しかし、lasso以外のペナルティ項の場合は、全係数が0でないにも関わらず全変数をペナルティ項なしで作成したmodelよりもdf(パラメータ数)は小さくなる、これは係数が0にshrinkされることが理由となっている。
ただし、AICは高次元モデルではオーバーフィットしやすい。
AICはn=∞を想定しているため、「 nが大きいとき に良い推定量となる」が、正確には「 n/df が大きいとき に良い推定量となる」。
AICは前述のように、 となるので、dfが大きいつまりパラメータ数が多い高次元の複雑なモデルほど大きくなりやすいので、例えばn数とdfが強い相関があるときはn数が大きくてもAICは悪くなる(大きくなる)。
AICはISのdevianceを使っているが、OOSの評価として機能している。これは(n=∞のとき)ISとOOSが同一分布からドローされていると仮定できたとき、理論的には2df = IS deviance - OOS devianceとなることが理由となっている。しかし、高次元データではパラメータに対してnが少なくなりがちなため、IS devianceは実際より小さく評価されている。このバイアスも込みで先程の式を書き直すと、真の標準偏差 とISから推定される標準偏差
を用いて、
となる、
これを用いてAICを修正した高次元にも使える評価方法として、AICcが提案されている。
なお、AICcはbig n/df のときにはAICと近似できるため、無難にAICcを使うほうがよい。また、AICcはCV-minのλの結果と近似できる。
BIC
Baysを使ったものとして、BICが存在する、
BIC = deviance + log(n) × df
これは、 を事前情報として想定している。AIC, AICcは予測に対して最適化されているが、BICは事前情報をもとに推定されたデータの真のモデルをどれくらい表しているかを評価している。
これは、保守的な使い方ができ、中くらいのデータ量ではcv-1seと近似できる。ただし、nが大きい場合は過小評価される。
## plot CV results and the various IC ll <- log(spender$lambda) ## the sequence of lambdas par(mfrow=c(2,1)) plot(cv.spender) plot(ll, AIC(spender)/n, xlab="log lambda", ylab="IC/n", pch=21, bg="orange") abline(v=ll[which.min(AIC(spender))], col="orange", lty=3) abline(v=ll[which.min(BIC(spender))], col="green", lty=3) abline(v=ll[which.min(AICc(spender))], col="black", lty=3) points(ll, BIC(spender)/n, pch=21, bg="green") points(ll, AICc(spender)/n, pch=21, bg="black") legend("topleft", bty="n", fill=c("black","orange","green"),legend=c("AICc","AIC","BIC"))
λを変えたときの評価指標、凡例が見えなくなっているが緑がBIC。
上図がCVでの結果となり、左の点線がCV-min、右がCV-1se。下図がCIでの評価でオレンジの線がAIC、黒がAICs、緑がBICでの最小λとなり、CV-minとAICc、CV-1seとBICがある程度同じ位置を示していることが確認できる。
また、各係数とλの関係でまとめると以下。
## all metrics, together in a path plot. plot(spender, col="grey") abline(v=ll[which.min(AICc(spender))], col="black", lty=2) abline(v=ll[which.min(AIC(spender))], col="orange", lty=2) abline(v=ll[which.min(BIC(spender))], col="green", lty=2) abline(v=log(cv.spender$lambda.min), col="blue", lty=2) abline(v=log(cv.spender$lambda.1se), col="purple", lty=2) legend("topright", bty="n", lwd=1, col=c("black","orange","green","blue","purple"), legend=c("AICc","AIC","BIC","CV.min","CV.1se"))
実データ
ホッケー選手のゴールをした履歴データを用いる。
library(gamlr) # loads Matrix as well data(hockey) head(goal) # homegoal season team.away team.home period differential playoffs # 1 0 20022003 DAL EDM 1 0 0 # 2 0 20022003 DAL EDM 1 -1 0 # 3 1 20022003 DAL EDM 2 -2 0 # 4 0 20022003 DAL EDM 2 -1 0 # 5 1 20022003 DAL EDM 3 -2 0 # 6 1 20022003 DAL EDM 3 -1 0
CVで評価する。
# Combine the covariates all together x <- cBind(config,team,player) # cBind binds together two sparse matrices # build 'y': home vs away, binary response y <- goal$homegoal ## to run cv, just use cv.gamlr instead of gamlr cv.nhlreg <- cv.gamlr(x, y, free=1:(ncol(config)+ncol(team)), family="binomial", verb=TRUE, standardize=FALSE) ## plot them together par(mfrow=c(1,2)) plot(cv.nhlreg) plot(cv.nhlreg$gamlr) ## cv.gamlr has included a gamlr object into cv.nhlreg
このとき、configとteam列はモデルに含めたいので、 free オプションにデザインマトリクスを与えることで、configとteamに対してペナルティを付けない処理をおこなっている。
また、 standardize=FALSE をつけることで全変数に対してスケールを整えることで共通のペナルティ項で制御をおこなっている。逆にいえば、 standardize=TRUE の場合はplayer毎に個別のペナルティ項を設定することになる。今回は、player毎に出場回数に差があることから標準偏差が大きく異なるためFALSEにしている。これは特殊ケースなので基本的にはTRUEを設定する。
# Without standardize=FALSE, you'd be multiplying the penalty for each # coefficient (player effect) by that player's standard deviation in onice. # # The players with big SD in onice are guys who play a lot. # Players with small SD are those who play little (almost all zeros). # # So weighting penalty by SD in this case is exactly what you don't want: # a bigger penalty for people with many minutes on ice, a smaller penalty # for those who seldom play. Indeed, running the regression without # standardize=FALSE leads to a bunch of farm teamers coming up tops. nhlreg.std <- gamlr(x, y, free=1:(ncol(config)+ncol(team)), family="binomial", standardize=FALSE) Bstd <- coef(nhlreg.std)[colnames(player),] Bstd[order(Bstd, decreasing=TRUE)[1:10]] # PETER_FORSBERG TYLER_TOFFOLI ONDREJ_PALAT ZIGMUND_PALFFY SIDNEY_CROSBY JOE_THORNTON # 0.7548254 0.6292577 0.6284040 0.4426997 0.4131174 0.3837632 # PAVEL_DATSYUK LOGAN_COUTURE ERIC_FEHR MARTIN_GELINAS # 0.3761981 0.3682103 0.3677283 0.3577613 # NOTE: this is an exceptional case! You almost always want standardize=TRUE.
FALSEのまま話をすすめると、それぞれのλは以下
log(nhlreg$lambda[which.min(AICc(nhlreg))]) # -9.165555 log(nhlreg$lambda[which.min(AIC(nhlreg))]) # -9.165555 log(nhlreg$lambda[which.min(BIC(nhlreg))]) # -6.653644 log(cv.nhlreg$lambda.min) # -8.979488 log(cv.nhlreg$lambda.1se) # -8.328251
AICとAICcは同様。それに対して、CV-minは近い値だが、BICは大きめに出ている。これは大規模データなため過小評価された結果となっている。