(cache)モデルで扱うデータの前処理をrecipesで行う

ドーモ。ホクソエムの @u_ribo です。本業ではモデリングとは離れたギョームをしています。寂しくなったので、Rのrecipesパッケージについて紹介します。

tidymodels.github.io

モデルに適用するデータの前処理

Rでのモデル式 (model formula) の記述って、利用時に不便を感じることや覚えるのが難しい面が時々ありませんか？例えば、y ~ .」は右辺のドットが、目的変数以外の全ての変数を説明変数として扱うことを示しますが、説明変数に対数変換などの変数変換を行うにはy ~ log(.)という記述はできず、結局、説明変数を「+」でつなげていくことになります。また、交互作用項の指定には「x1 * x2」や「(x1 + x2)^2」、「:」を使う表記が可能ですが、この表記には最初は混乱しませんか？（単に私が不勉強なだけということもあります）

加えて、多くのモデルでは欠損への処理が必要となったり、適用するモデルによって（例えば複数の説明変数を扱うK-NNやSVM）は、変数間の標準化が必要です。

ここで紹介するrecipesパッケージを使うことで、こうしたモデル構築に伴うデータの前処理を楽にすることが期待されます。

モデルに適用するデータの前処理
対数変換を行うモデルを例に
- recipesパッケージによるデータの前処理
交互作用を扱うモデルを例に
stepに指定可能な処理
まとめ

対数変換を行うモデルを例に

recipesの例を示すのに、まずは標準のモデル構築の手続きを見てみます。次のモデルは、「アヒル本」こと「StanとRでベイズ統計モデリング (松浦(2016.) 共立出版)」の第7章で取り上げられている関東近郊の架空物件についての費用(Y)を説明する変数として物件の広さ(Area)を適用したもので、変数変換の例として応答変数、説明変数ともに対数変換を施しています。

library(magrittr)
library(broom)

# サポートサイトから架空物件データを読み込む
df_rental <-
  readr::read_csv("https://raw.githubusercontent.com/MatsuuraKentaro/RStanBook/master/chap07/input/data-rental.txt")

lm(formula = log(Y) ~ log(Area), data = df_rental) %>% 
  # 推定された効果を出力
  tidy()

 ## # A tibble: 2 x 5
 ##   term        estimate std.error statistic  p.value
 ##   <chr>          <dbl>     <dbl>     <dbl>    <dbl>
 ## 1 (Intercept)     1.84    0.196       9.38 2.73e-15
 ## 2 log(Area)       1.11    0.0557     19.9  3.72e-36

線形回帰モデルを実行するlm()のformulaにモデル式を定義します。このモデル式を読むと、

Yを目的変数とする
Areaを説明変数として扱う
目的変数および説明変数に与えた2つの変数を対数変換する

という手順になっているのがわかります。

recipesでは、こうしたモデルの記述とデータの用意を、tidyverseでおなじみのパイプ演算子を使って段階に分けて記述していきます。それによりシームレスなモデル適用の手続きを可能にします。百聞は一見に如かずということで、recipesの働きを見ていきましょう。

ここでCMです

まだまだ絶賛発売中！皆様是非ッ！！

StanとRでベイズ統計モデリング (Wonderful R)

作者: 松浦健太郎,石田基広
出版社/メーカー: 共立出版
発売日: 2016/10/25
メディア: 単行本
この商品を含むブログ (10件) を見る

前処理の大事さやデータに応じた処理の適用方法についての詳細は

前処理大全[データ分析のためのSQL/R/Python実践テクニック]

作者: 本橋智光
出版社/メーカー: 技術評論社
発売日: 2018/04/13
メディア: 大型本
この商品を含むブログ (1件) を見る

を是非ッ！！！

Osaka.Stan #4 Chapter 7 回帰分析の悩みどころ (7.1–7.5) from Hiroyuki Muto

www.slideshare.net

recipesパッケージによるデータの前処理

recipesではモデルの記述と変数の処理を調理作業に例えて表します。まず食材を用意し、何を作るかを考え（recipes()でモデルを定義）、準備して(prep())、料理を作っていきます(juiceやbake())。

recipe() --> prep() --> bake() and juice()

ではrecipesを使った例を示します。

library(recipes)

mod_rec <- 
  df_rental %>% 
  recipe(formula = Y ~ Area)

まず、recipe()でモデル式の定義をします。ここでは変数変換するという記述はしません。データと、データに含まれるどの変数をモデルに用いるかだけを指定します。

summary(mod_rec)

## # A tibble: 2 x 4
##   variable type    role      source  
##   <chr>    <chr>   <chr>     <chr>   
## 1 Area     numeric predictor original
## 2 Y        numeric outcome   original

summary()でrecipe()の適用結果を見ると、各変数がモデル式の中でどのような関係を持っているか、数値なのかカテゴリなのかという情報がわかります。ここでは2変数が定義されていて、数値からなる説明変数が1つあることを示しています。

モデルの変数を対数にとる処理はstep_log()を使って指定します。レシピには、食材をどのように処理して料理を作るのかの記述がないといけませんよね。次のコードで、レシピに対数変換を行う手順 (step) を追加するという意味になります。

mod_rec <- 
  mod_rec %>% 
  # step_log()の引数に、対数変換する変数を指定する
  step_log(Y, Area)

step_log()の引数中で対象となる変数を指定しています。今回のモデルでは説明変数は1つだけですが、複数の変数を変換する必要がある場合にはall_predictors()という変数を使うと便利です。recipesにはこうした、モデル中で変数の役割や数値かカテゴリかというデータの種類に応じた変数選択のための関数が用意されていて、下記の処理は上記の処理と同じものです。

mod_rec %>% 
  step_log(all_predictors(), all_outcomes())

この他の変数選択の補助関数については?has_role()で確認できます。またstep_*()ではtidyselectパッケージを利用した効率的な列選択も可能です。すなわちdplyrパッケージ等でも使えるstarts_with()やcontains()です。

話題を元に戻します。ここまでで「目的変数と説明変数を対数に変換する」という単純なレシピができました。次はこのレシピを架空物件データに適用します。まずはprep()でデータを当てはめます。関数を実行すると、step_log()の対象となるArea, Yが対数変換される処理が確認できます。

rec_trained <- 
  prep(mod_rec, retain = TRUE, verbose = TRUE)

## oper 1 step log [training]

# df_rentalのY, Area列が対数をとった値になっている
rec_trained

## Data Recipe
## 
## Inputs:
## 
##       role #variables
##    outcome          1
##  predictor          1
## 
## Training data contained 100 data points and no missing data.
## 
## Operations:
## 
## Log transformation on Y, Area [trained]

一方でまだ料理は完成ではありません。変数変換した値がどうなっているかを見ていませんよね。レシピの最初に与えた架空物件データ df_rental をレシピのどうりに調理するにはjuice()を使います。その後、目的のlm()を実行します。途中、print()を使ってjuice()の結果を見ながら処理の流れを確認します。

rec_trained %>% 
  juice() %>%
  # Y, Areaがそれぞれ対数変換されているのを確認
  print() %>% 
  lm(Y ~ Area, data = .) %>% 
  tidy()

## # A tibble: 100 x 2
##    Area     Y
##   <dbl> <dbl>
## 1  3.71  5.57
## 2  4.13  6.25
## 3  4.01  6.51
## 4  4.04  6.16
## 5  4.04  6.12
## 6  2.74  5.03
## # ... with 94 more rows

## # A tibble: 2 x 5
##   term        estimate std.error statistic  p.value
##   <chr>          <dbl>     <dbl>     <dbl>    <dbl>
## 1 (Intercept)     1.84    0.196       9.38 2.73e-15
## 2 Area            1.11    0.0557     19.9  3.72e-36

注目して欲しいのは、lm()内のモデル式にはlog()が使われていない点です。これはrecipesによりデータの変数変換の処理がすでに行われているためです。実行結果は当然最初のlog(Y) ~ log(Area)と同じものです。

また今回は元のデータにレシピを与えましたが、同じモデルを異なるデータに適用するにはbake(newdata = )に対象となるデータを与えます。これは教師データでモデルを用意し、テストデータを適用する時などに役立ちます。

rec_test <- 
  rec_trained %>% 
  bake(newdata = df_rental)

rec_test %>% 
  lm(Y ~ Area, data = .) %>% 
  tidy()

## # A tibble: 2 x 5
##   term        estimate std.error statistic  p.value
##   <chr>          <dbl>     <dbl>     <dbl>    <dbl>
## 1 (Intercept)     1.84    0.196       9.38 2.73e-15
## 2 Area            1.11    0.0557     19.9  3.72e-36

交互作用を扱うモデルを例に

今度は交互作用を扱うモデルのデータをreciesで用意する例を示します。 Rに用意されている、1973年のニューヨーク州の大気データセット (airquality) から、Ozoneを説明するモデルを考えます。ここでは交互作用を見るために日射量Solar.Rと風力Windを指定します。モデル式は次のように表現できます。

Ozone ~ (Solar.R + Wind)^2

このデータをrecipesで作ります。対数変換の例と同じく、データセットとモデルに使う変数の定義をrecipe()で行います。

mod_rec_airq <- 
  airquality %>% 
  recipe(Ozone ~ Solar.R + Wind)

交互作用の処理はstep_interact()で指定します。その後の処理は先のモデルと同じです。データへのレシピの適用まで一気にやってしまいましょう。

rec_test_airq <- 
  mod_rec_airq %>% 
  # (Solar.R + Wind)^2 としても良い
  step_interact(terms = ~ Solar.R:Wind) %>% 
  prep(retain = TRUE) %>% 
  juice()

交互作用を含めた処理を適用したデータを見ると、Solar.R_x_Windなる列が追加されています。この列はstep_interact()に指定した記述により作成された値からなります。

rec_test_airq

## # A tibble: 153 x 4
##   Solar.R  Wind Ozone Solar.R_x_Wind
##     <int> <dbl> <int>          <dbl>
## 1     190   7.4    41          1406 
## 2     118   8      36           944 
## 3     149  12.6    12          1877.
## 4     313  11.5    18          3600.
## 5      NA  14.3    NA            NA 
## 6      NA  14.9    28            NA 
## # ... with 147 more rows

それでは回帰モデルを実行します。交互作用項はモデルで記述するのではなく説明変数として与えます。すなわち次のようになります。

lm(Ozone ~ Solar.R + Wind + Solar.R_x_Wind, data = rec_test_airq) %>% 
  tidy()

## # A tibble: 4 x 5
##   term           estimate std.error statistic  p.value
##   <chr>             <dbl>     <dbl>     <dbl>    <dbl>
## 1 (Intercept)     34.5     17.6          1.95 0.0532  
## 2 Solar.R          0.324    0.0839       3.86 0.000193
## 3 Wind            -1.59     1.51        -1.06 0.293   
## 4 Solar.R_x_Wind  -0.0203   0.00725     -2.80 0.00609

lm(Ozone ~ (Solar.R + Wind)^2, data = airquality) %>% 
  tidy()

 ## # A tibble: 4 x 5
 ##   term         estimate std.error statistic  p.value
 ##   <chr>           <dbl>     <dbl>     <dbl>    <dbl>
 ## 1 (Intercept)   34.5     17.6          1.95 0.0532  
 ## 2 Solar.R        0.324    0.0839       3.86 0.000193
 ## 3 Wind          -1.59     1.51        -1.06 0.293   
 ## 4 Solar.R:Wind  -0.0203   0.00725     -2.80 0.00609

モデル式の構文、(Solar.R + Wind)^2を使った場合と同じ結果が得られました。

このモデルから交互作用項を取り除くには、Solar.R_x_Windをモデル式から除くだけなので応用も簡単です。

lm(Ozone ~ Solar.R + Wind, data = rec_test_airq) %>% 
  tidy()

## # A tibble: 3 x 5
##   term        estimate std.error statistic  p.value
##   <chr>          <dbl>     <dbl>     <dbl>    <dbl>
## 1 (Intercept)   77.2      9.07        8.52 1.05e-13
## 2 Solar.R        0.100    0.0263      3.82 2.24e- 4
## 3 Wind          -5.40     0.673      -8.02 1.34e-12

stepに指定可能な処理

step_log()以外にも、モデル構築時の前処理に便利なstep_*()というプレフィックスの関数が豊富に用意されています。現在では以下の場面で使える関数があります。

Basic: logs, roots, polynomials, logits, hyperbolics
Encodings: dummy variables, “other” factor level collapsing, discretization
Date Features: Encodings for day/doy/month etc, holiday indicators
Filters: correlation, near-zero variables, linear dependencies
Imputation: K-nearest neighbors, bagged trees, mean/mode imputation,
Normalization/Transformations: center, scale, range, Box-Cox, Yeo-Johnson
Dimension Reduction: PCA, kernel PCA, ICA, Isomap, data depth features, class distances
Others: spline basis functions, interactions, spatial sign

recipesでは、モデルの実行に欠かせないデータの前処理をこうしたステップにわけ、パイプ演算子と組みわせることで複数の処理を適用することを可能にしています。

処理を適用する変数の選択には all_predictors()やall_numeric()といった変数に応じた選択に加えて、tidyselectの関数が利用でき、柔軟な処理ができるようになっています。

まとめ

recipesパッケージでは統計解析や機械学習の各種モデリングを適用する際に必要となってくるデータの前処理を行うのに役立つ
recipesによるデータ変換はrecipe(), prep(), juice() and bake() という手順をとる
step_*()による豊富な前処理関数が用意されている
変数選択には、全ての説明変数、全ての数値変数を指定したりtidyselectによる関数の適用が可能

興味を持たれた方は、パッケージのvignettesやuseR!での発表やRStudioのwebinar動画をご覧ください。

Enjoy!

株式会社ホクソエムのブログ