<< トップページへ

<< 目次へ

最終更新日: 2009.12.17

ごく簡単な説明(準備偏)



1 semパッケージのインストール

 Rで共分散構造分析を行うためにはsemパッケージをインストールする必要がある。また、使用の際にはライブラリからsemパッケージを呼び出さなければならない。

> install.packages("sem")   # semパッケージのダウンロード・インストール
> library(sem)              # semパッケージをライブラリから呼び出す


2 sem()の引数

 基本的にsem()に指定すべき引数は3つある。1つは解析しようとするモデルが格納されているオブジェクト、もう1つは分散・共分散行列(あるいは相関行列でもよい)が格納されているオブジェクト、そしてサンプルサイズである。

> sem(MODEL, COVMAT, N)

2.1 MODELオブジェクトの用意

 semパッケージに含まれているspecify.model()を使うことになるが、基本的にはテキストファイルを読み込むのが簡単である。例えば、重回帰モデル:Y = X + Zを解析するならば、次のようなテキストファイルを用意する。なお、最初の行にある「--- ここから ---」という部分と、最終行にある「--- ここまで ---」というのは私が分かりやすく明示するために書いたもので、テキストファイルにはこの2つの行を書く必要はない。

--- MODEL.txt ここから ---
Y <- X, beta1, NA
Y <- Z, beta2, NA
X <-> Z, gamma, NA
X <-> X, sita1, NA
Y <-> Y, sita2, NA
Z <-> Z, sita3, NA
--- MODEL.txt ここまで ---

 上の命令について、第1列にはモデル式を、第2列にはパラメータを、第3列にはパラメータの初期値を書いている。例えば、第1行目の命令は「Y <- X, beta1, NA」となっている。この行に含まれる「Y <- X」というのはY = Xという回帰モデルを指定していることになる。また、第2列目にあるbeta1というのはパラメータの名前で、これはY = βXという回帰モデルの説明変数Xの係数βであることを意味している。これは自分が分かりやすいような、好きな名前を付けてよい。そして第3列目の「NA」というのは、パラメータの初期値で特に指定する必要のないときはNAと書いておけばよいのである。もし初期値を1にするならば、NAの代わりに1を入力しておけばよい。

 ところで、Y <- Xという一方向矢印は回帰関係を表しているが、X <-> Zのように双方向矢印は相関関係を表している。ただし、X <-> Xのように同じ変数同士を両方向矢印で結んだ場合は、誤差項を表している。Rでは回帰関係、相関関係、誤差項の全てをMODELオブジェクトとして用意しなければならないので、注意が必要である(LISRELのSIMPLISでは単にモデルの指定(回帰関係)のみを指定すればよいが、Rではそうはいかないのが面倒である)。

2.2 COVMATオブジェクトの用意

 分散・共分散行列を作成するためにはvar()を用いればよいし、相関行列はcor()で作成することができる。出力結果の見易さを考慮して、分散・共分散行列、あるいは相関行列でも行ラベルと列ラベルをそれぞれrownames()とcolnames()を用いてラベリングしておいたほうがよいだろう。

> cor.mat <- matrix(c(
+ 1.0000000, 0         , 0,
+ 0.2505807, 1.0000000 , 0,
+ 0.4948428, 0.5980343 , 1.0000000), ncol=3, byrow=T)

> rownames(cor.mat) <- colnames(cor.mat) <- c("X", "Y", "Z")

2.3 実行例

 あらかじめ2.1で作成したMODEL.txtというテキストファイルがワークディレクトリに保存されているとする(getwd()で確認できる)。それから先ほど2.2で用意した相関行列がcor.matというオブジェクトに格納されているとする。このとき、sem()を用いて共分散構造分析を実行するためには次のようにコマンドすることになる。

> my.model <- specify.model("MODEL.txt")   # 解析するモデルをテキストファイルから読み込む
> result <- sem(my.model, cor.mat, N=15)   # サンプルサイズはN=15である
> summary(result)                          # summary()で結果を出力できる

 Model Chisquare =  -7.7716e-15   Df =  0 Pr(>Chisq) = NA
 Chisquare (null model) =  10.188   Df =  3
 Goodness-of-fit index =  1
 BIC =  -7.7716e-15 

 Normalized Residuals
     Min.   1st Qu.    Median      Mean   3rd Qu.      Max. 
-5.87e-16  0.00e+00  0.00e+00 -2.05e-17  0.00e+00  2.01e-16 

 Parameter Estimates
      Estimate  Std Error z value  Pr(>|z|)          
beta1 -0.060059 0.24597   -0.24417 0.8071013 Y <--- X
beta2  0.627754 0.24597    2.55211 0.0107072 Y <--- Z
gamma  0.494843 0.29815    1.65969 0.0969774 Z <--> X
sita1  1.000000 0.37811    2.64470 0.0081762 X <--> X
sita2  0.639631 0.24183    2.64492 0.0081709 Y <--> Y
sita3  1.000000 0.37811    2.64470 0.0081762 Z <--> Z

 Iterations =  0

3 パスダイアグラムの作図

 Graphvizというソフトウェアが必要で、これはhttp://www.graphviz.org/より無料でダウンロード・インストールして使用することができる。

3.1 R上でpath.diagram()を用いてdotファイルのソースコードを表示させる

# 数値出力の桁数は引数digitsで指定できる
> path.diagram(result, ignore.double=FALSE, edge.labels="values", digits=3)
digraph "sem.result" {
  rankdir=LR;
  size="8,8";
  node [fontname="Helvetica" fontsize=14 shape=box];
  edge [fontname="Helvetica" fontsize=10];
  center=1;
  "X" -> "Y" [label="-0.06"];
  "Z" -> "Y" [label="0.628"];
  "X" -> "Z" [label="0.495" dir=both];
  "X" -> "X" [label="1" dir=both];
  "Y" -> "Y" [label="0.64" dir=both];
  "Z" -> "Z" [label="1" dir=both];
}

3.2 Gvedit.exeでパスダイアグラムを出力する

 先ほど3.1で出力されたdotファイルのソースコードを適当なテキストエディタにコピペして ファイル名.dot という拡張子で保存する。それができたら、[スタート] -> [すべてのプログラム] -> [Graphviz 2.26] -> [Gvedit.exe]の順に選択していく。Gvedit.exeが起動したら、[File] -> [Open]でdotファイルを読み込む。するとコードが表示されるので[F5]キーを押すか、[Graphviz] -> [Run]の順に選択していく。あるいは画面右上の人が走るような絵が描いてあるボタンをクリックしてもよい。

 すると以下のようなウィンドウが現れる。「Output File Name」と「Output File Type」を適当に変更しても良い。今回はCドライブにpathdiagram.pngというファイルを出力せよ、としている。デフォルトでは出力先は読み込んだdotファイルが存在しているディレクトリで、ファイルタイプはgifになっている。