重回帰分析における AIC, BIC など

以下のデータ例をモデル式 y ~ x1 + x2 で重回帰分析する。

using DataFrames

df = DataFrame(
    x1 = [46, 49, 55, 55, 52, 51, 53, 50, 50, 57, 51, 37, 45, 45, 54, 47, 54, 55, 40, 54],
    x2 = [49, 36, 50, 60, 54, 63, 44, 53, 53, 53, 48, 41, 44, 45, 53, 46, 59, 60, 47, 42],
    y  = [49, 40, 53, 61, 59, 50, 50, 50, 51, 52, 47, 41, 42, 44, 51, 48, 56, 61, 45, 47]);

using GLM, StatsBase

results = lm(@formula(y ~ x1 + x2), df)

StatsModels.TableRegressionModel{LinearModel{GLM.LmResp{Vector{Float64}}, GLM.DensePredChol{Float64, LinearAlgebra.CholeskyPivoted{Float64, Matrix{Float64}}}}, Matrix{Float64}}

y ~ 1 + x1 + x2

Coefficients:
────────────────────────────────────────────────────────────────────────
                Coef.  Std. Error     t  Pr(>|t|)   Lower 95%  Upper 95%
────────────────────────────────────────────────────────────────────────
(Intercept)  1.64221     6.8959    0.24    0.8146  -12.9069    16.1913
x1           0.433718    0.154839  2.80    0.0123    0.107038   0.760399
x2           0.530437    0.11467   4.63    0.0002    0.288505   0.772369
────────────────────────────────────────────────────────────────────────

残差分散

サンプルサイズを $n$，独立変数の個数を $p$ とする。
残差平方和を $S_e$ とすれば，残差分散 $d$ は（1）式で定義される。

$\hspace{2cm}\displaystyle d = \sum_{i=1}^n (y_i-\hat{y_i})^2/n = S_e/n　　　　(1)$

n = nobs(results) # size(df, 1)
p = length(coefnames(results)) - 1
Se = sum((response(results) .- predict(results)) .^ 2)

159.30705544126175

d = Se / n

7.9653527720630874

対数尤度

対数尤度 $loglik$ は（2）式で定義される。

$\hspace{2cm}loglik = -0.5(n\log 2\pi + n\log S_e -n\log n + n) = -0.5(n\log 2\pi + n\log d + n)　　　　(2)$

loglik = -0.5*(n*log(2*pi)+n*log(Se)-n*log(n)+n) # loglik: 定義式により計算する

-49.129782990684234

-0.5*(n*log(2*pi)+n*log(d)+n) # loglik: 定義式により計算する

-49.129782990684234

loglikelihood(results) # loglikelihood メソッドで計算できる

-49.12978299068424

AIC

AIC の定義は（3）式で定義される。

$\hspace{2cm}AIC = -2loglik+2(p+2) = n\log 2\pi + n\log d + n + 2(p+2)　　　　(3)$

-2*loglik+2*(p+2) # AIC: 定義式により計算する

106.25956598136847

n*log(2*pi)+n*log(d)+n+2*(p+2) # AIC: 定義式により計算する

106.25956598136847

aic(results) # aic メソッドで計算できる

106.25956598136848

なお，（3）式から定数部分を除いた（4）式による値を AIC とする場合もあるので注意が必要である。

$\hspace{2cm}AIC = n\log d + 2(p+1)　　　　(4)$

n*log(d)+2*(p+1) # AIC: 定義式により計算する

47.50202465318156

BIC

一般的に，「$-2\ \times$ 対数尤度$+ k\ \times$ パラメータ数」とすれば，AIC では $k=2$ である。

$k=\log($サンプルサイズ$)$ としたものは BIC または SBC(Schwarz's Bayesian criterion) と呼ばれる。

n*log(2*pi)+n*log(d)+n+2*(p+2) # AIC: 定義式により計算する

106.25956598136847

n*log(2*pi)+n*log(d)+n+log(n)*(p+2) # BIC: 定義式により計算する

110.24249507558443

bic(results) # bic メソッドで計算できる

110.24249507558444