Rの機械学習パッケージmlrのチュートリアル2（前処理からチューニングまで）

[Rの機械学習パッケージmlrのチュートリアル（タスクの作成から予測まで） - Qiita](https://qiita.com/nozma/items/bedcb35cba925764247a)の続き。

- 次→[Rの機械学習パッケージmlrのチュートリアル3（ベンチマーク試験から可視化まで） - Qiita](https://qiita.com/nozma/items/01bca01c7832d760e579)
- まとめたやつ→[mlrパッケージチュートリアル - Quick Walkthrough編](https://nozma.github.io/mlr_tutorial/)

前処理
==============

学習アルゴリズムを適用する前にデータに施すあらゆる種類の変換を前処理と呼ぶ。この中には、データの矛盾の発見と解決、欠損値への代入、外れ値の特定・除去・置換、数値データの離散化、カテゴリカルデータからのダミー変数の生成、標準化やBox-Cox変換などのあらゆる種類の変換、次元削減、特徴量の抽出・選択などが含まれる。

`mlr`は前処理に関して幾つかの選択肢を用意している。データフレームやタスクを直接変更するものは比較的単純な部類に属する。以下にその例を示す。中にはこれまでに既に取り上げたものもある。

-   `capLargeValues`: 大きな値や無限大の値の変換。
-   `createDummyFeature`: 因子型特徴量からのダミー変数の生成。
-   `dropFeatures`: 特徴量の削除。
-   `joinClassLevels`: (分類のみ)複数のクラスを併合して、大きな1つのクラスにする。
-   `mergeSmallFactorLevels`: 因子型特徴量において、例数の少ない水準を併合する。
-   `normalizeFeatures`: 正規化には複数の異なったやり方がある。標準化や特定の範囲への再スケールなど。
-   `removeConstantFeatures`: 1つの値しか持っていない特徴量(=定数)を除去する。
-   `subsetTask`: 観測値や特徴量をタスクから除去する。

また、以下のものについては別途チュートリアルを用意してある。

-   特徴量選択
-   欠損値への代入

前処理と学習器を融合する
------------------------

`mlr`のラッパー機能により、学習器と前処理を組み合わせることができる。これは、前処理が学習器に属し、訓練や予測の度に実行されるということを意味する。

このようにすることで非常に便利な点がある。データやタスクの変更なしに、簡単に学習器と前処理の組合せを変えることができるのだ。

また、これは前処理を行ってから学習器のパフォーマンスを測定する際にありがちな一般的な間違いを避けることにもつながる。前処理は学習アルゴリズムとは完全に独立したものだと考えられがちだ。学習器のパフォーマンスを測定する場合を考えてみよう。例えば、クロスバリデーションで雨処理を事前にデータセット全体に対して行い、学習と予測は学習器だけで行うような場合だ。前処理として何が行われたかによっては、評価が楽観的になる危険性がある。例えば、(欠損値への)平均値の代入という前処理が学習器の性能評価前に、データ全体を対象に行われたとすると、これは楽観的なパフォーマンス評価につながる。

前処理にはデータ依存的なものとデータ非依存的なものがあることをはっきりさせておこう。データ依存的な前処理とは、前処理のやり方がデータに依存しており、データセットが異なれば結果も異なるというようなもののことだ。一方でデータ非依存的な前処理は常に結果が同じになる。

データの間違いを修正したり、ID列のような学習に使うべきではないデータ列の除去のような前処理は、明らかにデータ非依存的である。一方、先程例に挙げた欠損値への平均値の代入はデータ依存的である。代入を定数で行うのであれば違うが。

前処理と組み合わせた学習器の性能評価を正しく行うためには、全てのデータ依存的な前処理をリサンプリングに含める必要がある。学習器と前処理を融合させれば、これは自動的に可能になる。

この目的のために、`mlr`パッケージは2つのラッパーを用意している。

-   `makePreprocWrapperCaret`は`caret`パッケージの`preProcess`関数に対するインターフェースを提供するラッパー。
-   `makePreprocWrapper`を使えば、訓練と予測の前の動作を定義することで独自の前処理を作成できる。

これらを使用する前処理は、`normalizeFeatures`などを使う前処理とは異なり、ラップされた学習器に組み込まれる。

-   タスクそのものは変更されない。
-   前処理はデータ全体に対して予め行われるのではなく、リサンプリングなど、訓練とテストの対が発生する毎に実行される。
-   前処理に関わる制御可能なパラメータは、学習器のパラメータと一緒に調整できる。

まずは`makePreprocWrapperCaret`の例から見ていこう。

`makePreprocWrapperCaret`を使用した前処理
-----------------------------------------

`makePreprocWrapperCaret`は`caret`パッケージの`preProcess`関数へのインターフェースだ。`PreProcess`関数は、欠損値への代入やスケール変換やBox-Cox変換、独立主成分分析による次元削減など、様々な手法を提供する関数だ。具体的に何が可能かは`preProcess`関数のヘルプページ([preProcess function | R Documentation](https://www.rdocumentation.org/packages/caret/versions/6.0-78/topics/preProcess))を確認してもらいたい。

まず、`makePreprocWrapperCaret`と`preProcess`の違いを確認しておこう。

-   `makePreprocWrapperCaret`は`preProcess`とほぼ同じ仮引数を持つが、仮引数名に`ppc.`というプレフィックスが付く。
-   上記の例外は`method`引数だ。この引数は`makePreprocWrapperCaret`には無い。その代わりに、本来`method`に渡す前処理に関するオプションは、対応する仮引数に論理値を指定することで制御する。

例を見よう。`preProcess`では行列またはデータフレーム`x`に対して、次のように前処理を行う。

``` r
preProcess(x, method= c("knnInpute", "pca"), pcaComp = 10)
```

一方、`makePreporcWrapperCaret`では、`Learner`クラスのオブジェクトまたはクラスの名前(`"classif.lda"`など)を引数にとって、次のように前処理を指定する。

``` r
makePreprocWrapperCaret(learner, ppc.knnImpute = TRUE, ppc.pca = TRUE, ppc.pcaComp = 10)
```

この例のように複数の前処理(注: kNNを使った代入と主成分分析)を有効にした場合、それらは特定の順序で実行される。詳細は`preProcess`関数のヘルプを確認してほしい(訳注: Details後半の"The operations are applied in this order:..."以下。主成分分析は代入後に実施。)。

以下に主成分分析による次元削減の例を示そう。これは無闇に使用して良い手法ではないが、高次元のデータで問題が起こるような学習器や、データの回転が有用な学習器に対しては有効である。

例では`soner.task`を用いる。これは208の観測値と60の特徴量を持つ。

``` r
sonar.task
```

    $> Supervised task: Sonar-example
    $> Type: classif
    $> Target: Class
    $> Observations: 208
    $> Features:
    $> numerics  factors  ordered 
    $>       60        0        0 
    $> Missings: FALSE
    $> Has weights: FALSE
    $> Has blocking: FALSE
    $> Classes: 2
    $>   M   R 
    $> 111  97 
    $> Positive class: M

今回は、`MASS`パッケージによる二次判別分析と、主成分分析による前処理を組み合わせる。また、閾値として0.9を設定する。これはつまり、主成分が累積寄与率90%を保持しなければならないという指示になる。データは主成分分析の前に自動的に標準化される。

``` r
lrn = makePreprocWrapperCaret("classif.qda", ppc.pca = TRUE, ppc.thresh = 0.9)
lrn
```

    $> Learner classif.qda.preproc from package MASS
    $> Type: classif
    $> Name: ; Short name: 
    $> Class: PreprocWrapperCaret
    $> Properties: twoclass,multiclass,numerics,factors,prob
    $> Predict-Type: response
    $> Hyperparameters: ppc.BoxCox=FALSE,ppc.YeoJohnson=FALSE,ppc.expoTrans=FALSE,ppc.center=TRUE,ppc.scale=TRUE,ppc.range=FALSE,ppc.knnImpute=FALSE,ppc.bagImpute=FALSE,ppc.medianImpute=FALSE,ppc.pca=TRUE,ppc.ica=FALSE,ppc.spatialSign=FALSE,ppc.thresh=0.9,ppc.na.remove=TRUE,ppc.k=5,ppc.fudge=0.2,ppc.numUnique=3

ラップされた学習器を`soner.task`によって訓練する。訓練したモデルを確認することで、22の主成分が訓練に使われたことがわかるだろう。

``` r
mod = train(lrn, sonar.task)
mod
```

    $> Model for learner.id=classif.qda.preproc; learner.class=PreprocWrapperCaret
    $> Trained on: task.id = Sonar-example; obs = 208; features = 60
    $> Hyperparameters: ppc.BoxCox=FALSE,ppc.YeoJohnson=FALSE,ppc.expoTrans=FALSE,ppc.center=TRUE,ppc.scale=TRUE,ppc.range=FALSE,ppc.knnImpute=FALSE,ppc.bagImpute=FALSE,ppc.medianImpute=FALSE,ppc.pca=TRUE,ppc.ica=FALSE,ppc.spatialSign=FALSE,ppc.thresh=0.9,ppc.na.remove=TRUE,ppc.k=5,ppc.fudge=0.2,ppc.numUnique=3

``` r
getLearnerModel(mod)
```

    $> Model for learner.id=classif.qda; learner.class=classif.qda
    $> Trained on: task.id = Sonar-example; obs = 208; features = 22
    $> Hyperparameters:

``` r
getLearnerModel(mod, more.unwrap = TRUE)
```

    $> Call:
    $> qda(f, data = getTaskData(.task, .subset, recode.target = "drop.levels"))
    $> 
    $> Prior probabilities of groups:
    $>         M         R 
    $> 0.5336538 0.4663462 
    $> 
    $> Group means:
    $>          PC1        PC2        PC3         PC4         PC5         PC6
    $> M  0.5976122 -0.8058235  0.9773518  0.03794232 -0.04568166 -0.06721702
    $> R -0.6838655  0.9221279 -1.1184128 -0.04341853  0.05227489  0.07691845
    $>          PC7         PC8        PC9       PC10        PC11          PC12
    $> M  0.2278162 -0.01034406 -0.2530606 -0.1793157 -0.04084466 -0.0004789888
    $> R -0.2606969  0.01183702  0.2895848  0.2051963  0.04673977  0.0005481212
    $>          PC13       PC14        PC15        PC16        PC17        PC18
    $> M -0.06138758 -0.1057137  0.02808048  0.05215865 -0.07453265  0.03869042
    $> R  0.07024765  0.1209713 -0.03213333 -0.05968671  0.08528994 -0.04427460
    $>          PC19         PC20        PC21         PC22
    $> M -0.01192247  0.006098658  0.01263492 -0.001224809
    $> R  0.01364323 -0.006978877 -0.01445851  0.001401586

二次判別分析について、主成分分析を使う場合と使わない場合をベンチマーク試験により比較してみよう。今回の例では各クラスの例数が少ないので、二次判別分析の際のエラーを防ぐためにリサンプリングにおいて層別サンプリングを行っている点に注意してほしい。リサンプリング手法については後ほど解説する。

``` r
rin = makeResampleInstance("CV", iters = 3, stratify = TRUE, task = sonar.task)
res = benchmark(list("classif.qda", lrn), sonar.task, rin, show.info = FALSE)
res
```

    $>         task.id          learner.id mmce.test.mean
    $> 1 Sonar-example         classif.qda      0.3505176
    $> 2 Sonar-example classif.qda.preproc      0.2213251

今回の例では、二次判別分析に対して主成分分析による前処理が効果的だったことがわかる。

前処理オプションと学習器パラメータの連結チューニング
----------------------------------------------------

今の例をもう少し最適化できないか考えてみよう。今回、任意に設定した0.9という閾値によって、主成分は22になった。しかし、主成分の数はもっと多いほうが良いかもしれないし、少ないほうが良いかもしれない。また、`qda`関数にはクラス共分散行列やクラス確率の推定方法を制御するためのいくつかのオプションがある。

学習機と前処理のパラメータは、連結してチューニングすることができる。まずは、ラップされた学習器の全てのパラメータを`getParamSet`関数で確認してみよう。

``` r
getParamSet(lrn)
```

    $>                      Type len     Def                      Constr Req
    $> ppc.BoxCox        logical   -   FALSE                           -   -
    $> ppc.YeoJohnson    logical   -   FALSE                           -   -
    $> ppc.expoTrans     logical   -   FALSE                           -   -
    $> ppc.center        logical   -    TRUE                           -   -
    $> ppc.scale         logical   -    TRUE                           -   -
    $> ppc.range         logical   -   FALSE                           -   -
    $> ppc.knnImpute     logical   -   FALSE                           -   -
    $> ppc.bagImpute     logical   -   FALSE                           -   -
    $> ppc.medianImpute  logical   -   FALSE                           -   -
    $> ppc.pca           logical   -   FALSE                           -   -
    $> ppc.ica           logical   -   FALSE                           -   -
    $> ppc.spatialSign   logical   -   FALSE                           -   -
    $> ppc.thresh        numeric   -    0.95                    0 to Inf   -
    $> ppc.pcaComp       integer   -       -                    1 to Inf   -
    $> ppc.na.remove     logical   -    TRUE                           -   -
    $> ppc.k             integer   -       5                    1 to Inf   -
    $> ppc.fudge         numeric   -     0.2                    0 to Inf   -
    $> ppc.numUnique     integer   -       3                    1 to Inf   -
    $> ppc.n.comp        integer   -       -                    1 to Inf   -
    $> method           discrete   -  moment            moment,mle,mve,t   -
    $> nu                numeric   -       5                    2 to Inf   Y
    $> predict.method   discrete   - plug-in plug-in,predictive,debiased   -
    $>                  Tunable Trafo
    $> ppc.BoxCox          TRUE     -
    $> ppc.YeoJohnson      TRUE     -
    $> ppc.expoTrans       TRUE     -
    $> ppc.center          TRUE     -
    $> ppc.scale           TRUE     -
    $> ppc.range           TRUE     -
    $> ppc.knnImpute       TRUE     -
    $> ppc.bagImpute       TRUE     -
    $> ppc.medianImpute    TRUE     -
    $> ppc.pca             TRUE     -
    $> ppc.ica             TRUE     -
    $> ppc.spatialSign     TRUE     -
    $> ppc.thresh          TRUE     -
    $> ppc.pcaComp         TRUE     -
    $> ppc.na.remove       TRUE     -
    $> ppc.k               TRUE     -
    $> ppc.fudge           TRUE     -
    $> ppc.numUnique       TRUE     -
    $> ppc.n.comp          TRUE     -
    $> method              TRUE     -
    $> nu                  TRUE     -
    $> predict.method      TRUE     -

`ppc.`というプレフィックスのついたものが前処理のパラメータで、他が`qda`関数のパラメータだ。主成分分析の閾値を`ppc.thresh`を使って調整する代わりに、主成分の数そのものを`ppc.pcaComp`を使って調整できる。さらに、`qda`関数に対しては、2種類の事後確率推定法(通常のプラグイン推定と不偏推定)を試してみよう。

今回は解像度10でグリッドサーチを行ってみよう。もっと解像度を高くしたくなるかもしれないが、今回はあくまでデモだ。

``` r
ps = makeParamSet(
  makeIntegerParam("ppc.pcaComp", lower = 1, upper = getTaskNFeats(sonar.task)),
  makeDiscreteParam("predict.method", values = c("plug-in", "debiased"))
)
ctrl = makeTuneControlGrid(resolution = 10)
res = tuneParams(lrn, sonar.task, rin, par.set = ps, control = ctrl, show.info = FALSE)
res
```

    $> Tune result:
    $> Op. pars: ppc.pcaComp=21; predict.method=plug-in
    $> mmce.test.mean=0.212

``` r
as.data.frame(res$opt.path)[1:3]
```

    $>    ppc.pcaComp predict.method mmce.test.mean
    $> 1            1        plug-in      0.5284334
    $> 2            8        plug-in      0.2311939
    $> 3           14        plug-in      0.2118703
    $> 4           21        plug-in      0.2116632
    $> 5           27        plug-in      0.2309869
    $> 6           34        plug-in      0.2739821
    $> 7           40        plug-in      0.2933057
    $> 8           47        plug-in      0.3029676
    $> 9           53        plug-in      0.3222912
    $> 10          60        plug-in      0.3505176
    $> 11           1       debiased      0.5579020
    $> 12           8       debiased      0.2502415
    $> 13          14       debiased      0.2503796
    $> 14          21       debiased      0.2550725
    $> 15          27       debiased      0.2792271
    $> 16          34       debiased      0.3128364
    $> 17          40       debiased      0.2982747
    $> 18          47       debiased      0.2839199
    $> 19          53       debiased      0.3224983
    $> 20          60       debiased      0.3799172

`"plug-in"`と`"debiased"`のいずれでも少なめ(27以下)の主成分が有効で、`"plug-in"`の方が若干エラー率が低いようだ。

独自の前処理ラッパーを書く
--------------------------

`makePreprocWrapperCaret`で不満があれば、`makePreprocWrapper`関数で独自の前処理ラッパーを作成できる。

ラッパーに関するチュートリアルでも説明しているが、ラッパーは**訓練**と**予測**という2つのメソッドを使って実装される。前処理ラッパーの場合は、メソッドは学習と予測の前に何をするかを指定するものであり、これは完全にユーザーが指定する。

以下に例として、訓練と予測の前にデータの中心化とスケーリングを行うラッパーの作成方法を示そう。k最近傍法やサポートベクターマシン、ニューラルネットワークなどは通常スケーリングされた特徴量を必要とする。多くの組み込みスケーリング手法は、データセットを事前にスケーリングし、テストデータセットもそれに従ってスケーリングされる。以下では、学習器にスケーリングオプションを追加し、`scale`関数と組み合わせる方法を示す。

今回この単純な例を選んだのはあくまで説明のためだ。中心化とスケーリングは`makePreprocWrapperCaret`でも可能だということに注意してほしい。

### 訓練関数の指定

**訓練**(ステップで使う)関数は以下の引数を持つ関数でなければならない。

-   `data`: 全ての特徴量と目的変数を列として含むデータフレーム。
-   `target`: `data`に含まれる目的変数の名前。
-   `args`: 前処理に関わるその他の引数とパラメータのリスト。

この関数は`$data`と`$control`を要素として持つリストを返す必要がある。`$data`は前処理されたデータセットを、`$control`には予測のために必要な全ての情報を格納する。

スケーリングのための訓練関数の定義例を以下に示す。これは数値型の特徴量に対して`scale`関数を呼び出し、スケーリングされたデータと関連するスケーリングパラメータを返す。

`args`は`scale`関数の引数である`center`と`scale`引数を含み、予測で使用するためにこれを`$control`スロットに保持する。これらの引数は、論理値または数値型の列の数に等しい長さの数値型ベクトルで指定する必要がある。`center`引数は数値を渡された場合にはその値を各データから引くが、`TRUE`が指定された場合には平均値を引く。`scale`引数は数値を渡されるとその値で各データを割るが、`TRUE`の場合は標準偏差か二乗平均平方根を引く(いずれになるかは`center`引数に依存する)。2つの引数のいずれかor両方に`TRUE`が指定された場合には、この値を予測の段階で使用するためには返り値の`$control`スロットに保持しておく必要があるという点に注意しよう。

``` r
trainfun = function(data, target, args = list(center, scale)){
  ## 数値特徴量を特定する
  cns = colnames(data)
  nums = setdiff(cns[sapply(data, is.numeric)], target)
  ## 数値特徴量を抽出し、scale関数を呼び出す
  x = as.matrix(data[, nums, drop = FALSE])
  x = scale(x, center = args$center, scale = args$scale)
  ## スケーリングパラメータを後で予測に使うためにcontrolに保持する
  control = args
  if(is.logical(control$center) && control$center){
    control$center = attr(x, "scaled:center")
  }
  if(is.logical(control$scale) && control$scale){
    control$scale = attr(x, "scaled:scale")
  }
  ## 結果をdataにまとめる
  data = data[, setdiff(cns, nums), drop = FALSE]
  data = cbind(data, as.data.frame(x))
  return(list(data = data, control = control))
}
```

### 予測関数の指定

**予測**(ステップで使う)関数は以下の引数を持つ必要がある。

-   `data`: 特徴量のみをもつデータフレーム。(予測ステップでは目的変数の値は未知なのが普通だ。)
-   `target`: 目的変数の名前。
-   `args`: 訓練関数に渡された`args`。
-   `control`: 訓練関数が返したもの。

この関数は前処理済みのデータを返す。

今回の例では、予測関数は数値特徴量を訓練ステージで`control`に保持されたパラメータを使ってスケーリングする。

``` r
predictfun = function(data, target, args, control){
  ## 数値特徴量の特定
  cns = colnames(data)
  nums = cns[sapply(data, is.numeric)]
  ## データから数値特徴量を抽出してscale関数を適用する
  x = as.matrix(data[, nums, drop = FALSE])
  x = scale(x, center = control$center, scale = control$scale)
  ## dataにまとめて返す
  data = data[, setdiff(cns, nums), drop = FALSE]
  data = cbind(data, as.data.frame(x))
  return(data)
}
```

### 前処理ラッパーの作成

以下では、ニューラルネットワークによる回帰(これは自前のスケーリングオプションを持たない)をベースの学習器として前処理ラッパーを作成する。

先に定義した**訓練**および**予測**関数を`makePreprocWrapper`関数の`train`と`predict`引数に渡す。`par.vals`には、訓練関数の`args`に渡すパラメータをリストとして渡す。

``` r
lrn = makeLearner("regr.nnet", trace = FALSE, decay = 1e-02)
lrn = makePreprocWrapper(lrn, train = trainfun, predict = predictfun,
                         par.vals = list(center = TRUE, scale = TRUE))
```

データセット`BostonHousing`を対象にして、スケーリングの有無による平均二乗誤差の違いを確認してみよう。

``` r
rdesc = makeResampleDesc("CV", iters = 3)

## スケーリングあり(上で前処理を指定した)
r = resample(lrn, bh.task, resampling = rdesc, show.info = FALSE)
r
```

    $> Resample Result
    $> Task: BostonHousing-example
    $> Learner: regr.nnet.preproc
    $> Aggr perf: mse.test.mean=  18
    $> Runtime: 0.137429

``` r
## 前処理無しの学習器を再度作る
lrn = makeLearner("regr.nnet", trace = FALSE, decay = 1e-02)
r = resample(lrn, bh.task, resampling = rdesc, show.info = FALSE)
r
```

    $> Resample Result
    $> Task: BostonHousing-example
    $> Learner: regr.nnet
    $> Aggr perf: mse.test.mean=41.5
    $> Runtime: 0.101203

### 前処理と学習器のパラメータを連結してチューニングする

前処理のオプションをどのように設定すれば特定のアルゴリズムに対して上手くいくのかということは、明確には分からないことが多い。`makePreprocWrapperCaret`の例で、既に前処理と学習器のパラメータを両方ともチューニングする方法を既に見た。

スケーリングの例では、ニューラルネットに対してスケーリングと中心化の両方を行うのが良いのか、いずれか片方なのか、あるいは行わないほうが良いのかという点を確認することができる。`center`と`scale`をチューニングするためには、適切な種類の`LearnerParam`をパラメータセットに追加する必要がある。

前述のように、`center`と`scale`には数値型か論理値型のいずれかを指定できるが、今回は論理値型のパラメータとしてチューニングしよう。

``` r
lrn = makeLearner("regr.nnet", trace = FALSE)
lrn = makePreprocWrapper(lrn, train = trainfun, predict = predictfun,
                         par.set = makeParamSet(
                           makeLogicalLearnerParam("center"),
                           makeLogicalLearnerParam("scale")
                         ),
                         par.vals = list(center = TRUE, scale = TRUE))
lrn
```

    $> Learner regr.nnet.preproc from package nnet
    $> Type: regr
    $> Name: ; Short name: 
    $> Class: PreprocWrapper
    $> Properties: numerics,factors,weights
    $> Predict-Type: response
    $> Hyperparameters: size=3,trace=FALSE,center=TRUE,scale=TRUE

今回はグリッドサーチで`nnet`の`decay`パラメータと`scale`の`center`と`scale`パラメータをチューニングする。

``` r
rdesc = makeResampleDesc("Holdout")
ps = makeParamSet(
  makeDiscreteLearnerParam("decay", c(0, 0.05, 0.1)),
  makeLogicalLearnerParam("center"),
  makeLogicalLearnerParam("scale")
)
crrl = makeTuneControlGrid()
res = tuneParams(lrn, bh.task, rdesc, par.set = ps, control = ctrl, show.info = FALSE)
res
```

    $> Tune result:
    $> Op. pars: decay=0.05; center=TRUE; scale=TRUE
    $> mse.test.mean=11.2

``` r
as.data.frame(res$opt.path)
```

    $>    decay center scale mse.test.mean dob eol error.message exec.time
    $> 1      0   TRUE  TRUE      57.95746   1  NA          <NA>     0.039
    $> 2   0.05   TRUE  TRUE      11.23583   2  NA          <NA>     0.042
    $> 3    0.1   TRUE  TRUE      15.44886   3  NA          <NA>     0.043
    $> 4      0  FALSE  TRUE      84.89302   4  NA          <NA>     0.019
    $> 5   0.05  FALSE  TRUE      16.63278   5  NA          <NA>     0.041
    $> 6    0.1  FALSE  TRUE      13.80628   6  NA          <NA>     0.043
    $> 7      0   TRUE FALSE      64.98619   7  NA          <NA>     0.029
    $> 8   0.05   TRUE FALSE      55.94930   8  NA          <NA>     0.040
    $> 9    0.1   TRUE FALSE      26.67453   9  NA          <NA>     0.048
    $> 10     0  FALSE FALSE      63.27422  10  NA          <NA>     0.023
    $> 11  0.05  FALSE FALSE      34.35454  11  NA          <NA>     0.044
    $> 12   0.1  FALSE FALSE      42.57609  12  NA          <NA>     0.043

### 前処理ラッパー関数

よい前処理ラッパーを作成したのであれば、それを関数としてカプセル化するのは良いアイデアだ。他の人も使えると便利だろうから`mlr`に追加して欲しい、というのであれば[Issues · mlr-org/mlr](https://github.com/mlr-org/mlr/issues)からコンタクトして欲しい。

``` r
makePreprocWrapperScale = function(learner, center = TRUE, scale = TRUE) {
  trainfun = function(data, target, args = list(center, scale)) {
    cns = colnames(data)
    nums = setdiff(cns[sapply(data, is.numeric)], target)
    x = as.matrix(data[, nums, drop = FALSE])
    x = scale(x, center = args$center, scale = args$scale)
    control = args
    if (is.logical(control$center) && control$center)
      control$center = attr(x, "scaled:center")
    if (is.logical(control$scale) && control$scale)
      control$scale = attr(x, "scaled:scale")
    data = data[, setdiff(cns, nums), drop = FALSE]
    data = cbind(data, as.data.frame(x))
    return(list(data = data, control = control))
  }
  predictfun = function(data, target, args, control) {
    cns = colnames(data)
    nums = cns[sapply(data, is.numeric)]
    x = as.matrix(data[, nums, drop = FALSE])
    x = scale(x, center = control$center, scale = control$scale)
    data = data[, setdiff(cns, nums), drop = FALSE]
    data = cbind(data, as.data.frame(x))
    return(data)
  }
  makePreprocWrapper(
    learner,
    train = trainfun,
    predict = predictfun,
    par.set = makeParamSet(
      makeLogicalLearnerParam("center"),
      makeLogicalLearnerParam("scale")
    ),
    par.vals = list(center = center, scale = scale)
  )
}

lrn = makePreprocWrapperScale("classif.lda")
train(lrn, iris.task)
```

    $> Model for learner.id=classif.lda.preproc; learner.class=PreprocWrapper
    $> Trained on: task.id = iris-example; obs = 150; features = 4
    $> Hyperparameters: center=TRUE,scale=TRUE

学習器の性能を評価する
======================

`mlr`には学習機の予測性能について様々な側面から評価する方法が備えられている。性能指標は、`predict`の返すオブジェクトと目的の性能指標を指定して`performance`関数を呼び出すことで計算できる。

利用可能な性能指標
------------------

`mlr`はすべての種類の学習問題に対して多数の性能指標を提供している。分類問題に対する典型的な性能指標としては、平均誤分類率(mmce)、精度(acc)、ROC曲線などが使える。回帰問題に対しては、平均二乗偏差(mse)、平均絶対誤差(mae)などが一般に使用される。他にもクラスタリング問題では、Dunn Index(dunn)が、生存時間分析に対してはConcordance Index(cindex)が、コスト考慮型予測問題ではMisclassification Penalty(mcp)など、様々な指標が利用可能である。また、訓練にかかった時間(timetrain)、予測にかかった時間(timepredict)、その合計(timeboth)も性能指標の一つとしてアクセスできる。

どのような指標が実装されているかについては、[Implemented Performance Measures - mlr tutorial](https://mlr-org.github.io/mlr-tutorial/devel/html/measures/index.html)および[measures function | R Documentation](https://www.rdocumentation.org/packages/mlr/versions/2.10/topics/measures)を確認してもらいたい。

もし新たな指標を実装したり、標準的でない誤分類コストを指標に含めたいと思うのであれば、[Create Custom Measures - mlr tutorial](https://mlr-org.github.io/mlr-tutorial/devel/html/create_measure/index.html)を見てもらいたい。

指標の一覧
----------

各指標の詳細については上述のImplemented Performance Measuresを確認してもらうとして、特定のプロパティを持つ指標や、特定のタスクに利用可能な指標を確認したければ`listMeasures`関数を使うと良い。

``` r
## 多クラス問題に対する分類指標
listMeasures("classif", properties = "classif.multi")
```

    $>  [1] "kappa"            "multiclass.brier" "multiclass.aunp" 
    $>  [4] "multiclass.aunu"  "qsr"              "ber"             
    $>  [7] "logloss"          "wkappa"           "timeboth"        
    $> [10] "timepredict"      "acc"              "lsr"             
    $> [13] "featperc"         "multiclass.au1p"  "multiclass.au1u" 
    $> [16] "ssr"              "timetrain"        "mmce"

``` r
## iris.taskに対する分類指標
listMeasures(iris.task)
```

    $>  [1] "kappa"            "multiclass.brier" "multiclass.aunp" 
    $>  [4] "multiclass.aunu"  "qsr"              "ber"             
    $>  [7] "logloss"          "wkappa"           "timeboth"        
    $> [10] "timepredict"      "acc"              "lsr"             
    $> [13] "featperc"         "multiclass.au1p"  "multiclass.au1u" 
    $> [16] "ssr"              "timetrain"        "mmce"

簡便のため、それぞれの学習問題に対しては、よく使われる指標がデフォルトとして指定してある。例えば回帰では平均二乗偏差が、分類では平均誤分類率がデフォルトだ。何がデフォルトであるかは`getDefaultMeasure`関数を使うと確認できる。また、この関数のヘルプでデフォルトに使用される指標の一覧が確認できる。

``` r
## iris.taskのデフォルト指標
getDefaultMeasure(iris.task)
```

    $> Name: Mean misclassification error
    $> Performance measure: mmce
    $> Properties: classif,classif.multi,req.pred,req.truth
    $> Minimize: TRUE
    $> Best: 0; Worst: 1
    $> Aggregated by: test.mean
    $> Note: Defined as: mean(response != truth)

``` r
## 回帰のデフォルト指標
getDefaultMeasure(makeLearner("regr.lm"))
```

    $> Name: Mean of squared errors
    $> Performance measure: mse
    $> Properties: regr,req.pred,req.truth
    $> Minimize: TRUE
    $> Best: 0; Worst: Inf
    $> Aggregated by: test.mean
    $> Note: Defined as: mean((response - truth)^2)

性能指標を計算する
------------------

例として、勾配ブースティングマシンを`BostonHousing`データの一部に適用し、残りのデータから標準の性能指標である平均二乗偏差を計算してみよう。

``` r
n = getTaskSize(bh.task)
lrn = makeLearner("regr.gbm", n.trees = 1000)
mod = train(lrn, task = bh.task, subset = seq(1, n, 2))
pred = predict(mod, task = bh.task, subset = seq(2, n, 2))

performance(pred)
```

    $>      mse 
    $> 42.85008

他の指標の例として中央値二乗誤差(medse)を求めてみよう。

``` r
performance(pred, measures = medse)
```

    $>    medse 
    $> 8.930711

もちろん、独自に作成した指標も含めて、複数の指標を一度に計算することもできる。その場合、求めたい指標をリストにして渡す。

``` r
performance(pred, measures = list(mse, medse, mae))
```

    $>       mse     medse       mae 
    $> 42.850084  8.930711  4.547737

上記の方法は、学習問題や性能指標の種類が異なっても基本的には同じである。

指標計算に必要な情報
--------------------

一部の性能指標では、計算のために予測結果だけでなく、タスクやフィット済みモデルも必要とする。

一例は訓練にかかった時間(timetrain)だ。

``` r
performance(pred, measures = timetrain, model = mod)
```

    $> timetrain 
    $>     0.082

クラスター分析に関わる多くの性能指標はタスクを必要とする。

``` r
lrn = makeLearner("cluster.kmeans", centers = 3)
mod = train(lrn, mtcars.task)
pred = predict(mod, task = mtcars.task)

performance(pred, measures = dunn, task = mtcars.task)
```

    $>      dunn 
    $> 0.1462919

また、いくつかの指標は特定の種類の予測を必要とする。例えば2クラス分類におけるAUC(これはROC\[receiver operating characteristic\]曲線の下側の面積\[Area Under Curve\]である)を計算するためには、事後確率が必要である。ROC分析に関する詳細が必要であれば[ROC Analysis - mlr tutorial](https://mlr-org.github.io/mlr-tutorial/devel/html/roc_analysis/index.html)を確認してほしい。

``` r
lrn = makeLearner("classif.rpart", predict.type = "prob")
mod = train(lrn, task = sonar.task)
pred = predict(mod, task = sonar.task)

performance(pred, measures = auc)
```

    $>       auc 
    $> 0.9224018

また、分類問題に利用可能な性能指標(偽陽性率fprなど)の多くは、2クラス分類のみに利用可能であるという点に注意してもらいたい。

性能指標へのアクセス
--------------------

`mlr`における性能指標は`Measure`クラスのオブジェクトである。オブジェクトを通じて指標のプロパティ等には直接アクセスすることができる。各スロットに関する説明は[makeMeasure function | R Documentation](https://www.rdocumentation.org/packages/mlr/versions/2.10/topics/makeMeasure)を確認してもらいたい。

``` r
str(mmce)
```

    $> List of 10
    $>  $ id        : chr "mmce"
    $>  $ minimize  : logi TRUE
    $>  $ properties: chr [1:4] "classif" "classif.multi" "req.pred" "req.truth"
    $>  $ fun       :function (task, model, pred, feats, extra.args)  
    $>  $ extra.args: list()
    $>  $ best      : num 0
    $>  $ worst     : num 1
    $>  $ name      : chr "Mean misclassification error"
    $>  $ note      : chr "Defined as: mean(response != truth)"
    $>  $ aggr      :List of 4
    $>   ..$ id        : chr "test.mean"
    $>   ..$ name      : chr "Test mean"
    $>   ..$ fun       :function (task, perf.test, perf.train, measure, group, pred)  
    $>   ..$ properties: chr "req.test"
    $>   ..- attr(*, "class")= chr "Aggregation"
    $>  - attr(*, "class")= chr "Measure"

2クラス分類
-----------

### 性能と閾値の関係をプロットする

2クラス分類問題においては、予測された確率からクラスラベルへの割り当てを行う際の確率の閾値を設定できるということを思い出してもらいたい。`generateThreshVsPrefData`と`plotThreshVsPref`は、学習器のパフォーマンスと閾値の関係をプロットできる便利な関数だ。

パフォーマンスのプロットと閾値の自動チューニングに関して詳しい情報は[ROC Analysis - mlr tutorial](https://mlr-org.github.io/mlr-tutorial/devel/html/roc_analysis/index.html)を確認してほしい。

以下の例では、`Sonar`データセットを用い、偽陽性率(fpr)、偽陰性率(fnr)、平均誤分類率(mmce)を設定可能な範囲の閾値に対してプロットしている。

``` r
lrn = makeLearner("classif.lda", predict.type = "prob")
n = getTaskSize(sonar.task)
mod = train(lrn, task = sonar.task, subset = seq(1, n, by = 2))
pred = predict(mod, task = sonar.task, subset = seq(2, n, by = 2))

d = generateThreshVsPerfData(pred, measures = list(fpr, fnr, mmce))
plotThreshVsPerf(d)
```

![unnamed-chunk-12-1.png](https://qiita-image-store.s3.amazonaws.com/0/109348/d1957889-ee35-9061-7701-905de41f5fcf.png)

リサンプリング
==============

一般的に学習機の性能評価はリサンプリングを通じて行われる。リサンプリングの概要は次のようなものである。まず、データセット全体を*D*として、これを訓練セット*D*<sup>\**b*</sup>とテストセット*D* \\ *D*<sup>\**b*</sup>に分割する。この種の分割を*B*回行う(つまり、*b* = 1, ..., *B*とする)。そして、それぞれのテストセット、訓練セットの対を用いて訓練と予測を行い、性能指標*S*(*D*<sup>\**b*</sup>, *D* \\ *D*<sup>\**b*</sup>)を計算する。これにより*B*個の性能指標が得られるが、これを集約する(一般的には平均値が用いられる)。リサンプリングの方法には、クロスバリデーションやブートストラップなど様々な手法が存在する。

もしさらに詳しく知りたいのであれば、Simonによる論文([Resampling Strategies for Model Assessment and Selection | SpringerLink](https://link.springer.com/chapter/10.1007%2F978-0-387-47509-7_8))を読むのは悪い選択ではないだろう。また、Berndらによる論文、[Resampling methods for meta-model validation with recommendations for evolutionary computation](https://www.mitpressjournals.org/doi/pdf/10.1162/EVCO_a_00069)では、リサンプリング手法の統計的な背景に対して多くの説明がなされている。

リサンプリング手法を決める
--------------------------

`mlr`では`makeResampleDesc`関数を使ってリサンプリング手法を設定する。この関数にはリサンプリング手法の名前とともに、手法に応じてその他の情報(例えば繰り返し数など)を指定する。サポートしているサンプリング手法は以下のとおりである。

-   `CV`: クロスバリデーション(Cross-varidation)
-   `LOO`: 一つ抜き法(Leave-one-out cross-varidation)
-   `RepCV`: Repeatedクロスバリデーション(Repeated cross-varidation)
-   `Bootstrap`: out-of-bagブートストラップとそのバリエーション(b632等)
-   `Subsample`: サブサンプリング(モンテカルロクロスバリデーションとも呼ばれる)
-   `Holdout`: ホールドアウト法

3-fold(3分割)クロスバリデーションの場合は

``` r
rdesc = makeResampleDesc("CV", iters = 3)
rdesc
```

    $> Resample description: cross-validation with 3 iterations.
    $> Predict: test
    $> Stratification: FALSE

ホールドアウト法の場合は

``` r
rdesc = makeResampleDesc("Holdout")
rdesc
```

    $> Resample description: holdout with 0.67 split rate.
    $> Predict: test
    $> Stratification: FALSE

という具合だ。

これらのリサンプルdescriptionのうち、よく使うものは予め別名が用意してある。例えばホールドアウト法は`hout`、クロスバリデーションは`cv5`や`cv10`などよく使う分割数に対して定義してある。

``` r
hout
```

    $> Resample description: holdout with 0.67 split rate.
    $> Predict: test
    $> Stratification: FALSE

``` r
cv3
```

    $> Resample description: cross-validation with 3 iterations.
    $> Predict: test
    $> Stratification: FALSE

リサンプリングを実行する
------------------------

`resample`関数は指定されたリサンプリング手法により、学習機をタスク上で評価する。

最初の例として、`BostonHousing`データに対する線形回帰を3分割クロスバリデーションで評価してみよう。

*K*分割クロスバリデーションはデータセット*D*を*K*個の(ほぼ)等しいサイズのサブセットに分割する。*K*回の繰り返しの*b*番目では、*b*番目のサブセットがテストに、残りが訓練に使用される。

`resample`関数に学習器を指定する際には、`Learner`クラスのオブジェクトか学習器の名前(`regr.lm`など)のいずれを渡しても良い。性能指標は指定しなければ学習器に応じたデフォルトが使用される(回帰の場合は平均二乗誤差)。

``` r
rdesc = makeResampleDesc("CV", iters = 3)

r = resample("regr.lm", bh.task, rdesc)
```

    $> [Resample] cross-validation iter 1: mse.test.mean=25.1
    $> [Resample] cross-validation iter 2: mse.test.mean=23.1
    $> [Resample] cross-validation iter 3: mse.test.mean=21.9
    $> [Resample] Aggr. Result: mse.test.mean=23.4

``` r
r
```

    $> Resample Result
    $> Task: BostonHousing-example
    $> Learner: regr.lm
    $> Aggr perf: mse.test.mean=23.4
    $> Runtime: 0.047179

ここで`r`に格納したオブジェクトは`ResampleResult`クラスである。この中には評価結果の他に、実行時間や予測値、リサンプリング毎のフィット済みモデルなどが格納されている。

``` r
## 中身をざっと確認
names(r)
```

    $>  [1] "learner.id"     "task.id"        "task.desc"      "measures.train"
    $>  [5] "measures.test"  "aggr"           "pred"           "models"        
    $>  [9] "err.msgs"       "err.dumps"      "extract"        "runtime"

`r$measures.test`には各テストセットの性能指標が入っている。

``` r
## 各テストセットの性能指標
r$measures.test
```

    $>   iter      mse
    $> 1    1 25.13717
    $> 2    2 23.12795
    $> 3    3 21.91527

`r$aggr`には集約(aggrigate)後の性能指標が入っている。

``` r
## 集約後の性能指標
r$aggr
```

    $> mse.test.mean 
    $>      23.39346

名前`mse.test.mean`は、性能指標が`mse`であり、`test.mean`によりデータが集約されていることを表している。`test.mean`は多くの性能指標においてデフォルトの集約方法であり、その名前が示すようにテストデータの性能指標の平均値である。

`mlr`ではどのような種類の学習器も同じようにリサンプリングを行える。以下では、分類問題の例として`Sonar`データセットに対する分類木を5反復のサブサンプリングで評価してみよう。

サブサンプリングの各繰り返しでは、データセット*D*はランダムに訓練データとテストデータに分割される。このとき、テストデータには指定の割合のデータ数が割り当てられる。この反復が1の場合はホールドアウト法と同じである。

評価したい性能指標はリストとしてまとめて指定することもできる。以下の例では平均誤分類、偽陽性・偽陰性率、訓練時間を指定している。

``` r
rdesc = makeResampleDesc("Subsample", iter = 5, split = 4/5)
lrn = makeLearner("classif.rpart", parms = list(split = "information"))
r = resample(lrn, sonar.task, rdesc, measures = list(mmce, fpr, fnr, timetrain))
```

    $> [Resample] subsampling iter 1: mmce.test.mean=0.405,fpr.test.mean= 0.5,fnr.test.mean=0.318,timetrain.test.mean=0.019
    $> [Resample] subsampling iter 2: mmce.test.mean=0.262,fpr.test.mean=0.381,fnr.test.mean=0.143,timetrain.test.mean=0.016
    $> [Resample] subsampling iter 3: mmce.test.mean=0.19,fpr.test.mean=0.304,fnr.test.mean=0.0526,timetrain.test.mean=0.015
    $> [Resample] subsampling iter 4: mmce.test.mean=0.429,fpr.test.mean=0.35,fnr.test.mean= 0.5,timetrain.test.mean=0.013
    $> [Resample] subsampling iter 5: mmce.test.mean=0.333,fpr.test.mean=0.235,fnr.test.mean= 0.4,timetrain.test.mean=0.036
    $> [Resample] Aggr. Result: mmce.test.mean=0.324,fpr.test.mean=0.354,fnr.test.mean=0.283,timetrain.test.mean=0.0198

``` r
r
```

    $> Resample Result
    $> Task: Sonar-example
    $> Learner: classif.rpart
    $> Aggr perf: mmce.test.mean=0.324,fpr.test.mean=0.354,fnr.test.mean=0.283,timetrain.test.mean=0.0198
    $> Runtime: 0.160807

もし指標を後から追加したくなったら、`addRRMeasure`関数を使うと良い。

``` r
addRRMeasure(r, list(ber, timepredict))
```

    $> Resample Result
    $> Task: Sonar-example
    $> Learner: classif.rpart
    $> Aggr perf: mmce.test.mean=0.324,fpr.test.mean=0.354,fnr.test.mean=0.283,timetrain.test.mean=0.0198,ber.test.mean=0.318,timepredict.test.mean=0.005
    $> Runtime: 0.160807

デフォルトでは`resample`関数は進捗と途中結果を表示するが、`show.info=FALSE`で非表示にもできる。このようなメッセージを完全に制御したかったら、[Configuration - mlr tutorial](https://mlr-org.github.io/mlr-tutorial/devel/html/configureMlr/index.html)を確認してもらいたい。

上記例では学習器を明示的に作成してから`resample`に渡したが、代わりに学習器の名前を指定しても良い。その場合、学習器のパラメータは`...`引数を通じて渡すことができる。

``` r
resample("classif.rpart", parms = list(split = "information"), sonar.task, rdesc,
         measures = list(mmce, fpr, fnr, timetrain), show.info = FALSE)
```

    $> Resample Result
    $> Task: Sonar-example
    $> Learner: classif.rpart
    $> Aggr perf: mmce.test.mean=0.267,fpr.test.mean=0.246,fnr.test.mean=0.282,timetrain.test.mean=0.0146
    $> Runtime: 0.258963

リサンプル結果へのアクセス
--------------------------

学習器の性能以外にも、リサンプル結果から様々な情報を得ることが出来る。例えばリサンプリングの各繰り返しに対応する予測値やフィット済みモデル等だ。以下で情報の取得の仕方をみていこう。

### 予測値

デフォルトでは、`ResampleResult`はリサンプリングで得た予測値を含んでいる。メモリ節約などの目的でこれを止めさせたければ、`resample`関数に`keep.pred = FALSE`を指定する。

予測値は`$pred`スロットに格納されている。また、`getRRPredictions`関数を使ってアクセスすることもできる。

``` r
r$pred
```

    $> Resampled Prediction for:
    $> Resample description: subsampling with 5 iterations and 0.80 split rate.
    $> Predict: test
    $> Stratification: FALSE
    $> predict.type: response
    $> threshold: 
    $> time (mean): 0.00
    $>    id truth response iter  set
    $> 1 161     M        M    1 test
    $> 2 151     M        R    1 test
    $> 3 182     M        M    1 test
    $> 4 114     M        R    1 test
    $> 5  76     R        R    1 test
    $> 6 122     M        M    1 test
    $> ... (210 rows, 5 cols)

``` r
pred = getRRPredictions(r)
pred
```

    $> Resampled Prediction for:
    $> Resample description: subsampling with 5 iterations and 0.80 split rate.
    $> Predict: test
    $> Stratification: FALSE
    $> predict.type: response
    $> threshold: 
    $> time (mean): 0.00
    $>    id truth response iter  set
    $> 1 161     M        M    1 test
    $> 2 151     M        R    1 test
    $> 3 182     M        M    1 test
    $> 4 114     M        R    1 test
    $> 5  76     R        R    1 test
    $> 6 122     M        M    1 test
    $> ... (210 rows, 5 cols)

ここで作成した`pred`は`ResamplePrediction`クラスのオブジェクトである。これは`Prediction`オブジェクトのように`$data`にデータフレームとして予測値と真値(教師あり学習の場合)が格納されている。`as.data.frame`を使って直接`$data`スロットの中身を取得できる。さらに、`Prediction`オブジェクトに対するゲッター関数は全て利用可能である。

``` r
head(as.data.frame(pred))
```

    $>    id truth response iter  set
    $> 1 161     M        M    1 test
    $> 2 151     M        R    1 test
    $> 3 182     M        M    1 test
    $> 4 114     M        R    1 test
    $> 5  76     R        R    1 test
    $> 6 122     M        M    1 test

``` r
head(getPredictionTruth(pred))
```

    $> [1] M M M M R M
    $> Levels: M R

データフレームの`iter`と`set`は繰り返し回数とデータセットの種類(訓練なのかテストなのか)を示している。

デフォルトでは予測はテストセットだけに行われるが、`makeResampleDesc`に対し、`predict = "train"`を指定で訓練セットだけに、`predict = "both"`を指定で訓練セットとテストセットの両方に予測を行うことが出来る。後で例を見るが、*b632*や*b632+*のようなブートストラップ手法ではこれらの設定が必要となる。

以下では単純なホールドアウト法の例を見よう。つまり、テストセットと訓練セットへの分割は一度だけ行い、予測は両方のデータセットを用いて行う。

``` r
rdesc = makeResampleDesc("Holdout", predict = "both")

r = resample("classif.lda", iris.task, rdesc, show.info = FALSE)
r
```

    $> Resample Result
    $> Task: iris-example
    $> Learner: classif.lda
    $> Aggr perf: mmce.test.mean=0.02
    $> Runtime: 0.0246351

``` r
r$aggr
```

    $> mmce.test.mean 
    $>           0.02

(`predict="both"`の指定にかかわらず、`r$aggr`ではテストデータに対するmmceしか計算しないことに注意してもらいたい。訓練セットに対して計算する方法はこの後で説明する。)

リサンプリング結果から予測を取り出す方法として、`getRRPredictionList`を使う方法もある。これは、分割されたデータセット(訓練/テスト)それぞれと、リサンプリングの繰り返し毎に分割した単位でまとめた予測結果のリストを返す。

``` r
getRRPredictionList(r)
```

    $> $train
    $> $train$`1`
    $> Prediction: 100 observations
    $> predict.type: response
    $> threshold: 
    $> time: 0.00
    $>      id      truth   response
    $> 85   85 versicolor versicolor
    $> 13   13     setosa     setosa
    $> 140 140  virginica  virginica
    $> 109 109  virginica  virginica
    $> 70   70 versicolor versicolor
    $> 27   27     setosa     setosa
    $> ... (100 rows, 3 cols)
    $> 
    $> 
    $> 
    $> $test
    $> $test$`1`
    $> Prediction: 50 observations
    $> predict.type: response
    $> threshold: 
    $> time: 0.00
    $>      id      truth   response
    $> 82   82 versicolor versicolor
    $> 55   55 versicolor versicolor
    $> 29   29     setosa     setosa
    $> 147 147  virginica  virginica
    $> 44   44     setosa     setosa
    $> 83   83 versicolor versicolor
    $> ... (50 rows, 3 cols)

### 訓練済みモデルの抽出

リサンプリング毎に学習器は訓練セットにフィットさせられる。標準では、`WrappedModel`は`ResampleResult`オブジェクトには含まれておらず、`$models`スロットは空だ。これを保持するためには、`resample`関数を呼び出す際に引数`models = TRUE`を指定する必要がある。以下に生存時間分析の例を見よう。

``` r
## 3分割クロスバリデーション
rdesc = makeResampleDesc("CV", iters = 3)

r = resample("surv.coxph", lung.task, rdesc, show.info = FALSE, models = TRUE)
r$models
```

    $> [[1]]
    $> Model for learner.id=surv.coxph; learner.class=surv.coxph
    $> Trained on: task.id = lung-example; obs = 111; features = 8
    $> Hyperparameters: 
    $> 
    $> [[2]]
    $> Model for learner.id=surv.coxph; learner.class=surv.coxph
    $> Trained on: task.id = lung-example; obs = 111; features = 8
    $> Hyperparameters: 
    $> 
    $> [[3]]
    $> Model for learner.id=surv.coxph; learner.class=surv.coxph
    $> Trained on: task.id = lung-example; obs = 112; features = 8
    $> Hyperparameters:

### 他の抽出方法

完全なフィット済みモデルを保持しようとすると、リサンプリングの繰り返し数が多かったりオブジェクトが大きかったりする場合にメモリの消費量が大きくなってしまう。モデルの全ての情報を保持する代わりに、`resample`関数の`extract`引数に指定することで必要な情報だけを保持することができる。引数`extract`に対しては、リサンプリング毎の各`WrapedModel`オブジェクトに適用するための関数を渡す必要がある。

以下では、`mtcars`データセットをk=3のk-meansでクラスタリングし、クラスター中心だけを保持する例を紹介する。

``` r
rdesc = makeResampleDesc("CV", iter = 3)

r = resample("cluster.kmeans", mtcars.task, rdesc, show.info = FALSE,
             centers = 3, extract = function(x){getLearnerModel(x)$centers})
```

    $> 
    $> This is package 'modeest' written by P. PONCET.
    $> For a complete list of functions, use 'library(help = "modeest")' or 'help.start()'.

``` r
r$extract
```

    $> [[1]]
    $>        mpg      cyl      disp    hp     drat       wt     qsec        vs
    $> 1 26.96667 4.000000  99.08333  89.5 4.076667 2.087167 18.26833 0.8333333
    $> 2 20.61429 5.428571 166.81429 104.0 3.715714 3.167857 19.11429 0.7142857
    $> 3 15.26667 8.000000 356.64444 216.0 3.251111 3.956556 16.55556 0.0000000
    $>          am     gear     carb
    $> 1 1.0000000 4.333333 1.500000
    $> 2 0.2857143 3.857143 3.000000
    $> 3 0.2222222 3.444444 3.666667
    $> 
    $> [[2]]
    $>        mpg      cyl     disp       hp   drat       wt     qsec        vs
    $> 1 15.03750 8.000000 351.1750 205.0000 3.1825 4.128125 17.08375 0.0000000
    $> 2 21.22222 5.111111 157.6889 110.1111 3.7600 2.945000 18.75889 0.6666667
    $> 3 31.00000 4.000000  76.1250  62.2500 4.3275 1.896250 19.19750 1.0000000
    $>          am     gear     carb
    $> 1 0.0000000 3.000000 3.375000
    $> 2 0.4444444 3.888889 2.888889
    $> 3 1.0000000 4.000000 1.250000
    $> 
    $> [[3]]
    $>        mpg      cyl     disp       hp     drat       wt     qsec        vs
    $> 1 14.68571 8.000000 384.8571 230.5714 3.378571 4.077000 16.34000 0.0000000
    $> 2 25.30000 4.500000 113.8125 101.2500 4.037500 2.307875 18.00125 0.7500000
    $> 3 16.96667 7.333333 276.1000 145.8333 2.981667 3.580000 18.20500 0.3333333
    $>          am     gear  carb
    $> 1 0.2857143 3.571429 4.000
    $> 2 0.7500000 4.250000 2.375
    $> 3 0.0000000 3.000000 2.000

他の例として、フィット済みの回帰木から変数の重要度を`getFeatureImportance`を使って抽出してみよう(より詳しい内容は[Feature Selection - mlr tutorial](https://mlr-org.github.io/mlr-tutorial/devel/html/feature_selection/index.html)を確認してもらいたい)。

``` r
r = resample("regr.rpart", bh.task, rdesc, show.info = FALSE, extract = getFeatureImportance)
r$extract
```

    $> [[1]]
    $> FeatureImportance:
    $> Task: BostonHousing-example
    $> 
    $> Learner: regr.rpart
    $> Measure: NA
    $> Contrast: NA
    $> Aggregation: function (x)  x
    $> Replace: NA
    $> Number of Monte-Carlo iterations: NA
    $> Local: FALSE
    $>      crim       zn    indus chas      nox       rm      age      dis rad
    $> 1 2689.63 867.4877 3719.711    0 2103.622 16096.32 2574.183 3647.211   0
    $>        tax  ptratio       b    lstat
    $> 1 1972.207 3712.621 395.486 8608.757
    $> 
    $> [[2]]
    $> FeatureImportance:
    $> Task: BostonHousing-example
    $> 
    $> Learner: regr.rpart
    $> Measure: NA
    $> Contrast: NA
    $> Aggregation: function (x)  x
    $> Replace: NA
    $> Number of Monte-Carlo iterations: NA
    $> Local: FALSE
    $>       crim       zn  indus chas      nox       rm      age     dis
    $> 1 7491.707 5423.593 7295.2    0 7348.742 14014.78 1391.373 2309.92
    $>        rad      tax  ptratio b    lstat
    $> 1 340.3975 1871.451 938.0743 0 17618.49
    $> 
    $> [[3]]
    $> FeatureImportance:
    $> Task: BostonHousing-example
    $> 
    $> Learner: regr.rpart
    $> Measure: NA
    $> Contrast: NA
    $> Aggregation: function (x)  x
    $> Replace: NA
    $> Number of Monte-Carlo iterations: NA
    $> Local: FALSE
    $>       crim       zn    indus chas     nox       rm     age      dis
    $> 1 2532.084 4637.312 6150.854    0 6015.01 11330.75 6843.29 2049.772
    $>        rad      tax  ptratio        b    lstat
    $> 1 525.3815 747.8954 1925.899 62.58285 17336.77

階層化とブロック化
------------------

-   カテゴリー変数に対する階層化とは、訓練セットとテストセット内で各値の比率が変わらないようにすることを指す。階層化が可能なのは目的変数がカテゴリーである場合(教師あり学習における分類や生存時間分析)や、説明変数がカテゴリーである場合に限られる。
-   ブロック化とは、観測値の一部分をブロックとして扱い、リサンプリングの間にブロックが分割されないように扱うことを指す。つまり、ブロック全体は訓練セットかテストセットのいずれかにまとまって属すことになる。

### 目的変数の階層化

分類においては、元のデータと同じ比率で各クラスの値が含まれていることが望ましい。これはクラス間の観測数が不均衡であったり、データセットの大きさが小さい場合に有効である。さもなければ、観測数が少ないクラスのデータが訓練セットに含まれないということが起こりうる。これは分類性能の低下やモデルのクラッシュにつながる。階層化リサンプリングを行うためには、`makeResampleDesc`実行時に`stratify = TRUE`を指定する。

``` r
rdesc = makeResampleDesc("CV", iters = 3, stratify = TRUE)

r = resample("classif.lda", iris.task, rdesc, show.info = FALSE)
r
```

    $> Resample Result
    $> Task: iris-example
    $> Learner: classif.lda
    $> Aggr perf: mmce.test.mean=0.02
    $> Runtime: 0.027998

階層化を生存時間分析に対して行う場合は、打ち切りの割合が制御される。

### 説明変数の階層化

説明変数の階層化が必要な場合もある。この場合は、`stratify.cols`引数に対して階層化したい因子型変数を指定する。

``` r
rdesc = makeResampleDesc("CV", iter = 3, stratify.cols = "chas")

r = resample("regr.rpart", bh.task, rdesc, show.info = FALSE)
r
```

    $> Resample Result
    $> Task: BostonHousing-example
    $> Learner: regr.rpart
    $> Aggr perf: mse.test.mean=23.7
    $> Runtime: 0.0576711

### ブロック化

いくつかの観測値が互いに関連しており、これらが訓練データとテストデータに分割されるのが望ましくない場合には、タスク作成時にその情報を`blocking`引数に因子型ベクトルを与えることで指定する。

``` r
## それぞれ30の観測値からなる5つのブロックを指定する例
task = makeClassifTask(data = iris, target = "Species", blocking = factor(rep(1:5, each = 30)))
task
```

    $> Supervised task: iris
    $> Type: classif
    $> Target: Species
    $> Observations: 150
    $> Features:
    $> numerics  factors  ordered 
    $>        4        0        0 
    $> Missings: FALSE
    $> Has weights: FALSE
    $> Has blocking: TRUE
    $> Classes: 3
    $>     setosa versicolor  virginica 
    $>         50         50         50 
    $> Positive class: NA

リサンプリングの詳細とリサンプルのインスタンス
----------------------------------------------

既に説明したように、リサンプリング手法は`makeResampleDesc`関数を使って指定する。

``` r
rdesc = makeResampleDesc("CV", iter = 3)
rdesc
```

    $> Resample description: cross-validation with 3 iterations.
    $> Predict: test
    $> Stratification: FALSE

``` r
str(rdesc)
```

    $> List of 4
    $>  $ id      : chr "cross-validation"
    $>  $ iters   : int 3
    $>  $ predict : chr "test"
    $>  $ stratify: logi FALSE
    $>  - attr(*, "class")= chr [1:2] "CVDesc" "ResampleDesc"

上記`rdesc`は`ResampleDesc`クラス(resample descriptionの略)を継承しており、原則として、リサンプリング手法に関する必要な情報(繰り返し数、訓練セットとテストセットの比率、階層化したい変数など)を全て含んでいる。

`makeResampleInstance`関数は、データセットに含まれるデータ数を直接指定するか、タスクを指定することで、`ResampleDesc`に従って訓練セットとテストセットの概要を生成する。

``` r
## taskに基づくリサンプルインスタンスの生成
rin = makeResampleInstance(rdesc, iris.task)
rin
```

    $> Resample instance for 150 cases.
    $> Resample description: cross-validation with 3 iterations.
    $> Predict: test
    $> Stratification: FALSE

``` r
str(rin)
```

    $> List of 5
    $>  $ desc      :List of 4
    $>   ..$ id      : chr "cross-validation"
    $>   ..$ iters   : int 3
    $>   ..$ predict : chr "test"
    $>   ..$ stratify: logi FALSE
    $>   ..- attr(*, "class")= chr [1:2] "CVDesc" "ResampleDesc"
    $>  $ size      : int 150
    $>  $ train.inds:List of 3
    $>   ..$ : int [1:100] 11 111 2 125 49 71 82 16 12 121 ...
    $>   ..$ : int [1:100] 18 11 111 20 2 125 16 121 70 68 ...
    $>   ..$ : int [1:100] 18 20 49 71 82 12 68 102 5 25 ...
    $>  $ test.inds :List of 3
    $>   ..$ : int [1:50] 1 5 8 10 14 15 18 20 23 32 ...
    $>   ..$ : int [1:50] 3 4 7 12 17 22 25 27 29 31 ...
    $>   ..$ : int [1:50] 2 6 9 11 13 16 19 21 24 26 ...
    $>  $ group     : Factor w/ 0 levels: 
    $>  - attr(*, "class")= chr "ResampleInstance"

``` r
## データセットのサイズを指定してリサンプルインスタンスを生成する例
rin = makeResampleInstance(rdesc, size = nrow(iris))
str(rin)
```

    $> List of 5
    $>  $ desc      :List of 4
    $>   ..$ id      : chr "cross-validation"
    $>   ..$ iters   : int 3
    $>   ..$ predict : chr "test"
    $>   ..$ stratify: logi FALSE
    $>   ..- attr(*, "class")= chr [1:2] "CVDesc" "ResampleDesc"
    $>  $ size      : int 150
    $>  $ train.inds:List of 3
    $>   ..$ : int [1:100] 23 15 5 81 143 38 102 145 85 132 ...
    $>   ..$ : int [1:100] 99 23 78 41 15 81 108 128 102 145 ...
    $>   ..$ : int [1:100] 99 78 41 5 108 128 143 38 132 84 ...
    $>  $ test.inds :List of 3
    $>   ..$ : int [1:50] 1 3 7 8 9 10 11 13 16 19 ...
    $>   ..$ : int [1:50] 2 5 17 21 24 30 34 36 38 39 ...
    $>   ..$ : int [1:50] 4 6 12 14 15 18 20 22 23 28 ...
    $>  $ group     : Factor w/ 0 levels: 
    $>  - attr(*, "class")= chr "ResampleInstance"

ここで`rin`は`ResampleInstance`クラスを継承しており、訓練セットとテストセットのインデックスをリストとして含んでいる。

`ResampleDesc`が`resample`に渡されると、インスタンスの生成は内部的に行われる。もちろん、`ResampleInstance`を直接渡すこともできる。

リサンプルの詳細(resample description)とリサンプルのインスタンス、そしてリサンプル関数と分割するのは、複雑にしすぎているのではと感じるかもしれないが、幾つかの利点がある。

-   リサンプルインスタンスを用いると、同じ訓練セットとテストセットを用いて学習器の性能比較を行うことが容易になる。これは、既に実施した性能比較試験に対し、他の手法を追加したい場合などに特に便利である。また、後で結果を再現するためにデータとリサンプルインスタンスをセットで保管しておくこともできる。

``` r
rdesc = makeResampleDesc("CV", iter = 3)
rin = makeResampleInstance(rdesc, task = iris.task)

## 同じインスタンスを使い、2種類の学習器で性能指標を計算する
r.lda = resample("classif.lda", iris.task, rin, show.info = FALSE)
r.rpart = resample("classif.rpart", iris.task, rin, show.info = FALSE)
c("lda" = r.lda$aggr, "rpart" = r.rpart$aggr)
```

    $>   lda.mmce.test.mean rpart.mmce.test.mean 
    $>                 0.02                 0.06

-   新しいリサンプリング手法を追加したければ、`ResampleDesc`および`ResampleInstance`クラスのインスタンスを作成すればよく、`resample`関数やそれ以上のメソッドに触る必要はない。

通常、`makeResampleInstance`を呼び出したときの訓練セットとテストセットのインデックスはランダムに割り当てられる。主にホールドアウト法においては、これを完全にマニュアルで行わなければならない場面がある。これは`makeFixedHoldoutInstance`関数を使うと実現できる。

``` r
rin = makeFixedHoldoutInstance(train.inds = 1:100, test.inds = 101:150, size = 150)
rin
```

    $> Resample instance for 150 cases.
    $> Resample description: holdout with 0.67 split rate.
    $> Predict: test
    $> Stratification: FALSE

性能指標の集約
--------------

リサンプリングそれぞれに対して性能指標を計算したら、それを集計する必要がある。

大半のリサンプリング手法(ホールドアウト法、クロスバリデーション、サブサンプリングなど)では、性能指標はテストデータのみで計算され、平均によって集約される。

`mlr`における性能指標を表現する`Measure`クラスのオブジェクトは、`$aggr`スロットに対応するデフォルトの集約手法を格納している。大半は`test.mean`である。例外の一つは平均二乗誤差平方根(rmse)である。

``` r
## 一般的な集約手法
mmce$aggr
```

    $> Aggregation function: test.mean

``` r
## 具体的な計算方法
mmce$aggr$fun
```

    $> function (task, perf.test, perf.train, measure, group, pred) 
    $> mean(perf.test)
    $> <bytecode: 0x7fd6dfea4428>
    $> <environment: namespace:mlr>

``` r
## rmseの場合
rmse$aggr
```

    $> Aggregation function: test.rmse

``` r
## test.rmseの具体的な計算方法
rmse$aggr$fun
```

    $> function (task, perf.test, perf.train, measure, group, pred) 
    $> sqrt(mean(perf.test^2))
    $> <bytecode: 0x7fd6e28f8978>
    $> <environment: namespace:mlr>

`setAggrigation`関数を使うと、集約方法を変更することも出来る。利用可能な集約手法の一覧は[aggregations function | R Documentation](https://www.rdocumentation.org/packages/mlr/versions/2.10/topics/aggregations)を確認してほしい。

### 例: 一つの指標に複数の集約方法

`test.median`、`test.min`、`test.max`はそれぞれテストセットから求めた性能指標を中央値、最小値、最大値で集約する。

``` r
mseTestMedian = setAggregation(mse, test.median)
mseTestMin = setAggregation(mse, test.min)
mseTestMax = setAggregation(mse, test.max)
rdesc = makeResampleDesc("CV", iter = 3)
r = resample("regr.lm", bh.task, rdesc, show.info = FALSE, 
             measures = list(mse, mseTestMedian, mseTestMin, mseTestMax))
r
```

    $> Resample Result
    $> Task: BostonHousing-example
    $> Learner: regr.lm
    $> Aggr perf: mse.test.mean=24.2,mse.test.median=23.9,mse.test.min=20.8,mse.test.max=27.8
    $> Runtime: 0.0312719

``` r
r$aggr
```

    $>   mse.test.mean mse.test.median    mse.test.min    mse.test.max 
    $>        24.17104        23.92198        20.82022        27.77090

### 例: 訓練セットの誤差を計算する

平均誤分類率を訓練セットとテストセットに対して計算する例を示す。`makeResampleDesc`実行時に`predict = "both"`を指定しておく必要があることに注意してもらいたい。

``` r
mmceTrainMean = setAggregation(mmce, train.mean)
rdesc = makeResampleDesc("CV", iters = 3, predict = "both")
r = resample("classif.rpart", iris.task, rdesc, measures = list(mmce, mmceTrainMean))
```

    $> [Resample] cross-validation iter 1: mmce.train.mean=0.03,mmce.test.mean=0.08
    $> [Resample] cross-validation iter 2: mmce.train.mean=0.05,mmce.test.mean=0.04
    $> [Resample] cross-validation iter 3: mmce.train.mean=0.01,mmce.test.mean= 0.1
    $> [Resample] Aggr. Result: mmce.test.mean=0.0733,mmce.train.mean=0.03

### 例: ブートストラップ

out-of-bagブートストラップ推定では、まず元のデータセット*D*から重複ありの抽出によって*D*<sup>\*1</sup>, ..., *D*<sup>\**B*</sup>と*B*個の新しいデータセット(要素数は元のデータセットと同じ)を作成する。そして、*b*回目の繰り返しでは、*D*<sup>\**b*</sup>を訓練セットに使い、使われなかった要素*D* \\ *D*<sup>\**b*</sup>をテストセットに用いて各繰り返しに対する推定値を計算し、最終的に*B*個の推定値を得る。

out-of-bagブートストラップの変種である*b632*と*b632+*では、訓練セットのパフォーマンスとOOBサンプルのパフォーマンスの凸結合を計算するため、訓練セットに対する予測と適切な集計方法を必要とする。

``` r
## ブートストラップをリサンプリング手法に選び、予測は訓練セットとテストセットの両方に行う
rdesc = makeResampleDesc("Bootstrap", predict = "both", iters = 10)

## b632およびb632+専用の集計手法を設定する
mmceB632 = setAggregation(mmce, b632)
mmceB632plus = setAggregation(mmce, b632plus)

r = resample("classif.rpart", iris.task, rdesc, measures = list(mmce, mmceB632, mmceB632plus),
             show.info = FALSE)
r$measures.train
```

    $>    iter        mmce        mmce        mmce
    $> 1     1 0.006666667 0.006666667 0.006666667
    $> 2     2 0.033333333 0.033333333 0.033333333
    $> 3     3 0.026666667 0.026666667 0.026666667
    $> 4     4 0.026666667 0.026666667 0.026666667
    $> 5     5 0.020000000 0.020000000 0.020000000
    $> 6     6 0.026666667 0.026666667 0.026666667
    $> 7     7 0.033333333 0.033333333 0.033333333
    $> 8     8 0.026666667 0.026666667 0.026666667
    $> 9     9 0.026666667 0.026666667 0.026666667
    $> 10   10 0.006666667 0.006666667 0.006666667

``` r
r$aggr
```

    $> mmce.test.mean      mmce.b632  mmce.b632plus 
    $>     0.04612359     0.03773677     0.03841055

便利な関数
----------

これまでに説明した方法は柔軟ではあるが、学習器を少し試してみたい場合にはタイプ数が多くて面倒だ。`mlr`には様々な略記法が用意してあるが、リサンプリング手法についても同様にショートカットが用意されている。ホールドアウトやクロスバリデーション、ブートストラップ(b632)等のよく使うリサンプリング手法にはそれぞれ特有の関数が用意してある。

``` r
crossval("classif.lda", iris.task, iters = 3, measures = list(mmce, ber))
```

    $> [Resample] cross-validation iter 1: mmce.test.mean=0.04,ber.test.mean=0.0303
    $> [Resample] cross-validation iter 2: mmce.test.mean=0.02,ber.test.mean=0.0167
    $> [Resample] cross-validation iter 3: mmce.test.mean=   0,ber.test.mean=   0
    $> [Resample] Aggr. Result: mmce.test.mean=0.02,ber.test.mean=0.0157

    $> Resample Result
    $> Task: iris-example
    $> Learner: classif.lda
    $> Aggr perf: mmce.test.mean=0.02,ber.test.mean=0.0157
    $> Runtime: 0.0347431

``` r
bootstrapB632plus("regr.lm", bh.task, iters = 3, measures = list(mse, mae))
```

    $> [Resample] OOB bootstrapping iter 1: mse.b632plus=21.1,mae.b632plus=3.15,mse.b632plus=23.7,mae.b632plus=3.59
    $> [Resample] OOB bootstrapping iter 2: mse.b632plus=18.4,mae.b632plus=3.04,mse.b632plus=28.8,mae.b632plus=3.98
    $> [Resample] OOB bootstrapping iter 3: mse.b632plus=23.1,mae.b632plus=3.35,mse.b632plus=16.1,mae.b632plus=2.99
    $> [Resample] Aggr. Result: mse.b632plus=22.2,mae.b632plus=3.41

    $> Resample Result
    $> Task: BostonHousing-example
    $> Learner: regr.lm
    $> Aggr perf: mse.b632plus=22.2,mae.b632plus=3.41
    $> Runtime: 0.0514081

ハイパーパラメータのチューニング
================================

多くの機械学習アルゴリズムはハイパーパラメータを持っている。学習器のチュートリアルでも説明したが、ハイパーパラメータとして特定の値を設定したければその値を`makeLearner`に渡すだけで良い。しかし、ハイパーパラメータの最適な値というのは大抵の場合は自明ではなく、できれば自動的に調整する手法が欲しい。

機械学習アルゴリズムをチューニングするためには、以下の点を指定する必要がある。

-   パラメータの探索範囲
-   最適化アルゴリズム(チューニングメソッドとも呼ぶ)
-   評価手法(すなわち、リサンプリング手法と性能指標)

パラメータの探索範囲: 例としてサポートベクターマシン(SVM)におけるパラメータ`C`の探索範囲を指定してみよう。

``` r
ps = makeParamSet(
  makeNumericParam("C", lower = 0.01, upper = 0.1)
)
```

最適化アルゴリズム: 例としてランダムサーチを指定してみよう。

``` r
ctrl = makeTuneControlRandom(maxit = 100L)
```

評価手法: リサンプリング手法として3分割クロスバリデーションを、性能指標として精度を指定してみよう。

``` r
rdesc = makeResampleDesc("CV", iter = 3L)
measure = acc
```

評価手法の指定方法については既に説明したところであるので、ここから先は探索範囲と最適化アルゴリズムの指定方法と、チューニングをどのように行い、結果にどのようにアクセスし、さらにチューニング結果を可視化する方法について幾つかの例を通して説明していこう。

このセクションを通して、例としては分類問題を取り上げるが、他の学習問題についても同様の手順で実行できるはずだ。

このさき、irisの分類タスクを使用して、SVMのハイパーパラメータを放射基底関数(RBF)カーネルを使ってチューニングする例を説明する。以下の例では、コストパラメータ`C`と、RBFカーネルのパラメータ`sigma`をチューニングする。

パラメータ探索空間の指定
------------------------

チューニングの際にまず指定しなければならないのは値の探索範囲である。これは例えば"いくつかの値の中のどれか"かもしれないし、"10<sup>−10</sup>から10<sup>10</sup>までの間の中のどこか"かもしれない。

探索空間の指定に際して、パラメータの探索範囲についての情報を含む`ParamSet`オブジェクトを作成する。これには`makeParamSet`関数を用いる。

例として、パメータ`C`と`sigma`の探索範囲を両方共0.5, 1.0, 1.5, 2.0という離散値に設定する例を見よう。それぞれのパラメータにどのような名前が使われているのかは、`kernlab`パッケージで定義されている(cf. [kernlab package | R Documentation](https://www.rdocumentation.org/packages/kernlab/versions/0.9-25))。

``` r
discrete_ps = makeParamSet(
  makeDiscreteParam("C", values = c(0.5, 1.0, 1.5, 2.0)),
  makeDiscreteParam("sigma", values = c(0.5, 1.0, 1.5, 2.0))
)
discrete_ps
```

    $>           Type len Def      Constr Req Tunable Trafo
    $> C     discrete   -   - 0.5,1,1.5,2   -    TRUE     -
    $> sigma discrete   -   - 0.5,1,1.5,2   -    TRUE     -

連続値の探索範囲を指定する際には`makeDiscreteParam`の代わりに`makeNumericParam`を使用する。また、探索範囲として10<sup>−10</sup>から10<sup>10</sup>のような範囲を指定する際には、`trafo`引数に変換用の関数を指定できる(trafoはtransformationの略)。変換用の関数を指定した場合、変換前のスケールで行われ、学習アルゴリズムに値を渡す前に変換が行われる。

``` r
num_ps = makeParamSet(
  makeNumericParam("C", lower = -10, upper = 10, trafo = function(x) 10^x),
  makeNumericParam("sigma", lower = -10, upper = 10, trafo = function(x) 10^x)
)
```

他にも数多くのパラメータが利用できるが、詳しくは[makeParamSet function | R Documentation](https://www.rdocumentation.org/packages/ParamHelpers/versions/1.10/topics/makeParamSet)を確認してもらいたい。

パラメータをリストの形で指定しなければならない関数もあるが、`mlr`を通じてその関数を扱う場合、`mlr`はできるかぎりリスト構造を除去し、パラメータを直接指定できるように試みる。例えばSVMを実行する関数の`ksvm`は、`kpar`引数に`sigma`のようなカーネルパラメータをリストで渡す必要がある。今例を見たように、`mlr`は`sigma`を直接扱うことができる。この仕組みのおかげで、`mlr`は様々なパッケージの学習器を統一したインターフェースで扱うことができるのだ。

最適化アルゴリズムの指定
------------------------

パラメータの探索範囲を決めたら次は最適化アルゴリズムを指定する。`mlr`において最適化アルゴリズムは`TuneControl`クラスのオブジェクトとして扱われる。

グリッドサーチは適切なパラメータを見つけるための標準的な(しかし遅い)方法の一つだ。

先に例を挙げた`discrete_ps`の場合、グリッドサーチは単純に値の全ての組合せを探索する。

``` r
ctrl = makeTuneControlGrid()
```

`num_ps`の場合は、グリッドサーチは探索範囲をまず均等なサイズのステップに分割する。標準では分割数は10だが、これは`resolution`引数で変更できる。ここでは`resolution`に15を指定してみよう。なお、ここで言う均等な15分割というのは、`10^seq(-10, 10, length.out = 15)`という意味である。

``` r
ctrl = makeTuneControlGrid(resolution = 15L)
```

クロスバリデーション以外にも多くの最適化アルゴリズムが利用可能であるが、詳しくは[TuneControl function | R Documentation](https://www.rdocumentation.org/packages/mlr/versions/2.10/topics/TuneControl)を確認してもらいたい。

グリッドサーチは一般的には遅すぎるので、ランダムサーチについても検討してみよう。ランダムサーチはその名の通り値をランダムに選択する。`maxit`引数に試行回数を指定できる。

``` r
ctrl = makeTuneControlRandom(maxit = 200L)
```

チューニングの実行
------------------

パラメータの探索範囲と最適化アルゴリズムを決めたら、いよいよチューニングの実行の時だ。あとは、リサンプリング手法と評価尺度を設定する必要がある。

今回は3分割クロスバリデーションをパラメータ設定の評価に使用する。まずはリサンプリングdescriptionを生成する。

``` r
rdesc = makeResampleDesc("CV", iters = 3L)
```

では、今まで作成したものを組合せて、`tuneParams`関数によりパラメータチューニングを実行しよう。今回は`discrete_ps`に対してグリッドサーチを行う。

``` r
ctrl = makeTuneControlGrid()
res = tuneParams("classif.ksvm", task = iris.task, resampling = rdesc,
                 par.set = discrete_ps, control = ctrl)
```

    $> [Tune] Started tuning learner classif.ksvm for parameter set:

    $>           Type len Def      Constr Req Tunable Trafo
    $> C     discrete   -   - 0.5,1,1.5,2   -    TRUE     -
    $> sigma discrete   -   - 0.5,1,1.5,2   -    TRUE     -

    $> With control class: TuneControlGrid

    $> Imputation value: 1

    $> [Tune-x] 1: C=0.5; sigma=0.5

    $> [Tune-y] 1: mmce.test.mean=0.0533; time: 0.0 min

    $> [Tune-x] 2: C=1; sigma=0.5

    $> [Tune-y] 2: mmce.test.mean=0.06; time: 0.0 min

    $> [Tune-x] 3: C=1.5; sigma=0.5

    $> [Tune-y] 3: mmce.test.mean=0.0533; time: 0.0 min

    $> [Tune-x] 4: C=2; sigma=0.5

    $> [Tune-y] 4: mmce.test.mean=0.0533; time: 0.0 min

    $> [Tune-x] 5: C=0.5; sigma=1

    $> [Tune-y] 5: mmce.test.mean=0.0667; time: 0.0 min

    $> [Tune-x] 6: C=1; sigma=1

    $> [Tune-y] 6: mmce.test.mean=0.0667; time: 0.0 min

    $> [Tune-x] 7: C=1.5; sigma=1

    $> [Tune-y] 7: mmce.test.mean=0.0667; time: 0.0 min

    $> [Tune-x] 8: C=2; sigma=1

    $> [Tune-y] 8: mmce.test.mean=0.0733; time: 0.0 min

    $> [Tune-x] 9: C=0.5; sigma=1.5

    $> [Tune-y] 9: mmce.test.mean=0.0733; time: 0.0 min

    $> [Tune-x] 10: C=1; sigma=1.5

    $> [Tune-y] 10: mmce.test.mean=0.0733; time: 0.0 min

    $> [Tune-x] 11: C=1.5; sigma=1.5

    $> [Tune-y] 11: mmce.test.mean=0.0733; time: 0.0 min

    $> [Tune-x] 12: C=2; sigma=1.5

    $> [Tune-y] 12: mmce.test.mean=0.0733; time: 0.0 min

    $> [Tune-x] 13: C=0.5; sigma=2

    $> [Tune-y] 13: mmce.test.mean=0.0867; time: 0.0 min

    $> [Tune-x] 14: C=1; sigma=2

    $> [Tune-y] 14: mmce.test.mean=0.08; time: 0.0 min

    $> [Tune-x] 15: C=1.5; sigma=2

    $> [Tune-y] 15: mmce.test.mean=0.08; time: 0.0 min

    $> [Tune-x] 16: C=2; sigma=2

    $> [Tune-y] 16: mmce.test.mean=0.0733; time: 0.0 min

    $> [Tune] Result: C=2; sigma=0.5 : mmce.test.mean=0.0533

``` r
res
```

    $> Tune result:
    $> Op. pars: C=2; sigma=0.5
    $> mmce.test.mean=0.0533

`tuneParams`はパラメータの全ての組み合わせに対してクロスバリデーションによる性能評価を行い、最も良い値を出した組合せをパラメータとして採用する。性能指標を指定しなかった場合は誤分類率(mmce)が使用される。

それぞれの`measure`は、その値を最大化すべきか最小化すべきかを**知って**いる。

``` r
mmce$minimize
```

    $> [1] TRUE

``` r
acc$minimize
```

    $> [1] FALSE

もちろん、他の指標をリストとして同時に`tuneParams`に渡すこともできる。この場合、最初の指標が最適化に使われ、残りの指標は単に計算されるだけとなる。もし複数の指標を同時に最適化したいと考えるのであれば、[Advanced Tuning - mlr tutorial](https://mlr-org.github.io/mlr-tutorial/devel/html/advanced_tune/index.html)を参照してほしい。

誤分類率の代わりに精度(acc)を計算する例を示そう。同時に、他の性能指標として精度の標準偏差を求めるため、`setAggregation`関数を使用している。また、今回は探索範囲`num_set`に対して100回のランダムサーチを行う。100回分の出力は長くなるので、`show.info = FALSE`を指定している。

``` r
ctrl = makeTuneControlRandom(maxit = 100L)
res = tuneParams("classif.ksvm", task = iris.task, resampling = rdesc, par.set = num_ps,
                 control = ctrl, measures = list(acc, setAggregation(acc, test.sd)), show.info = FALSE)
res
```

    $> Tune result:
    $> Op. pars: C=1.79e+04; sigma=5.07e-06
    $> acc.test.mean=0.973,acc.test.sd=0.0115

チューニング結果へのアクセス
----------------------------

チューニングの結果は`TuneResult`クラスのオブジェクトである。見つかった最適値は`$x`スロット、性能指標については`$y`スロットを通じてアクセスできる。

``` r
res$x
```

    $> $C
    $> [1] 17936.23
    $> 
    $> $sigma
    $> [1] 5.074312e-06

``` r
res$y
```

    $> acc.test.mean   acc.test.sd 
    $>    0.97333333    0.01154701

最適化されたパラメータをセットした学習器は次のように作成できる。

``` r
lrn = setHyperPars(makeLearner("classif.ksvm"), par.vals = res$x)
lrn
```

    $> Learner classif.ksvm from package kernlab
    $> Type: classif
    $> Name: Support Vector Machines; Short name: ksvm
    $> Class: classif.ksvm
    $> Properties: twoclass,multiclass,numerics,factors,prob,class.weights
    $> Predict-Type: response
    $> Hyperparameters: fit=FALSE,C=1.79e+04,sigma=5.07e-06

あとはこれまでと同じだ。irisデータセットに対して再度学習と予測を行ってみよう。

``` r
m = train(lrn, iris.task)
predict(m, task = iris.task)
```

    $> Prediction: 150 observations
    $> predict.type: response
    $> threshold: 
    $> time: 0.00
    $>   id  truth response
    $> 1  1 setosa   setosa
    $> 2  2 setosa   setosa
    $> 3  3 setosa   setosa
    $> 4  4 setosa   setosa
    $> 5  5 setosa   setosa
    $> 6  6 setosa   setosa
    $> ... (150 rows, 3 cols)

しかし、この方法だと最適化された状態のハイパーパラメータの影響しか見ることができない。検索時に生成された他の値を使った場合の影響はどのように確認すれば良いだろうか？

ハイパーパラメータチューニングの影響を調査する
----------------------------------------------

`generateHyperParsEffectData`を使うと、サーチ中に生成された全ての値について調査を行うことができる。

``` r
generateHyperParsEffectData(res)
```

    $> HyperParsEffectData:
    $> Hyperparameters: C,sigma
    $> Measures: acc.test.mean,acc.test.sd
    $> Optimizer: TuneControlRandom
    $> Nested CV Used: FALSE
    $> Snapshot of data:
    $>           C     sigma acc.test.mean acc.test.sd iteration exec.time
    $> 1  3.729807 -5.483632     0.9200000  0.05291503         1     0.050
    $> 2 -3.630108 -6.520324     0.2933333  0.02309401         2     0.050
    $> 3  6.028592 -7.074359     0.9600000  0.02000000         3     0.046
    $> 4  6.454348 -3.380043     0.9400000  0.03464102         4     0.046
    $> 5 -2.516612  2.594908     0.2933333  0.02309401         5     0.051
    $> 6 -8.067325 -9.560126     0.2933333  0.02309401         6     0.046

この中に含まれているパラメータの値は**オリジナルのスケール**であることに注意しよう。`trafo`に指定した関数で変換後の値が欲しければ、`trafo`引数に`TRUE`を指定する必要がある。

``` r
generateHyperParsEffectData(res, trafo = TRUE)
```

    $> HyperParsEffectData:
    $> Hyperparameters: C,sigma
    $> Measures: acc.test.mean,acc.test.sd
    $> Optimizer: TuneControlRandom
    $> Nested CV Used: FALSE
    $> Snapshot of data:
    $>              C        sigma acc.test.mean acc.test.sd iteration exec.time
    $> 1 5.367932e+03 3.283738e-06     0.9200000  0.05291503         1     0.050
    $> 2 2.343645e-04 3.017702e-07     0.2933333  0.02309401         2     0.050
    $> 3 1.068050e+06 8.426382e-08     0.9600000  0.02000000         3     0.046
    $> 4 2.846740e+06 4.168277e-04     0.9400000  0.03464102         4     0.046
    $> 5 3.043601e-03 3.934671e+02     0.2933333  0.02309401         5     0.051
    $> 6 8.563965e-09 2.753432e-10     0.2933333  0.02309401         6     0.046

また、リサンプリングの部分で説明したように、テストデータに加えて訓練データに対しても性能指標を求められることに注意してもらいたい。

``` r
rdesc2 = makeResampleDesc("Holdout", predict = "both")
res2 = tuneParams("classif.ksvm", task = iris.task, resampling = rdesc2, par.set = num_ps,
                  control = ctrl, measures = list(acc, setAggregation(acc, train.mean)), show.info = FALSE)
generateHyperParsEffectData(res2)
```

    $> HyperParsEffectData:
    $> Hyperparameters: C,sigma
    $> Measures: acc.test.mean,acc.train.mean
    $> Optimizer: TuneControlRandom
    $> Nested CV Used: FALSE
    $> Snapshot of data:
    $>            C     sigma acc.test.mean acc.train.mean iteration exec.time
    $> 1  6.5343098  7.837925          0.30           1.00         1     0.026
    $> 2  0.7882147  1.106272          0.80           1.00         2     0.023
    $> 3 -4.4712023  3.183959          0.28           0.36         3     0.027
    $> 4 -5.1312952 -4.234676          0.28           0.36         4     0.025
    $> 5  4.3156192  4.550673          0.30           1.00         5     0.026
    $> 6 -6.5391392 -3.873775          0.28           0.36         6     0.023

パラメータ値の評価結果は`plotHyperParsEffect`関数を使うと簡単に可視化できる。例を示そう。以下では、繰り返し毎の性能指標の変化をプロットしている。ここで`res`は先に示したものとほぼ同じだが、2つの性能指標を使用している。

``` r
res = tuneParams("classif.ksvm", task = iris.task, resampling = rdesc, par.set = num_ps,
                 control = ctrl, measures = list(acc, mmce), show.info = FALSE)
data = generateHyperParsEffectData(res)
plotHyperParsEffect(data, x = "iteration", y = "acc.test.mean", plot.type = "line")
```

![unnamed-chunk-22-1.png](https://qiita-image-store.s3.amazonaws.com/0/109348/b51227bb-dda8-d0d6-904a-8c42015d5426.png)