以下はWindows10、R version 4.2.1で実行しています。
データを読み込むにあたって、作業場所の指定とそこにデータが置いてあることが前提になります。 R本体を操作している場合は、「ファイル」→「ディレクトリの変更」でデータの置いてある場所を指定するのが簡単です。 RStudioを操作している場合は、例えばデスクトップの「R」という名前のフォルダにデータがあるとすると
setwd("C:/Users/ユーザー名/Desktop/R")を最初に実行するのがよいでしょう。ユーザー名は各自異なるので注意です。
使用するデータを読み込み、表示します。なお、R version 4.1以前で読み込む際には
data1 <- read.csv("ファイル名.csv", fileEncoding = "UTF-8-BOM")のように、エンコードのオプションを指定する必要があります。version 4.2以降では必要ありませんので、以下ではオプションを指定していません。
data1 <- read.csv("ch10-12ex.csv")使用するパッケージを読み込みますが、パッケージはインストールしておく必要があります。 以下のようにコマンドでインストールするか、RやRStudioからクリックでインストールすることもできます。
install.packages(“stargazer”, dependencies = TRUE)パッケージの機能が使えるように読み込みます。
library(dplyr)##
## 次のパッケージを付け加えます: 'dplyr'## 以下のオブジェクトは 'package:stats' からマスクされています:
##
## filter, lag## 以下のオブジェクトは 'package:base' からマスクされています:
##
## intersect, setdiff, setequal, unionlibrary(psych)
library(car)## 要求されたパッケージ carData をロード中です##
## 次のパッケージを付け加えます: 'car'## 以下のオブジェクトは 'package:psych' からマスクされています:
##
## logit## 以下のオブジェクトは 'package:dplyr' からマスクされています:
##
## recodelibrary(stargazer)##
## Please cite as:## Hlavac, Marek (2022). stargazer: Well-Formatted Regression and Summary Statistics Tables.## R package version 5.2.3. https://CRAN.R-project.org/package=stargazer従属変数の準備をします。
data1 <- data1 %>% mutate(
q651=recode(q06_5_1,'99=NA'),
q652=recode(q06_5_2,'99=NA'),
q653=recode(q06_5_3,'99=NA'),
q654=recode(q06_5_4,'99=NA'),
q655=recode(q06_5_5,'99=NA'),
q656=recode(q06_5_6,'99=NA'),
q657=recode(q06_5_7,'99=NA'))欠損値へのリコードがうまくできているか確認します。
table(data1$q06_5_1, data1$q651, useNA = "ifany")##
## 1 2 3 4 5 6 7 8 9 10 <NA>
## 1 408 0 0 0 0 0 0 0 0 0 0
## 2 0 170 0 0 0 0 0 0 0 0 0
## 3 0 0 254 0 0 0 0 0 0 0 0
## 4 0 0 0 236 0 0 0 0 0 0 0
## 5 0 0 0 0 428 0 0 0 0 0 0
## 6 0 0 0 0 0 282 0 0 0 0 0
## 7 0 0 0 0 0 0 139 0 0 0 0
## 8 0 0 0 0 0 0 0 168 0 0 0
## 9 0 0 0 0 0 0 0 0 95 0 0
## 10 0 0 0 0 0 0 0 0 0 294 0
## 99 0 0 0 0 0 0 0 0 0 0 273table(data1$q06_5_2, data1$q652, useNA = "ifany")##
## 1 2 3 4 5 6 7 8 9 10 <NA>
## 1 107 0 0 0 0 0 0 0 0 0 0
## 2 0 33 0 0 0 0 0 0 0 0 0
## 3 0 0 85 0 0 0 0 0 0 0 0
## 4 0 0 0 127 0 0 0 0 0 0 0
## 5 0 0 0 0 291 0 0 0 0 0 0
## 6 0 0 0 0 0 197 0 0 0 0 0
## 7 0 0 0 0 0 0 163 0 0 0 0
## 8 0 0 0 0 0 0 0 285 0 0 0
## 9 0 0 0 0 0 0 0 0 226 0 0
## 10 0 0 0 0 0 0 0 0 0 1078 0
## 99 0 0 0 0 0 0 0 0 0 0 155table(data1$q06_5_3, data1$q653, useNA = "ifany")##
## 1 2 3 4 5 6 7 8 9 10 <NA>
## 1 189 0 0 0 0 0 0 0 0 0 0
## 2 0 113 0 0 0 0 0 0 0 0 0
## 3 0 0 211 0 0 0 0 0 0 0 0
## 4 0 0 0 228 0 0 0 0 0 0 0
## 5 0 0 0 0 537 0 0 0 0 0 0
## 6 0 0 0 0 0 389 0 0 0 0 0
## 7 0 0 0 0 0 0 240 0 0 0 0
## 8 0 0 0 0 0 0 0 285 0 0 0
## 9 0 0 0 0 0 0 0 0 117 0 0
## 10 0 0 0 0 0 0 0 0 0 237 0
## 99 0 0 0 0 0 0 0 0 0 0 201table(data1$q06_5_4, data1$q654, useNA = "ifany")##
## 1 2 3 4 5 6 7 8 9 10 <NA>
## 1 457 0 0 0 0 0 0 0 0 0 0
## 2 0 200 0 0 0 0 0 0 0 0 0
## 3 0 0 311 0 0 0 0 0 0 0 0
## 4 0 0 0 246 0 0 0 0 0 0 0
## 5 0 0 0 0 534 0 0 0 0 0 0
## 6 0 0 0 0 0 349 0 0 0 0 0
## 7 0 0 0 0 0 0 199 0 0 0 0
## 8 0 0 0 0 0 0 0 157 0 0 0
## 9 0 0 0 0 0 0 0 0 50 0 0
## 10 0 0 0 0 0 0 0 0 0 81 0
## 99 0 0 0 0 0 0 0 0 0 0 163table(data1$q06_5_5, data1$q655, useNA = "ifany")##
## 1 2 3 4 5 6 7 8 9 10 <NA>
## 1 600 0 0 0 0 0 0 0 0 0 0
## 2 0 344 0 0 0 0 0 0 0 0 0
## 3 0 0 426 0 0 0 0 0 0 0 0
## 4 0 0 0 312 0 0 0 0 0 0 0
## 5 0 0 0 0 409 0 0 0 0 0 0
## 6 0 0 0 0 0 193 0 0 0 0 0
## 7 0 0 0 0 0 0 67 0 0 0 0
## 8 0 0 0 0 0 0 0 83 0 0 0
## 9 0 0 0 0 0 0 0 0 47 0 0
## 10 0 0 0 0 0 0 0 0 0 85 0
## 99 0 0 0 0 0 0 0 0 0 0 181table(data1$q06_5_6, data1$q656, useNA = "ifany")##
## 1 2 3 4 5 6 7 8 9 10 <NA>
## 1 87 0 0 0 0 0 0 0 0 0 0
## 2 0 46 0 0 0 0 0 0 0 0 0
## 3 0 0 89 0 0 0 0 0 0 0 0
## 4 0 0 0 97 0 0 0 0 0 0 0
## 5 0 0 0 0 303 0 0 0 0 0 0
## 6 0 0 0 0 0 226 0 0 0 0 0
## 7 0 0 0 0 0 0 157 0 0 0 0
## 8 0 0 0 0 0 0 0 278 0 0 0
## 9 0 0 0 0 0 0 0 0 272 0 0
## 10 0 0 0 0 0 0 0 0 0 999 0
## 99 0 0 0 0 0 0 0 0 0 0 193table(data1$q06_5_7, data1$q657, useNA = "ifany")##
## 1 2 3 4 5 6 7 8 9 10 <NA>
## 1 270 0 0 0 0 0 0 0 0 0 0
## 2 0 133 0 0 0 0 0 0 0 0 0
## 3 0 0 249 0 0 0 0 0 0 0 0
## 4 0 0 0 262 0 0 0 0 0 0 0
## 5 0 0 0 0 487 0 0 0 0 0 0
## 6 0 0 0 0 0 316 0 0 0 0 0
## 7 0 0 0 0 0 0 226 0 0 0 0
## 8 0 0 0 0 0 0 0 248 0 0 0
## 9 0 0 0 0 0 0 0 0 122 0 0
## 10 0 0 0 0 0 0 0 0 0 257 0
## 99 0 0 0 0 0 0 0 0 0 0 177合計して価値観のスコアを作成します。
attach(data1)
data1$value=q651+q652+q653+q654+q655+q656+q657独立変数の準備をします。
data1 <- data1 %>% mutate(
age=q10_5,
fem=recode(gender,'2=1; 1=0'),
mar=recode(q10_7,'1=0; 2=1; 3:9=NA'),
edu=as.factor(recode(q10_12,'1:2=1; 3:4=2; 5:6=3; 7:9=NA')),
wrk=as.factor(recode(q10_11,'1:3=1; 4:7=2; 9:10=3; 8=NA; 11:99=NA')))リコードがうまくできているか確認します。
table(data1$gender, data1$fem, useNA = "ifany")##
## 0 1
## 1 1556 0
## 2 0 1191table(data1$q10_7, data1$mar, useNA = "ifany")##
## 0 1 <NA>
## 1 892 0 0
## 2 0 1734 0
## 3 0 0 93
## 9 0 0 28table(data1$q10_12, data1$edu, useNA = "ifany")##
## 1 2 3 <NA>
## 1 25 0 0 0
## 2 597 0 0 0
## 3 0 177 0 0
## 4 0 327 0 0
## 5 0 0 1391 0
## 6 0 0 184 0
## 7 0 0 0 12
## 9 0 0 0 34table(data1$q10_11, data1$wrk, useNA = "ifany")##
## 1 2 3 <NA>
## 1 77 0 0 0
## 2 329 0 0 0
## 3 773 0 0 0
## 4 0 303 0 0
## 5 0 65 0 0
## 6 0 207 0 0
## 7 0 20 0 0
## 8 0 0 0 113
## 9 0 0 435 0
## 10 0 0 335 0
## 11 0 0 0 49
## 99 0 0 0 41reg1 <- lm(value~age+fem+mar+edu+wrk, data=data1)
summary(reg1)##
## Call:
## lm(formula = value ~ age + fem + mar + edu + wrk, data = data1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -35.508 -7.064 0.404 7.500 34.984
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 31.93944 1.15728 27.599 < 2e-16 ***
## age 0.12878 0.02248 5.729 1.16e-08 ***
## fem 0.45071 0.57581 0.783 0.43387
## mar 2.74634 0.61545 4.462 8.56e-06 ***
## edu2 -0.58234 0.79319 -0.734 0.46293
## edu3 -2.09145 0.64087 -3.263 0.00112 **
## wrk2 0.31086 0.67057 0.464 0.64300
## wrk3 1.92950 0.67605 2.854 0.00436 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 11.22 on 2012 degrees of freedom
## ( 727 個の観測値が欠損のため削除されました )
## Multiple R-squared: 0.07397, Adjusted R-squared: 0.07075
## F-statistic: 22.96 on 7 and 2012 DF, p-value: < 2.2e-16年齢が高い人、既婚者、無職の人は、いずれも1％水準で有意にスコアが高い（より保守的）ことがわかります。 年齢が1歳上がると、スコアが0.13上がり、既婚者は未婚者に比べてスコアが2.75高く、無職の人は正規職の人に比べて、スコアが1.93高くなっています。 また、学歴が高い人は、1％水準で有意にスコアが低い（よりリベラル）であることもわかります。高校以下に比べて、大学・大学院は2.09低くなっています。 性別による有意差はないようです。
data1 <- data1 %>% mutate(
agecl=as.factor(recode(age,'20:29=1;30:39=2;40:49=3;50:59=4;60:69=5')))
table(data1$age, data1$agecl, useNA="ifany")##
## 1 2 3 4 5
## 20 20 0 0 0 0
## 21 24 0 0 0 0
## 22 35 0 0 0 0
## 23 33 0 0 0 0
## 24 35 0 0 0 0
## 25 48 0 0 0 0
## 26 61 0 0 0 0
## 27 62 0 0 0 0
## 28 69 0 0 0 0
## 29 69 0 0 0 0
## 30 0 42 0 0 0
## 31 0 42 0 0 0
## 32 0 48 0 0 0
## 33 0 53 0 0 0
## 34 0 40 0 0 0
## 35 0 56 0 0 0
## 36 0 60 0 0 0
## 37 0 90 0 0 0
## 38 0 66 0 0 0
## 39 0 88 0 0 0
## 40 0 0 55 0 0
## 41 0 0 48 0 0
## 42 0 0 46 0 0
## 43 0 0 44 0 0
## 44 0 0 41 0 0
## 45 0 0 40 0 0
## 46 0 0 52 0 0
## 47 0 0 46 0 0
## 48 0 0 47 0 0
## 49 0 0 38 0 0
## 50 0 0 0 55 0
## 51 0 0 0 68 0
## 52 0 0 0 67 0
## 53 0 0 0 52 0
## 54 0 0 0 58 0
## 55 0 0 0 49 0
## 56 0 0 0 46 0
## 57 0 0 0 38 0
## 58 0 0 0 54 0
## 59 0 0 0 36 0
## 60 0 0 0 0 83
## 61 0 0 0 0 88
## 62 0 0 0 0 80
## 63 0 0 0 0 78
## 64 0 0 0 0 94
## 65 0 0 0 0 85
## 66 0 0 0 0 79
## 67 0 0 0 0 62
## 68 0 0 0 0 49
## 69 0 0 0 0 28データフレームに変数を追加せずに、２乗した年齢をI()で独立変数に加えます。
reg2 <- lm(value~age+I(age^2)+fem+mar+edu+wrk, data=data1)
summary(reg2)##
## Call:
## lm(formula = value ~ age + I(age^2) + fem + mar + edu + wrk,
## data = data1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -35.482 -6.894 0.389 7.624 33.621
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 46.326559 3.559080 13.016 < 2e-16 ***
## age -0.550426 0.160540 -3.429 0.000619 ***
## I(age^2) 0.007288 0.001706 4.273 2.02e-05 ***
## fem 0.575986 0.574107 1.003 0.315849
## mar 3.189671 0.621554 5.132 3.15e-07 ***
## edu2 -0.346384 0.791740 -0.437 0.661798
## edu3 -2.153229 0.638305 -3.373 0.000757 ***
## wrk2 0.098465 0.669556 0.147 0.883100
## wrk3 1.231315 0.692717 1.778 0.075634 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 11.18 on 2011 degrees of freedom
## ( 727 個の観測値が欠損のため削除されました )
## Multiple R-squared: 0.0823, Adjusted R-squared: 0.07865
## F-statistic: 22.54 on 8 and 2011 DF, p-value: < 2.2e-16reg3 <- lm(value~agecl+fem+mar+edu+wrk, data=data1)
summary(reg3)##
## Call:
## lm(formula = value ~ agecl + fem + mar + edu + wrk, data = data1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -36.534 -6.937 0.339 7.626 34.035
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 37.1161 0.9375 39.591 < 2e-16 ***
## agecl2 -1.4005 0.8924 -1.569 0.116721
## agecl3 -0.4357 0.9610 -0.453 0.650346
## agecl4 0.4037 0.9614 0.420 0.674563
## agecl5 3.8007 0.9814 3.873 0.000111 ***
## fem 0.3932 0.5725 0.687 0.492306
## mar 3.2846 0.6161 5.331 1.08e-07 ***
## edu2 -0.4076 0.7925 -0.514 0.607037
## edu3 -2.1432 0.6400 -3.349 0.000826 ***
## wrk2 0.1651 0.6697 0.246 0.805336
## wrk3 1.4759 0.6875 2.147 0.031934 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 11.19 on 2009 degrees of freedom
## ( 727 個の観測値が欠損のため削除されました )
## Multiple R-squared: 0.08123, Adjusted R-squared: 0.07665
## F-statistic: 17.76 on 10 and 2009 DF, p-value: < 2.2e-16stargazer(reg1, reg2, reg3, type="text", star.cutoffs = c(0.1, 0.05, 0.01), keep.stat=c("n","rsq","adj.rsq"))##
## ==========================================
## Dependent variable:
## -----------------------------
## value
## (1) (2) (3)
## ------------------------------------------
## age 0.129*** -0.550***
## (0.022) (0.161)
##
## I(age2) 0.007***
## (0.002)
##
## agecl2 -1.401
## (0.892)
##
## agecl3 -0.436
## (0.961)
##
## agecl4 0.404
## (0.961)
##
## agecl5 3.801***
## (0.981)
##
## fem 0.451 0.576 0.393
## (0.576) (0.574) (0.572)
##
## mar 2.746*** 3.190*** 3.285***
## (0.615) (0.622) (0.616)
##
## edu2 -0.582 -0.346 -0.408
## (0.793) (0.792) (0.792)
##
## edu3 -2.091*** -2.153*** -2.143***
## (0.641) (0.638) (0.640)
##
## wrk2 0.311 0.098 0.165
## (0.671) (0.670) (0.670)
##
## wrk3 1.930*** 1.231* 1.476**
## (0.676) (0.693) (0.688)
##
## Constant 31.939*** 46.327*** 37.116***
## (1.157) (3.559) (0.937)
##
## ------------------------------------------
## Observations 2,020 2,020 2,020
## R2 0.074 0.082 0.081
## Adjusted R2 0.071 0.079 0.077
## ==========================================
## Note: *p<0.1; **p<0.05; ***p<0.01年齢の2次項の係数が正ですので、年齢と価値観の関係は、下に凸の形状になっていることがわかります。頂点の年齢を求めると、39.3歳くらいになることがわかります。 年齢階級ダミーでみると、20代に比べて30代、40代はややスコアが低く、50代でやや高く、60代は有意に高いことからも、その形状であることが確認できます。
reg1m <- lm(value~age+mar+edu+wrk, data=data1, subset=fem==0)
summary(reg1m)##
## Call:
## lm(formula = value ~ age + mar + edu + wrk, data = data1, subset = fem ==
## 0)
##
## Residuals:
## Min 1Q Median 3Q Max
## -36.245 -7.308 0.352 7.803 34.254
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 30.80016 1.56290 19.707 < 2e-16 ***
## age 0.12701 0.03377 3.761 0.000178 ***
## mar 3.64978 0.87772 4.158 3.44e-05 ***
## edu2 0.86390 1.27645 0.677 0.498670
## edu3 -1.77824 0.86440 -2.057 0.039891 *
## wrk2 0.77950 0.91426 0.853 0.394054
## wrk3 2.69877 1.03934 2.597 0.009534 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 11.66 on 1164 degrees of freedom
## ( 385 個の観測値が欠損のため削除されました )
## Multiple R-squared: 0.07997, Adjusted R-squared: 0.07523
## F-statistic: 16.86 on 6 and 1164 DF, p-value: < 2.2e-16reg1f <- lm(value~age+mar+edu+wrk, data=data1, subset=fem==1)
summary(reg1f)##
## Call:
## lm(formula = value ~ age + mar + edu + wrk, data = data1, subset = fem ==
## 1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -33.334 -6.973 0.638 6.947 33.658
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 34.30946 1.60653 21.356 < 2e-16 ***
## age 0.11871 0.03104 3.824 0.000141 ***
## mar 1.81017 0.94710 1.911 0.056310 .
## edu2 -1.50699 1.00793 -1.495 0.135254
## edu3 -2.57994 0.95973 -2.688 0.007326 **
## wrk2 -0.11743 1.04259 -0.113 0.910351
## wrk3 1.56673 0.99268 1.578 0.114876
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 10.59 on 842 degrees of freedom
## ( 342 個の観測値が欠損のため削除されました )
## Multiple R-squared: 0.0677, Adjusted R-squared: 0.06106
## F-statistic: 10.19 on 6 and 842 DF, p-value: 6.654e-11stargazer(reg1m, reg1f, type="text", star.cutoffs = c(0.1, 0.05, 0.01), keep.stat=c("n","rsq","adj.rsq"))##
## =========================================
## Dependent variable:
## ----------------------------
## value
## (1) (2)
## -----------------------------------------
## age 0.127*** 0.119***
## (0.034) (0.031)
##
## mar 3.650*** 1.810*
## (0.878) (0.947)
##
## edu2 0.864 -1.507
## (1.276) (1.008)
##
## edu3 -1.778** -2.580***
## (0.864) (0.960)
##
## wrk2 0.779 -0.117
## (0.914) (1.043)
##
## wrk3 2.699*** 1.567
## (1.039) (0.993)
##
## Constant 30.800*** 34.309***
## (1.563) (1.607)
##
## -----------------------------------------
## Observations 1,171 849
## R2 0.080 0.068
## Adjusted R2 0.075 0.061
## =========================================
## Note: *p<0.1; **p<0.05; ***p<0.01男女別に推定したところ、男性のほうが、既婚のスコア上昇効果が大きいように見えます。 また、無職であることの効果が有意なのは男性だけであることもわかります。
reg4 <- lm(value~age+fem*mar+edu+wrk, data=data1)
summary(reg4)##
## Call:
## lm(formula = value ~ age + fem * mar + edu + wrk, data = data1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -35.935 -7.124 0.534 7.539 35.292
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 31.52316 1.17728 26.776 < 2e-16 ***
## age 0.12290 0.02268 5.419 6.70e-08 ***
## fem 1.87191 0.94622 1.978 0.048030 *
## mar 3.65229 0.77946 4.686 2.98e-06 ***
## edu2 -0.53197 0.79313 -0.671 0.502476
## edu3 -2.13009 0.64079 -3.324 0.000903 ***
## wrk2 0.52149 0.67932 0.768 0.442779
## wrk3 2.26964 0.69912 3.246 0.001188 **
## fem:mar -2.20256 1.16409 -1.892 0.058623 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 11.22 on 2011 degrees of freedom
## ( 727 個の観測値が欠損のため削除されました )
## Multiple R-squared: 0.07562, Adjusted R-squared: 0.07194
## F-statistic: 20.56 on 8 and 2011 DF, p-value: < 2.2e-16女性ダミーと既婚ダミーの交差項の係数の推定が負で有意となっています。これは11-4で見たように、女性の場合、男性に比べて既婚のスコア上昇効果が小さいことが確認できました。
女性未婚（femsin）、女性既婚（femmar）、男性未婚（malsin）、男性既婚（malmar）の各ダミーを作成します。
data1 <- data1 %>% mutate(
femsin=case_when(fem==1 & mar==0 ~ 1, fem!=1 | mar!=0 ~ 0),
femmar=case_when(fem==1 & mar==1 ~ 1, fem!=1 | mar!=1 ~ 0),
malsin=case_when(fem==0 & mar==0 ~ 1, fem!=0 | mar!=0 ~ 0),
malmar=case_when(fem==0 & mar==1 ~ 1, fem!=0 | mar!=1 ~ 0),)table(data1$femsin,data1$fem,useNA="ifany")##
## 0 1
## 0 1556 778
## 1 0 349
## <NA> 0 64table(data1$femsin,data1$mar,useNA="ifany")##
## 0 1 <NA>
## 0 543 1734 57
## 1 349 0 0
## <NA> 0 0 64table(data1$femmar,data1$fem,useNA="ifany")##
## 0 1
## 0 1556 349
## 1 0 778
## <NA> 0 64table(data1$femmar,data1$mar,useNA="ifany")##
## 0 1 <NA>
## 0 892 956 57
## 1 0 778 0
## <NA> 0 0 64table(data1$malsin,data1$fem,useNA="ifany")##
## 0 1
## 0 956 1191
## 1 543 0
## <NA> 57 0table(data1$malsin,data1$mar,useNA="ifany")##
## 0 1 <NA>
## 0 349 1734 64
## 1 543 0 0
## <NA> 0 0 57table(data1$malmar,data1$fem,useNA="ifany")##
## 0 1
## 0 543 1191
## 1 956 0
## <NA> 57 0table(data1$malmar,data1$mar,useNA="ifany")##
## 0 1 <NA>
## 0 892 778 64
## 1 0 956 0
## <NA> 0 0 57attach(data1)## 以下のオブジェクトは data1 (pos = 3) からマスクされています:
##
## gender, q04_1, q06_5_1, q06_5_2, q06_5_3, q06_5_4, q06_5_5,
## q06_5_6, q06_5_7, q10_1, q10_11, q10_12, q10_5, q10_7, q651, q652,
## q653, q654, q655, q656, q657data2 <- data.frame(fem,mar,femsin,femmar,malsin,malmar)
head(data2,n=10)## fem mar femsin femmar malsin malmar
## 1 0 1 0 0 0 1
## 2 0 1 0 0 0 1
## 3 1 0 1 0 0 0
## 4 1 1 0 1 0 0
## 5 0 1 0 0 0 1
## 6 0 0 0 0 1 0
## 7 1 0 1 0 0 0
## 8 1 0 1 0 0 0
## 9 1 0 1 0 0 0
## 10 0 1 0 0 0 1reg5 <- lm(value~age+femmar+femsin+malmar+edu+wrk,data=data1)
stargazer(reg4, reg5, type="text",star.cutoffs=c(0.1,0.05,0.01), keep.stat=c("n","rsq","adj.rsq"))##
## =========================================
## Dependent variable:
## ----------------------------
## value
## (1) (2)
## -----------------------------------------
## age 0.123*** 0.123***
## (0.023) (0.023)
##
## fem 1.872**
## (0.946)
##
## mar 3.652***
## (0.779)
##
## femmar 3.322***
## (0.796)
##
## femsin 1.872**
## (0.946)
##
## malmar 3.652***
## (0.779)
##
## edu2 -0.532 -0.532
## (0.793) (0.793)
##
## edu3 -2.130*** -2.130***
## (0.641) (0.641)
##
## wrk2 0.521 0.521
## (0.679) (0.679)
##
## wrk3 2.270*** 2.270***
## (0.699) (0.699)
##
## fem:mar -2.203*
## (1.164)
##
## Constant 31.523*** 31.523***
## (1.177) (1.177)
##
## -----------------------------------------
## Observations 2,020 2,020
## R2 0.076 0.076
## Adjusted R2 0.072 0.072
## =========================================
## Note: *p<0.1; **p<0.05; ***p<0.01組み合わせダミーから以下のように、交差項の係数の推定値が得られます（四捨五入されているので値は少しずれています）。
(malmar-malsin)-(femmar-femsin)=(3.652-0)-(3.322-1.872)=2.202
また、性別と既婚ダミーの交差項を使用した推定結果と見比べると、いくつか共通の推定値があることがわかります。
組み合わせダミーのベースカテゴリーを色々変えてみた結果が以下です。
reg6 <- lm(value~age+femmar+femsin+malsin+edu+wrk,data=data1)
reg7 <- lm(value~age+femmar+malmar+malsin+edu+wrk,data=data1)
reg8 <- lm(value~age+femsin+malmar+malsin+edu+wrk,data=data1)
stargazer(reg5, reg6, reg7, reg8, type="text",star.cutoffs=c(0.1,0.05,0.01), keep.stat=c("n","rsq","adj.rsq"))##
## ====================================================
## Dependent variable:
## ---------------------------------------
## value
## (1) (2) (3) (4)
## ----------------------------------------------------
## age 0.123*** 0.123*** 0.123*** 0.123***
## (0.023) (0.023) (0.023) (0.023)
##
## femmar 3.322*** -0.331 1.450
## (0.796) (0.708) (0.921)
##
## femsin 1.872** -1.780* -1.450
## (0.946) (0.930) (0.921)
##
## malmar 3.652*** 1.780* 0.331
## (0.779) (0.930) (0.708)
##
## malsin -3.652*** -1.872** -3.322***
## (0.779) (0.946) (0.796)
##
## edu2 -0.532 -0.532 -0.532 -0.532
## (0.793) (0.793) (0.793) (0.793)
##
## edu3 -2.130*** -2.130*** -2.130*** -2.130***
## (0.641) (0.641) (0.641) (0.641)
##
## wrk2 0.521 0.521 0.521 0.521
## (0.679) (0.679) (0.679) (0.679)
##
## wrk3 2.270*** 2.270*** 2.270*** 2.270***
## (0.699) (0.699) (0.699) (0.699)
##
## Constant 31.523*** 35.175*** 33.395*** 34.845***
## (1.177) (1.367) (1.226) (1.276)
##
## ----------------------------------------------------
## Observations 2,020 2,020 2,020 2,020
## R2 0.076 0.076 0.076 0.076
## Adjusted R2 0.072 0.072 0.072 0.072
## ====================================================
## Note: *p<0.1; **p<0.05; ***p<0.01