共分散行列を利用

lavaanパッケージ　gawk graphviz

(参考)blue backs「原因をさぐる統計学」

OSはlinux(zorinOS)
windowsではlavaan:::sem使用時に変数名が日本語だとエラーになりました。

(準備)「因子分析＆確証的因子分析」にあるmakeawk関数を読み込んでおく。

「原因をさぐる統計学」p.127

アメリカ　オハイオ州で行われたヘッドスタートの評価データの分析

library(lavaan)
#Sample Size は 1800
mat<-'
1
0.484 1
0.224 0.342 1
0.268 0.215 0.387 1
0.23 0.215 0.196 0.115 1
0.265 0.297 0.234 0.162 0.635 1
'
wheaton.cov <- getCov(mat,names=c("母親の学歴","父親の学歴","父親の職業","家庭の収入","MRT","ITPA"))
# 多重指標モデル
model<-'
社会的地位=~母親の学歴+父親の学歴+父親の職業+家庭の収入
子供の知能=~MRT+ITPA
子供の知能~~社会的地位
'
fit <- sem(model, sample.cov=wheaton.cov, sample.nobs=1800)
summary(fit, standardized=TRUE)
latent<-c("社会的地位","子供の知能")
makeawk(latent)
pars <- parameterEstimates(fit, standardized = TRUE)
g<-data.frame(pars$lhs,pars$op,pars$rhs," [label=",signif(pars$std.all,digits=2),"];")
write.table(g,"out",sep=" ",col.names=FALSE,row.names=FALSE,quote=FALSE, na="")
system( 'gawk -f "lavaan_dot.awk" "out">out.dot')

summary(fit, standardized=TRUE)の出力は省略

out.dotを編集

system('dot -Tpng out.dot -o out.png')

mat01.dot

# 外生的潜在変数が２つのモデル
model<-'
両親の学歴=~母親の学歴+父親の学歴
家庭の経済状態=~父親の職業+家庭の収入
子供の知能=~MRT+ITPA
子供の知能~両親の学歴+家庭の経済状態
両親の学歴~~家庭の経済状態
'
fit <- sem(model, sample.cov=wheaton.cov, sample.nobs=1800)
summary(fit, standardized=TRUE)
latent<-c("両親の学歴","家庭の経済状態","子供の知能")
makeawk(latent)
pars <- parameterEstimates(fit, standardized = TRUE)
g<-data.frame(pars$lhs,pars$op,pars$rhs," [label=",signif(pars$std.all,digits=2),"];")
write.table(g,"out",sep=" ",col.names=FALSE,row.names=FALSE,quote=FALSE, na="")
system( 'gawk -f "lavaan_dot.awk" "out">out.dot')

summary(fit, standardized=TRUE)の出力は省略

絶対位置を指定する。

pngへの変換にはdotではなくneatoを使う。

system('neato -Tpng out.dot -o out.png')

mat02.dot

「原因をさぐる統計学」p.151

カップルデータの因果モデルの分析

ノードのラベルが長いのでdotファイルにしてからラベルをつけます。（ノードに改行を入れる必要あり）

# 外生的潜在変数が２つのモデル
# データ数は340
mat<-'
1
0.47 1
0.46 0.27 1
0.312 0.223 0.495 1
0.628 0.421 0.498 0.381 1
0.596 0.347 0.586 0.422 0.816 1
'
wheaton.cov <- getCov(mat,names=c("X1","X2","X3","X4","X5","X6"))
model<-'
F1=~X1+X2
F2=~X3+X4
F3=~X5+X6
F3~F1+F2
F1~~F2
'
fit <- sem(model, sample.cov=wheaton.cov, sample.nobs=340)
summary(fit, standardized=TRUE)
latent<-c("F1","F2","F3")
makeawk(latent)
pars <- parameterEstimates(fit, standardized = TRUE)
g<-data.frame(pars$lhs,pars$op,pars$rhs," [label=",signif(pars$std.all,digits=2),"];")
write.table(g,"out",sep=" ",col.names=FALSE,row.names=FALSE,quote=FALSE, na="")
system( 'gawk -f "lavaan_dot.awk" "out">out.dot')

summary(fit, standardized=TRUE)の出力は省略

絶対位置を指定する。
ノードにラベルをつける

system('neato -Tpng out.dot -o out.png')

mat03.dot