［R］貨幣賃金率変化率の回帰モデルの推定（「線形回帰分析」（朝倉書店）pp.119-120）

表3.7の値をCSV形式で入力したtable3_7.csvを、カレントディレクトリに置いておく。

> dtf <- read.csv("table3_7.csv", header = TRUE)
> xx1i <- 1 / as.double(dtf$RU)
> xx2i <- as.double(dtf$CPIDOT)
> yyi <- as.double(dtf$WDOT)
> n <- length(xx1i)
> k <- 3
> degf <- n - k
> mxy <- matrix(yyi, ncol = 1)
> mxxx <- matrix(c(rep(1.0, n), xx1i, xx2i), ncol = 3)
> mxxxtxxi <- solve(t(mxxx) %*% mxxx)
> mxbh <- mxxxtxxi %*% t(mxxx) %*% mxy
> mxyh <- mxxx %*% mxbh
> mxe <- mxy - mxyh
> sse2 <- t(mxe) %*% mxe
> s2 <- as.vector(sse2 / degf)
> s <- sqrt(s2)
> mxvv <- s2 * mxxxtxxi
> mxs <- matrix(sqrt(diag(mxvv)), ncol = 1)
> tj <- as.vector(mxbh / mxs)
> pj <- 2 * pt(abs(tj), degf, lower.tail = FALSE)
> bh <- as.vector(mxbh)
> yyhi <- as.vector(mxyh)
> ei <- yyi - yyhi
> yym <- mean(yyi)
> rr2 <- sum((yyhi - yym) ^ 2) / sum((yyi - yym) ^ 2)
> rrb2 <- 1 - (n - 1) / (n - k) * (1 - rr2)
> print(bh)  # 推定した回帰モデルのパラメーター
[1] -4.6884718 16.2320984  0.8915398
> print(tj)  # t値
[1] -8.938476 12.501712 12.623158
> print(pj)  # p値
[1] 6.096213e-12 5.248937e-17 3.628214e-17
> print(rr2)  # 決定係数
[1] 0.9509842
> print(rrb2)  # 自由度修正済み決定係数
[1] 0.9490236
> print(s)  # 誤差項の分散の不偏推定量の標準誤差
[1] 1.580779