Rのplmパッケージのpgmmクラスで算出したモデルの情報基準(AIC, BIC)を求める関数の作成

とある案件で、R言語を用いて動学的パネル分析を行おうと思ったのですが、モデル選択の際、どうやらplmパッケージのpgmmクラスに情報基準(AIC,BIC)が存在しないようなので(R)自作してみたというお話。もちろん、chat-gpt4のサポートありです。

・まずは、適当にモデルの定義(pdataはpdata.frameクラスに変換済)。

#モデルの定義
z <- pgmm(y ~  lag(y,1) + lag(x,0:2) |
                   unit + time,
           data = pdata,
           effect = "twoways",
           model = "twosteps",
           transformation = 'ld')

・二つの情報基準を求める関数を作成(zはpgmmクラスを引き数とする)

#aicを求める関数
calculate_aic <- function(z) {
        residuals_z <- as.data.frame(residuals(z))　#モデルの残差の抽出
        log_lik_z <- -0.5 * sum(log(2 * pi * (residuals_z^2)) + 1)　#対数尤度（近似値）の計算
        k <- length(z$coefficients[[1]])　#次数の算出
        AIC_z <- -2 * log_lik_z + 2 * k #AICの算出
        return(AIC_z)　#リターンとしてAICを設定
}

#bicを求める関数
calculate_bic <- function(z) {
        residuals_z <- as.data.frame(residuals(z))　#モデルの残差の抽出
        log_lik_z <- -0.5 * sum(log(2 * pi * (residuals_z^2)) + 1)　#対数尤度（近似値）の計算
        k <- length(z$coefficients[[1]]) #次数の算出
        n <- nrow(residuals_z) * ncol(residuals_z) #観測数の算出
        BIC_z <- -2 * log_lik_z + k * log(n) #BICの算出
        return(BIC_z)　#リターンとしてAICを設定
}

・最後に関数への代入と、書き出し代入するのは、先ほど作成したモデルz

#モデルzのAICの代入と書き出し
AIC_z <- calculate_aic(z)
print(AIC_z)

#モデルzのBICの算出と書き出し
BIC_z <- calculate_bic(z)
print(BIC_z)

もう少し、オブジェクト指向っぽい書き方をすればよかったのかもしれないです。あと恥ずかしながらこれまで、動学パネルやGMMによる解析を扱ってこなかったので、適切な評価指標については、要学習といった感じです。