R

sem_PLSモデル

sem,qgraph,Rgraphviz

(参考)
はやしさんのブログ Rev.3
名城大学人間学部神谷研究室

(準備)　graphvizをインストールする

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
library(sem)
#library(semPlot)
library(qgraph)
library(knitr)
#blue backs「原因をさぐる統計学」p.88 データ
data <- read.table("http://wwwhum.meijo-u.ac.jp/labs/hh002/r/rfile/gan.csv", header=TRUE, sep=",", na.strings="NA", dec=".", 
strip.white=TRUE)
names(data)<-c("総熱量","肉類","乳製品","酒類","大腸ガン","直腸ガン")
soukan<-cor(data)
 
#covsオプションに"変数1, 変数2"という形で共分散を指定
#endog.variances=F (洋食傾向に誤差項をつけない)
model.c1 <- specifyEquations(endog.variances=F,
 covs="総熱量, 肉類,乳製品,酒類")
洋食傾向 = gam1*総熱量 + gam2*肉類+ gam3*乳製品 + gam4*酒類
model.c2 <- specifyModel()
消化管ガン -> 大腸ガン, lam61, NA
消化管ガン -> 直腸ガン, lam62, NA
消化管ガン <-> 消化管ガン, NA, 1
洋食傾向 -> 消化管ガン, NA, 1
洋食傾向 <-> 洋食傾向,NA,0
#モデルを結合
model.c<- combineModels(model.c1, model.c2)
#標準的なsemmodクラスに変換
model <- removeRedundantPaths(model.c)
fit <- sem:::sem(model,soukan, N=nrow(data))
kable(stdCoef(fit))
Std. Estimate
gam1 0.3449044 洋食傾向 <—- 総熱量
gam2 0.6735523 洋食傾向 <—- 肉類
gam3 -0.1391173 洋食傾向 <—- 乳製品
gam4 0.2119839 洋食傾向 <—- 酒類
V[総熱量] 1.0000000 総熱量 <—> 総熱量
C[総熱量,肉類] 0.7731110 肉類 <—> 総熱量
C[総熱量,乳製品] 0.7155246 乳製品 <—> 総熱量
C[総熱量,酒類] 0.6341239 酒類 <—> 総熱量
V[肉類] 1.0000000 肉類 <—> 肉類
C[肉類,乳製品] 0.7227996 乳製品 <—> 肉類
C[肉類,酒類] 0.4515932 酒類 <—> 肉類
V[乳製品] 1.0000000 乳製品 <—> 乳製品
C[乳製品,酒類] 0.2321328 酒類 <—> 乳製品
V[酒類] 1.0000000 酒類 <—> 酒類
lam61 0.9094601 大腸ガン <—- 消化管ガン
lam62 0.8272908 直腸ガン <—- 消化管ガン
0.0775069 消化管ガン <—> 消化管ガン
0.9604650 消化管ガン <—- 洋食傾向
0.0000000 洋食傾向 <—> 洋食傾向
V[大腸ガン] 0.1728824 大腸ガン <—> 大腸ガン
V[直腸ガン] 0.3155899 直腸ガン <—> 直腸ガン 
1
summary(fit)

Model Chisquare = 17.34207 Df = 3 Pr(>Chisq) = 0.0006010001
AIC = 53.34207
BIC = 5.791629

Normalized Residuals
Min. 1st Qu. Median Mean 3rd Qu. Max.
-0.3885000 -0.0000015 0.0000000 0.0068180 0.0000014 0.5215000

R-square for Endogenous Variables
洋食傾向 消化管ガン 大腸ガン 直腸ガン
1.0000 0.9225 0.8271 0.6844

Parameter Estimates
Estimate Std Error z value Pr(>|z|)
gam1 1.1898995 0.69655784 1.708257 8.758874e-02 洋食傾向 <—- 総熱量
gam2 2.3237115 0.98419115 2.361037 1.822393e-02 洋食傾向 <—- 肉類
gam3 -0.4799463 0.46171648 -1.039483 2.985803e-01 洋食傾向 <—- 乳製品
gam4 0.7313307 0.45144880 1.619964 1.052401e-01 洋食傾向 <—- 酒類
V[総熱量] 0.9999997 0.20851434 4.795832 1.620014e-06 総熱量 <—> 総熱量
C[総熱量,肉類] 0.7731107 0.18636684 4.148327 3.349131e-05 肉類 <—> 総熱量
C[総熱量,乳製品] 0.7155249 0.18129821 3.946674 7.924433e-05 乳製品 <—> 総熱量
C[総熱量,酒類] 0.6341241 0.17458730 3.632132 2.810895e-04 酒類 <—> 総熱量
V[肉類] 0.9999998 0.20851438 4.795832 1.620014e-06 肉類 <—> 肉類
C[肉類,乳製品] 0.7227997 0.18192446 3.973076 7.095045e-05 乳製品 <—> 肉類
C[肉類,酒類] 0.4515934 0.16177929 2.791417 5.247789e-03 酒類 <—> 肉類
V[乳製品] 1.0000003 0.20851447 4.795832 1.620014e-06 乳製品 <—> 乳製品
C[乳製品,酒類] 0.2321335 0.15136241 1.533627 1.251214e-01 酒類 <—> 乳製品
V[酒類] 1.0000004 0.20851449 4.795832 1.620014e-06 酒類 <—> 酒類
lam61 0.2531944 0.09773497 2.590623 9.580246e-03 大腸ガン <—- 消化管ガン
lam62 0.2303184 0.08890465 2.590623 9.580246e-03 直腸ガン <—- 消化管ガン
V[大腸ガン] 0.1728823 0.05975160 2.893350 3.811566e-03 大腸ガン <—> 大腸ガン
V[直腸ガン] 0.3155902 0.07671436 4.113834 3.891411e-05 直腸ガン <—> 直腸ガン

Iterations = 51 

1
qgraph(fit,edge.labels=T,layout ="tree",vsize.man=6,vsize.lat=8,fade=F)

1
2
3
4
library(Rgraphviz)
pathDiagram(fit,edge.labels="values",ignore.double=T,standardize=T,digits=2,file="sem3",graphics.fmt ="png",
node.font = c("serif", 14), edge.font = c("sans", 10),
min.rank="総熱量, 肉類, 乳製品, 酒類",max.rank="大腸ガン, 直腸ガン",same.rank=NULL)