broomを使って回帰分析の結果を出力する

いわずもがなであるが、『効果検証入門』は大変勉強になる。因果推論は大学・大学院のときに学んでいるが、その復習に使える。また、Rを使ってどのように因果推論を行うかについても有益な情報を与える。

作者:安井翔太
技術評論社

その中で、回帰分析の結果の出力方法について私が知らないやり方があったので、ここで書いてみる。それがbroomというパッケージを使用した方法である。

通常Rで線形回帰分析を行うときはlm関数を用いる。例のごとく、irisを用いて適当な線形回帰モデルを推定し、その結果をsummaryで出力する。

reg <- 
  lm(Sepal.Length ~ Sepal.Width + Petal.Length + Petal.Width, data = iris)
summary(reg)

Call:
lm(formula = Sepal.Length ~ Sepal.Width + Petal.Length + Petal.Width, 
    data = iris)

Residuals:
     Min       1Q   Median       3Q      Max 
-0.82816 -0.21989  0.01875  0.19709  0.84570 

Coefficients:
             Estimate Std. Error t value Pr(>|t|)    
(Intercept)   1.85600    0.25078   7.401 9.85e-12 ***
Sepal.Width   0.65084    0.06665   9.765  < 2e-16 ***
Petal.Length  0.70913    0.05672  12.502  < 2e-16 ***
Petal.Width  -0.55648    0.12755  -4.363 2.41e-05 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 0.3145 on 146 degrees of freedom
Multiple R-squared:  0.8586,    Adjusted R-squared:  0.8557 
F-statistic: 295.5 on 3 and 146 DF,  p-value: < 2.2e-16

この方法では、下に余計な情報が出てくるなど見づらい、kableExtraなどで出力できないなど不便なことが多い。

そこでbroom::tidyの出番である。推定された係数の情報のみを抽出し、tibbleとして出力してくれる。

reg <- 
  lm(Sepal.Length ~ Sepal.Width + Petal.Length + Petal.Width, data = iris) %>%
  broom::tidy()

# A tibble: 4 x 5
  term         estimate std.error statistic  p.value
  <chr>           <dbl>     <dbl>     <dbl>    <dbl>
1 (Intercept)     1.86     0.251       7.40 9.85e-12
2 Sepal.Width     0.651    0.0666      9.77 1.20e-17
3 Petal.Length    0.709    0.0567     12.5  7.66e-25
4 Petal.Width    -0.556    0.128      -4.36 2.41e- 5

大概の場合、関心のあるのは係数の推定値のみであるから、これでほとんど事足りる。加えてtibbleとして出力してくれるので、例えばkableExtra::kable(reg)として出力することができる。

論文への出力のためのtableはstargazerを用いることが多いが、書き方によってはこのbroomを用いることで置き換えることが可能かもしれない。