井出草平の研究ノート

memiscパッケージで分析結果を出力

memiscパッケージのmtableコマンドを使うと分析結果などをきれいに出力することができるようだ。

Package ‘memisc’
https://cran.r-project.org/web/packages/memisc/memisc.pdf

事前準備

回帰分析の結果を出すため、サンプルデータから分析を行っておく。以前のエントリーと同様にAERパッケージに入っているCPS1985というデータを使う。

library(AER) #AERパッケージの読み込み
data(CPS1985) #CPS1985の読み込み
lm0 <- lm(wage ~ age, data=CPS1985) #年齢が賃金を決めるモデル
lm1 <- update(lm0, ~. + gender) # 性別をモデルに加える
lm2 <- update(lm1, ~. + education) # 学歴をモデルに加える
lm3 <- update(lm2, ~. + experience) # 仕事の経験年数をモデルに加える

これで4つの回帰分析が実行された。この4つの結果をmtableコマンドでまとめる。

library(memisc) #memiscパッケージの読み込み
mt1234 <- mtable("Model1"=lm0, "Model2"=lm1, "Model3"=lm2,
                 "Model4"=lm3, #モデル名を指定
          summary.stats=c("sigma","R-squared","F","p","BIC","N")
          #出力する統計量を指定
)
print(mt1234) #テキストで出力する。

出力結果は以下のようになる。

==============================================================
                   Model1     Model2     Model3     Model4    
--------------------------------------------------------------
  (Intercept)      6.167***   6.929***  -4.843***  -1.957     
                  (0.723)    (0.720)    (1.244)    (6.835)    
  age              0.078***   0.085***   0.113***  -0.367     
                  (0.019)    (0.018)    (0.017)    (1.120)    
  gender: female             -2.275***  -2.335***  -2.344***  
                             (0.430)    (0.388)    (0.389)    
  education                              0.827***   1.307     
                                        (0.075)    (1.120)    
  experience                                        0.481     
                                                   (1.121)    
--------------------------------------------------------------
  sigma               5.1        4.9        4.5        4.5    
  R-squared           0.0        0.1        0.3        0.3    
  F                  17.2       23.0       59.9       44.9    
  p                   0.0        0.0        0.0        0.0    
  BIC              3264.5     3243.4     3138.2     3144.3    
  N                 534        534        534        534      
==============================================================

オプション

さまざまなオプションが用意されている。

HTMLで出力する

show_html(mt1234)

HTML形式のファイルを書き出すのは下記のコマンド。

write_html(mt1234,file="mt1234.html")

Rのワーキング・ディレクトリに出力される。

TeX形式で出力する。

toLatex(mt1234)

TeXファイルとして書き出すのは以下のコマンド。

write.mtable(mt1234,format="LaTeX", file="mt123.tex")

テキスト形式で出力

テキスト形式でも出力できる。下記のコマンドはタブ区切り。

write.mtable(mt1234,file="mt1234.txt")

CSV形式で出力

CSV形式だとExcelなどの加工が楽になる。colsep=","を入れてコンマ区切りを指定する。

write.mtable(mt1234,file="mt1234.csv",colsep=",")

独立変数をリラベル

独立変数の変数名が英語であったり略語であったりする場合リラベルできる。

mt1234 <- relabel(mt1234,
                 "(Intercept)" = "(切片)", #独立変数のリラベル
                age = "年齢",
               "gender: female" = "性別",
                education = "教育年数",
               experience = "就労経験"
)

独立変数が変更されていることが確認できる。

=========================================================
              Model1     Model2     Model3     Model4    
---------------------------------------------------------
  (切片)      6.167***   6.929***  -4.843***  -1.957     
             (0.723)    (0.720)    (1.244)    (6.835)    
  年齢        0.078***   0.085***   0.113***  -0.367     
             (0.019)    (0.018)    (0.017)    (1.120)    
  性別                  -2.275***  -2.335***  -2.344***  
                        (0.430)    (0.388)    (0.389)    
  教育年数                          0.827***   1.307     
                                   (0.075)    (1.120)    
  就労経験                                     0.481     
                                              (1.121)    
---------------------------------------------------------

統計量の出力を選ぶ

mt1234 <- mtable("Model1"=lm0, "Model2"=lm1, "Model3"=lm2,
                 "Model4"=lm3, #モデル名を指定
          summary.stats=c("sigma","R-squared","F","p","BIC","N")
          #出力する統計量を指定
)

summary.statsの所の指定は適宜変える。指定できるのは下記の統計量。

R-squared
adj. R-squared sigma
F
p
Log-likelihood
Deviance AIC BIC N

有意のアスタリスクを他の記号に変更する。

mt1234 <- mtable("Model1"=lm0, "Model2"=lm1, "Model3"=lm2, "Model4"=lm3,
           signif.symbols=c("a"=.05,
                            "b"=.01,
                            "c"=.001)
)
======================================================
                  Model1   Model2   Model3   Model4   
------------------------------------------------------
  (Intercept)      6.167c   6.929c  -4.843c  -1.957   
                  (0.723)  (0.720)  (1.244)  (6.835)  
  age              0.078c   0.085c   0.113c  -0.367   
                  (0.019)  (0.018)  (0.017)  (1.120)  
  gender: female           -2.275c  -2.335c  -2.344c  
                           (0.430)  (0.388)  (0.389)  
  education                          0.827c   1.307   
                                    (0.075)  (1.120)  
  experience                                  0.481   
                                             (1.121)  
======================================================