1 はじめに

本資料は、lavaan 公式サイトhttps://lavaan.ugent.be/ )の内容を主に参照しながら、 lavaanを使った分析の操作について簡潔に解説したものである。lavaanは現在も開発中であり、最新の情報については上記URLを参照されたい。また、拙著「縦断データ分析:構造方程式モデリングによるアプローチ」(宇佐美, 2026)の補足資料としても使用しており、本資料内における「本書」はこれを指す。 lavaanについての引用は、以下の通り:

  • Rosseel, Y. (2012). lavaan: An R Package for Structural Equation Modeling. Journal of Statistical Software, 48(2), 1–36.

2 準備

2.1 RStudioとR Markdownのインストール

lavaan は R パッケージとして実装されている。したがって、lavaan をインストールする前に、まず新しい版(4.0.0 以上)の R をインストールしておく必要がある。RStudio および R は、ここ からダウンロードできる。

R Markdown は、本資料のように、本文R コードを一つのファイルにまとめ、HTML / PDF / Word に出力できる形式である。このR Markdownファイル(lavaan_tutorial_ja.Rmd)を Knit する(=文章と Rコードをまとめて実行・変換する)には、RStudio に加えて、rmarkdownknitr を利用する必要がある。入っていない場合は、R コンソールで一度だけ次を実行する。

options(repos = c(CRAN = "https://cloud.r-project.org"))
install.packages(c("rmarkdown", "knitr"))

2.2 この資料の実行方法

このファイルを例えば lavaan_tutorial_ja.Rmd として保存したあと、RStudio 上で Knit を押すと HTML / PDF / Word に出力できる。

あるいは、R コンソールから次のように実行する。

rmarkdown::render("lavaan_tutorial_ja.Rmd")

HTML に限定すれば次のようにも書ける。

rmarkdown::render("lavaan_tutorial_ja.Rmd", output_format = "html_document")

2.3 lavaan の準備

lavaan パッケージがRにインストールされていない場合は、以下からインストールする。

if (!requireNamespace("lavaan", quietly = TRUE)) {
  install.packages("lavaan", dependencies = TRUE)
}
library(lavaan) #lavaanパッケージのダウンロード(Rを開いたら逐一必要)

2.4 組み込みデータの確認

lavaan には、以下のような学習用の組み込みデータがいくつか含まれている:

data("HolzingerSwineford1939")
data("PoliticalDemocracy")
head(HolzingerSwineford1939)
  id sex ageyr agemo  school grade       x1   x2    x3       x4   x5        x6
1  1   1    13     1 Pasteur     7 3.333333 7.75 0.375 2.333333 5.75 1.2857143
2  2   2    13     7 Pasteur     7 5.333333 5.25 2.125 1.666667 3.00 1.2857143
3  3   2    13     1 Pasteur     7 4.500000 5.25 1.875 1.000000 1.75 0.4285714
4  4   1    13     2 Pasteur     7 5.333333 7.75 3.000 2.666667 4.50 2.4285714
5  5   2    12     2 Pasteur     7 4.833333 4.75 0.875 2.666667 4.00 2.5714286
6  6   2    14     1 Pasteur     7 5.333333 5.00 2.250 1.000000 3.00 0.8571429
        x7   x8       x9
1 3.391304 5.75 6.361111
2 3.782609 6.25 7.916667
3 3.260870 3.90 4.416667
4 3.000000 5.30 4.861111
5 3.695652 6.30 5.916667
6 4.347826 6.65 7.500000
head(PoliticalDemocracy)
     y1       y2       y3       y4       y5       y6       y7       y8       x1
1  2.50 0.000000 3.333333 0.000000 1.250000 0.000000 3.726360 3.333333 4.442651
2  1.25 0.000000 3.333333 0.000000 6.250000 1.100000 6.666666 0.736999 5.384495
3  7.50 8.800000 9.999998 9.199991 8.750000 8.094061 9.999998 8.211809 5.961005
4  8.90 8.800000 9.999998 9.199991 8.907948 8.127979 9.999998 4.615086 6.285998
5 10.00 3.333333 9.999998 6.666666 7.500000 3.333333 9.999998 6.666666 5.863631
6  7.50 3.333333 6.666666 6.666666 6.250000 1.100000 6.666666 0.368500 5.533389
        x2       x3
1 3.637586 2.557615
2 5.062595 3.568079
3 6.255750 5.224433
4 7.567863 6.267495
5 6.818924 4.573679
6 5.135798 3.892270

3 lavaan を使う前に

lavaan は、潜在変数モデルのための free open-source の R パッケージである。パス解析、確認的因子分析(CFA)、構造方程式モデリング(SEM)、潜在成長モデル(LGM)、交差遅延モデル群、媒介モデルなど、幅広いモデルを設定・推定できる。

一方、公式の tutorial では現状、以下のような制約が示されている:

  • マルチレベル SEM におけるランダム傾き(random slopes)は未対応(lavaanでは、共分散構造に基づくモデルの推定が標準である[本書の2章])。
  • 連続変数・2値変数・順序カテゴリカル変数[本書の9章]以外の変数型への対応は限定的。
  • 離散型潜在変数(混合モデル、潜在クラス[本書の10章])についてのサポートはない。

4 lavaan のモデル構文 1

lavaanでは、モデルを文字列として記述する。

4.1 回帰式の設定

y ~ x1 + x2 + x3 + x4

これは、y を従属変数、x1x4 を独立変数とする回帰式である。 すなわち、記号チルダ(“~’)が回帰のオペレータであり、このオペレータの左には従属変数(y)を与え,右には独立変数(x)を“+”オペレータで区切りながら与える。以下のように、潜在変数(f1,f2,f3)を独立変数として含む場合も基本的に同様である。

y  ~ f1 + f2 + x1 + x2
f1 ~ f2 + f3
f2 ~ f3 + x1 + x2

4.2 潜在変数の定義

潜在変数はオペレータ =~ を使って定義する。以下は、例えばy1からy3は潜在変数f1を反映する変数(指標)であることを意味する。

f1 =~ y1 + y2 + y3
f2 =~ y4 + y5 + y6
f3 =~ y7 + y8 + y9 + y10

4.3 分散・共分散

分散・共分散は2重のチルダのオペレータ ~~ で表す。同じ変数同士については分散を表す。

y1 ~~ y1# 分散
y1 ~~ y2# 共分散
f1 ~~ f2# 共分散

4.4 切片・平均

切片や平均は ~ 1 で表す。

y1 ~ 1

5 確認的因子分析(CFA)の例

HolzingerSwineford1939 データ(2つの学校[PasteurとGrant-White]の7年生と8年生の生徒の知能検査のスコア[x1-x9])に対する3因子の CFA を行う。

  • visual(視覚):x1, x2, x3
  • textual(言語):x4, x5, x6
  • speed(速度):x7, x8, x9

5.1 3因子の CFAモデルの設定

モデル構文を一重引用符(’)を用いて入力する。

HS.model <- '
  visual  =~ x1 + x2 + x3
  textual =~ x4 + x5 + x6
  speed   =~ x7 + x8 + x9
'

5.2 モデルの当てはめ

次のようにして,CFAモデル(HS.model)をデータ(HolzingerSwineford1939)に当てはめる。最初の引数はユーザーが設定したモデルであり、2番目の引数は分析に用いるデータセットである。モデルを当てはめた後、関数summary()を用いて、分析結果の要約情報(推定値、標準誤差、z 値、p 値、適合度指標;本書の2.4-2.6節)を表示できる。

fit.cfa <- cfa(HS.model, data = HolzingerSwineford1939)
summary(fit.cfa, fit.measures = TRUE)
lavaan 0.6-21 ended normally after 35 iterations

  Estimator                                         ML
  Optimization method                           NLMINB
  Number of model parameters                        21

  Number of observations                           301

Model Test User Model:
                                                      
  Test statistic                                85.306
  Degrees of freedom                                24
  P-value (Chi-square)                           0.000

Model Test Baseline Model:

  Test statistic                               918.852
  Degrees of freedom                                36
  P-value                                        0.000

User Model versus Baseline Model:

  Comparative Fit Index (CFI)                    0.931
  Tucker-Lewis Index (TLI)                       0.896

Loglikelihood and Information Criteria:

  Loglikelihood user model (H0)              -3737.745
  Loglikelihood unrestricted model (H1)      -3695.092
                                                      
  Akaike (AIC)                                7517.490
  Bayesian (BIC)                              7595.339
  Sample-size adjusted Bayesian (SABIC)       7528.739

Root Mean Square Error of Approximation:

  RMSEA                                          0.092
  90 Percent confidence interval - lower         0.071
  90 Percent confidence interval - upper         0.114
  P-value H_0: RMSEA <= 0.050                    0.001
  P-value H_0: RMSEA >= 0.080                    0.840

Standardized Root Mean Square Residual:

  SRMR                                           0.065

Parameter Estimates:

  Standard errors                             Standard
  Information                                 Expected
  Information saturated (h1) model          Structured

Latent Variables:
                   Estimate  Std.Err  z-value  P(>|z|)
  visual =~                                           
    x1                1.000                           
    x2                0.554    0.100    5.554    0.000
    x3                0.729    0.109    6.685    0.000
  textual =~                                          
    x4                1.000                           
    x5                1.113    0.065   17.014    0.000
    x6                0.926    0.055   16.703    0.000
  speed =~                                            
    x7                1.000                           
    x8                1.180    0.165    7.152    0.000
    x9                1.082    0.151    7.155    0.000

Covariances:
                   Estimate  Std.Err  z-value  P(>|z|)
  visual ~~                                           
    textual           0.408    0.074    5.552    0.000
    speed             0.262    0.056    4.660    0.000
  textual ~~                                          
    speed             0.173    0.049    3.518    0.000

Variances:
                   Estimate  Std.Err  z-value  P(>|z|)
   .x1                0.549    0.114    4.833    0.000
   .x2                1.134    0.102   11.146    0.000
   .x3                0.844    0.091    9.317    0.000
   .x4                0.371    0.048    7.779    0.000
   .x5                0.446    0.058    7.642    0.000
   .x6                0.356    0.043    8.277    0.000
   .x7                0.799    0.081    9.823    0.000
   .x8                0.488    0.074    6.573    0.000
   .x9                0.566    0.071    8.003    0.000
    visual            0.809    0.145    5.564    0.000
    textual           0.979    0.112    8.737    0.000
    speed             0.384    0.086    4.451    0.000

lavaan の出力は、大きく 3つの部分に分けて理解できる。第1に、ヘッダー部分である。 ここには、lavaan のバージョン、最適化が正常終了したかどうか、反復回数、使用した推定法(例:ML)、最適化アルゴリズム(例:NLMINB)、分析モデル内の母数の数、分析に実際に用いられたサンプルサイズなどが表示される。さらに、Model Test User Model: の箇所では、分析モデルに対する検定統計量、自由度、p値 が示される。

第2に、適合度指標の部分である。これは fit.measures = TRUE を付けたときに表示され、ベースラインモデルの検定結果、CFI、TLI、対数尤度、AIC、BIC、SABIC、RMSEA(とその信頼区間)、SRMRなどが出力される。したがって、この部分を見ることで、分析モデルがデータにどの程度適合しているかを総合的に確認できる。

第3に、母数推定値の部分である。ここには Latent Variables:, Covariances:, Variances: などの見出しの下に、分析モデルに含まれる自由母数と固定母数が表形式で示される。通常は、まず因子負荷、次に共分散、最後に分散や残差分散が表示される。この表のうち、特に重要なのは次の4つの列である。

Estimate 各母数の推定値、または固定した値を表す列である。たとえば因子負荷、因子間共分散、分散などの大きさがここに示される。

Std.Err 各母数の推定値の標準誤差を表す列である。推定値の不確実性の大きさを示している。

z-value Wald 統計量を表す列であり、推定値を標準誤差で割った値である。母数が 0 からどれだけ離れているかを標準誤差単位で見た指標といえる。

P(>|z|) 母集団においてその母数が 0 であるという帰無仮説を検定したときの p 値を表す列である。

加えて、分散の出力では若干の注意が必要である。Variances: の欄で、観測変数名の前に .(ドット) が付いている場合、それはその変数が内生変数であり、表示されているのは 残差分散である。一方、外生的な潜在変数にはドットが付かず、その場合は潜在変数の全分散が表示される。

6 共通因子間の回帰モデルの例

PoliticalDemocracy データを用いた、共通因子間の回帰モデル(本書2.2.5節も参照)を示す。このデータは、Bollen(1989)以来広く利用されてきたデータセットでもある。

  • y1 1960年の報道の自由に関する専門家の評価
  • y2 1960年の政治的反対の自由
  • y3 1960年の選挙の公平性
  • y4 1960年に選出された立法府の有効性
  • y5 1965年の報道の自由に関する専門家による評価
  • y6 1965年の政治的反対の自由
  • y7 1965年の選挙の公平性
  • y8 1965年に選出された立法府の有効性
  • x1 1960年の一人当たりの国民総生産(GNP)
  • x2 1960年の一人当たりの無生物のエネルギー消費量
  • x3 1960年の産業における労働力の割合

ここで示している変数の意味(訳出)については、こちらのサイトを参考にさせて頂いた。

model.sem <- '
  ind60 =~ x1 + x2 + x3
  dem60 =~ y1 + y2 + y3 + y4
  dem65 =~ y5 + y6 + y7 + y8

  dem60 ~ ind60
  dem65 ~ ind60 + dem60

  y1 ~~ y5 #独自因子間の相関
  y2 ~~ y4 + y6 #独自因子間の相関
  y3 ~~ y7 #独自因子間の相関
  y4 ~~ y8 #独自因子間の相関
  y6 ~~ y8 #独自因子間の相関
'

fit.sem <- sem(model.sem, data = PoliticalDemocracy)
summary(fit.sem, standardized = TRUE, fit.measures = TRUE)
lavaan 0.6-21 ended normally after 68 iterations

  Estimator                                         ML
  Optimization method                           NLMINB
  Number of model parameters                        31

  Number of observations                            75

Model Test User Model:
                                                      
  Test statistic                                38.125
  Degrees of freedom                                35
  P-value (Chi-square)                           0.329

Model Test Baseline Model:

  Test statistic                               730.654
  Degrees of freedom                                55
  P-value                                        0.000

User Model versus Baseline Model:

  Comparative Fit Index (CFI)                    0.995
  Tucker-Lewis Index (TLI)                       0.993

Loglikelihood and Information Criteria:

  Loglikelihood user model (H0)              -1547.791
  Loglikelihood unrestricted model (H1)      -1528.728
                                                      
  Akaike (AIC)                                3157.582
  Bayesian (BIC)                              3229.424
  Sample-size adjusted Bayesian (SABIC)       3131.720

Root Mean Square Error of Approximation:

  RMSEA                                          0.035
  90 Percent confidence interval - lower         0.000
  90 Percent confidence interval - upper         0.092
  P-value H_0: RMSEA <= 0.050                    0.611
  P-value H_0: RMSEA >= 0.080                    0.114

Standardized Root Mean Square Residual:

  SRMR                                           0.044

Parameter Estimates:

  Standard errors                             Standard
  Information                                 Expected
  Information saturated (h1) model          Structured

Latent Variables:
                   Estimate  Std.Err  z-value  P(>|z|)   Std.lv  Std.all
  ind60 =~                                                              
    x1                1.000                               0.670    0.920
    x2                2.180    0.139   15.742    0.000    1.460    0.973
    x3                1.819    0.152   11.967    0.000    1.218    0.872
  dem60 =~                                                              
    y1                1.000                               2.223    0.850
    y2                1.257    0.182    6.889    0.000    2.794    0.717
    y3                1.058    0.151    6.987    0.000    2.351    0.722
    y4                1.265    0.145    8.722    0.000    2.812    0.846
  dem65 =~                                                              
    y5                1.000                               2.103    0.808
    y6                1.186    0.169    7.024    0.000    2.493    0.746
    y7                1.280    0.160    8.002    0.000    2.691    0.824
    y8                1.266    0.158    8.007    0.000    2.662    0.828

Regressions:
                   Estimate  Std.Err  z-value  P(>|z|)   Std.lv  Std.all
  dem60 ~                                                               
    ind60             1.483    0.399    3.715    0.000    0.447    0.447
  dem65 ~                                                               
    ind60             0.572    0.221    2.586    0.010    0.182    0.182
    dem60             0.837    0.098    8.514    0.000    0.885    0.885

Covariances:
                   Estimate  Std.Err  z-value  P(>|z|)   Std.lv  Std.all
 .y1 ~~                                                                 
   .y5                0.624    0.358    1.741    0.082    0.624    0.296
 .y2 ~~                                                                 
   .y4                1.313    0.702    1.871    0.061    1.313    0.273
   .y6                2.153    0.734    2.934    0.003    2.153    0.356
 .y3 ~~                                                                 
   .y7                0.795    0.608    1.308    0.191    0.795    0.191
 .y4 ~~                                                                 
   .y8                0.348    0.442    0.787    0.431    0.348    0.109
 .y6 ~~                                                                 
   .y8                1.356    0.568    2.386    0.017    1.356    0.338

Variances:
                   Estimate  Std.Err  z-value  P(>|z|)   Std.lv  Std.all
   .x1                0.082    0.019    4.184    0.000    0.082    0.154
   .x2                0.120    0.070    1.718    0.086    0.120    0.053
   .x3                0.467    0.090    5.177    0.000    0.467    0.239
   .y1                1.891    0.444    4.256    0.000    1.891    0.277
   .y2                7.373    1.374    5.366    0.000    7.373    0.486
   .y3                5.067    0.952    5.324    0.000    5.067    0.478
   .y4                3.148    0.739    4.261    0.000    3.148    0.285
   .y5                2.351    0.480    4.895    0.000    2.351    0.347
   .y6                4.954    0.914    5.419    0.000    4.954    0.443
   .y7                3.431    0.713    4.814    0.000    3.431    0.322
   .y8                3.254    0.695    4.685    0.000    3.254    0.315
    ind60             0.448    0.087    5.173    0.000    1.000    1.000
   .dem60             3.956    0.921    4.295    0.000    0.800    0.800
   .dem65             0.172    0.215    0.803    0.422    0.039    0.039

この例では、

  • =~:因子負荷
  • ~:回帰係数
  • ~~:分散・共分散

が全て含まれている。y1~~y5は、2つの観測変数の残差が相関することを意味する。これは,各変数が対応する潜在変数だけからでは説明できない何か別の共通原因があると考えられるときにしばしば仮定される。今の場合、2つの変数は同じ種類の変数であり異なる時期(1960年と1965年)に測定されていることを反映している。y2~~y4+y6は、2つの表現y2~~y4y2~~y6を結合したものであり、このような略記が可能である。

そして、引数standardized=TRUEは、標準化解についての出力を追加する.そのため、対応する2つの列が結果に追加されている(5,6列目)。第5列(Std.lv)は潜在変数だけが標準化され、第6列(Std.all)では潜在変数と観測変数がともに標準化された場合の結果である。

7 lavaan のモデル構文 2

母数の固定、初期値、ラベル、制約の付与について説明する。

7.1 母数の固定

因子負荷部分に、*とともに値を明示的に与えると、対応する負荷(係数・重み)がその値に固定される。

f =~ y1 + 1*y2 + 1*y3 + 1*y4

7.2 直交因子・分散固定・NA*

HS.model2 <- '
  visual  =~ x1 + x2 + x3 #識別の観点から、最初の変数(x1)の負荷が1に固定される。
  textual =~ x4 + x5 + x6 #識別の観点から、最初の変数(x4)の負荷が1に固定される。
  speed   =~ NA*x7 + x8 + x9 #最初の変数(x7)の負荷が1に固定されず自由推定される。

  visual  ~~ 0*speed #無相関を仮定
  textual ~~ 0*speed #無相関を仮定
  speed   ~~ 1*speed #(識別のため、)speedの共通因子の分散が1であることを仮定。
'

fit.cfa2 <- cfa(HS.model2, data = HolzingerSwineford1939)
summary(fit.cfa2)

7.3 orthogonal = TRUE (直交解)と std.lv = TRUE(標準化解)

fit.ortho <- cfa(HS.model, data = HolzingerSwineford1939, orthogonal = TRUE)
fit.stdlv <- cfa(HS.model, data = HolzingerSwineford1939, std.lv = TRUE)

7.4 初期値

start()を用いて、()内の数値を初期値として設定できる。

model.start <- '
  visual  =~ x1 + start(0.8)*x2 + start(1.2)*x3
  textual =~ x4 + start(0.5)*x5 + start(1.0)*x6
  speed   =~ x7 + start(0.7)*x8 + start(1.8)*x9
'

7.5 母数のラベリング

母数について名称を付与することで出力時に参照しやすくなる。また、後述する等値制約の際にもラベリングが利用できる。以下は、x3に対応する因子負荷の名称をmyLabelとしている。このように、ラベルは記号a-z, A-Zで始める必要がある。

model.label <- '
  ind60 =~ x1 + x2 + myLabel*x3
  dem60 =~ y1 + y2 + y3 + y4
  dem65 =~ y5 + y6 + y7 + y8
  dem60 ~ ind60
  dem65 ~ ind60 + dem60
  y1 ~~ y5
  y2 ~~ y4 + y6
  y3 ~~ y7
  y4 ~~ y8
  y6 ~~ y8
'

7.6 等値制約

以下のように同じラベルv2を用いることで、x2とx3に対する因子負荷を等値に制約できる。

model.eq <- '
  visual  =~ x1 + v2*x2 + v2*x3
  textual =~ x4 + x5 + x6
  speed   =~ x7 + x8 + x9
'

7.7 母数の非線形制約や不等式制約

以下では、因子負荷b1,b2,b3の間に非線形の制約(b2とb3の和の二乗がb1に等しい)や不等式制約(b1 > exp(b2 + b3))を設定している。

model.constr <- '
  y ~ b1*x1 + b2*x2 + b3*x3
  b1 == (b2 + b3)^2
  b1 > exp(b2 + b3)
'

8 潜在成長モデル(LGM)

潜在成長モデル(latent growth model; 本書の4章-6章)は縦断データの分析に用いられ、測定値の時間的変化(軌跡)やその個人差を一次や二次などの関数によって記述することを目的とする。成長因子(growth factor)と呼ばれる共通因子を介して軌跡の個人差を記述できる。下記の例では,Demo.growthという人工データを利用している。これは、読解力尺度に関する標準化された得点のデータで、4時点で測定されている。このデータに対して線形LGM(本書の4章;切片因子iと傾き因子sからなる2つの成長因子を含むLGM)を当てはめる:

8.1 growth() を用いた基本例

model.growth <- '
  i =~ 1*t1 + 1*t2 + 1*t3 + 1*t4
  s =~ 0*t1 + 1*t2 + 2*t3 + 3*t4
'

fit.growth <- growth(model.growth, data = Demo.growth)
summary(fit.growth)
lavaan 0.6-21 ended normally after 29 iterations

  Estimator                                         ML
  Optimization method                           NLMINB
  Number of model parameters                         9

  Number of observations                           400

Model Test User Model:
                                                      
  Test statistic                                 8.069
  Degrees of freedom                                 5
  P-value (Chi-square)                           0.152

Parameter Estimates:

  Standard errors                             Standard
  Information                                 Expected
  Information saturated (h1) model          Structured

Latent Variables:
                   Estimate  Std.Err  z-value  P(>|z|)
  i =~                                                
    t1                1.000                           
    t2                1.000                           
    t3                1.000                           
    t4                1.000                           
  s =~                                                
    t1                0.000                           
    t2                1.000                           
    t3                2.000                           
    t4                3.000                           

Covariances:
                   Estimate  Std.Err  z-value  P(>|z|)
  i ~~                                                
    s                 0.618    0.071    8.686    0.000

Intercepts:
                   Estimate  Std.Err  z-value  P(>|z|)
    i                 0.615    0.077    8.007    0.000
    s                 1.006    0.042   24.076    0.000

Variances:
                   Estimate  Std.Err  z-value  P(>|z|)
   .t1                0.595    0.086    6.944    0.000
   .t2                0.676    0.061   11.061    0.000
   .t3                0.635    0.072    8.761    0.000
   .t4                0.508    0.124    4.090    0.000
    i                 1.932    0.173   11.194    0.000
    s                 0.587    0.052   11.336    0.000

8.2 cfa() を用いた明示的な指定

cfa()を用いた設定も可能である。一見冗長にはなるが、本書の中で提供されているコードは、表現の統一およびモデル上の仮定の明示の観点から、基本的にLGMにおいてもcfa()を利用している。

model.growth2 <- '
  i =~ 1*t1 + 1*t2 + 1*t3 + 1*t4
  s =~ 0*t1 + 1*t2 + 2*t3 + 3*t4

  t1 ~ 0*1
  t2 ~ 0*1
  t3 ~ 0*1
  t4 ~ 0*1

  i ~ mui*1
  s ~ mus*1

  i ~~ phii*i
  s ~~ phis*s
  i ~~ phiis*s

  t1 ~~ psi1*t1
  t2 ~~ psi2*t2
  t3 ~~ psi3*t3
  t4 ~~ psi4*t4
'

fit.growth2 <- cfa(model.growth2, data = Demo.growth)
summary(fit.growth2)
lavaan 0.6-21 ended normally after 29 iterations

  Estimator                                         ML
  Optimization method                           NLMINB
  Number of model parameters                         9

  Number of observations                           400

Model Test User Model:
                                                      
  Test statistic                                 8.069
  Degrees of freedom                                 5
  P-value (Chi-square)                           0.152

Parameter Estimates:

  Standard errors                             Standard
  Information                                 Expected
  Information saturated (h1) model          Structured

Latent Variables:
                   Estimate  Std.Err  z-value  P(>|z|)
  i =~                                                
    t1                1.000                           
    t2                1.000                           
    t3                1.000                           
    t4                1.000                           
  s =~                                                
    t1                0.000                           
    t2                1.000                           
    t3                2.000                           
    t4                3.000                           

Covariances:
                   Estimate  Std.Err  z-value  P(>|z|)
  i ~~                                                
    s       (phis)    0.618    0.071    8.686    0.000

Intercepts:
                   Estimate  Std.Err  z-value  P(>|z|)
   .t1                0.000                           
   .t2                0.000                           
   .t3                0.000                           
   .t4                0.000                           
    i        (mui)    0.615    0.077    8.007    0.000
    s        (mus)    1.006    0.042   24.076    0.000

Variances:
                   Estimate  Std.Err  z-value  P(>|z|)
    i       (phii)    1.932    0.173   11.194    0.000
    s       (phis)    0.587    0.052   11.336    0.000
   .t1      (psi1)    0.595    0.086    6.944    0.000
   .t2      (psi2)    0.676    0.061   11.061    0.000
   .t3      (psi3)    0.635    0.072    8.761    0.000
   .t4      (psi4)    0.508    0.124    4.090    0.000

8.3 条件付潜在成長モデル(本書の4.6節)

時不変的な説明変数(x1,x2)の他、時変的な説明変数(c1,c2,c3,c4)を含む場合の条件付LGM。

model.cg <- '
  i =~ 1*t1 + 1*t2 + 1*t3 + 1*t4
  s =~ 0*t1 + 1*t2 + 2*t3 + 3*t4

  i ~ x1 + x2
  s ~ x1 + x2

  t1 ~ c1
  t2 ~ c2
  t3 ~ c3
  t4 ~ c4

  t1 ~ 0*1
  t2 ~ 0*1
  t3 ~ 0*1
  t4 ~ 0*1

  i ~ inti*1
  s ~ ints*1

  i ~~ phii*i
  s ~~ phis*s
  i ~~ phiis*s

  t1 ~~ psi1*t1
  t2 ~~ psi2*t2
  t3 ~~ psi3*t3
  t4 ~~ psi4*t4
'
fit.conditionalgrowth <- cfa(model.cg, data = Demo.growth)
summary(fit.conditionalgrowth)
lavaan 0.6-21 ended normally after 31 iterations

  Estimator                                         ML
  Optimization method                           NLMINB
  Number of model parameters                        17

  Number of observations                           400

Model Test User Model:
                                                      
  Test statistic                                26.059
  Degrees of freedom                                21
  P-value (Chi-square)                           0.204

Parameter Estimates:

  Standard errors                             Standard
  Information                                 Expected
  Information saturated (h1) model          Structured

Latent Variables:
                   Estimate  Std.Err  z-value  P(>|z|)
  i =~                                                
    t1                1.000                           
    t2                1.000                           
    t3                1.000                           
    t4                1.000                           
  s =~                                                
    t1                0.000                           
    t2                1.000                           
    t3                2.000                           
    t4                3.000                           

Regressions:
                   Estimate  Std.Err  z-value  P(>|z|)
  i ~                                                 
    x1                0.608    0.060   10.134    0.000
    x2                0.604    0.064    9.412    0.000
  s ~                                                 
    x1                0.262    0.029    9.198    0.000
    x2                0.522    0.031   17.083    0.000
  t1 ~                                                
    c1                0.143    0.050    2.883    0.004
  t2 ~                                                
    c2                0.289    0.046    6.295    0.000
  t3 ~                                                
    c3                0.328    0.044    7.361    0.000
  t4 ~                                                
    c4                0.330    0.058    5.655    0.000

Covariances:
                   Estimate  Std.Err  z-value  P(>|z|)
 .i ~~                                                
   .s       (phis)    0.075    0.040    1.855    0.064

Intercepts:
                   Estimate  Std.Err  z-value  P(>|z|)
   .t1                0.000                           
   .t2                0.000                           
   .t3                0.000                           
   .t4                0.000                           
   .i       (inti)    0.580    0.062    9.368    0.000
   .s       (ints)    0.958    0.029   32.552    0.000

Variances:
                   Estimate  Std.Err  z-value  P(>|z|)
   .i       (phii)    1.079    0.112    9.609    0.000
   .s       (phis)    0.224    0.027    8.429    0.000
   .t1      (psi1)    0.580    0.080    7.230    0.000
   .t2      (psi2)    0.596    0.054   10.969    0.000
   .t3      (psi3)    0.481    0.055    8.745    0.000
   .t4      (psi4)    0.535    0.098    5.466    0.000

9 交差遅延パネルモデル(CLPM;本書の7章)

別添の(欠測値NAを含む、人工データである) semexample.csvデータを利用した、交差遅延パネルモデルによる分析例を示す。利用者は、まずこの semexample.csvデータを ダウンロードし、lavaan_tutorial_ja.Rmd と同じフォルダに保存しておくことを前提とする。

SLEEP <- read.csv("C:/Users/Satoshi Usami/Desktop/semexample.csv", header = TRUE) #分析者のPCの環境に合わせて変更する必要あり。データの読み込みの方法は、本資料の最後(21節)を参照のこと。
head(SLEEP)
  ID Sex GHQ1 GHQ2 GHQ3 GHQ4 GHQ5 GHQ6 sleep1 sleep2 sleep3 sleep4 sleep5
1  1  NA    0   NA   NA   NA   NA   NA     NA     NA     NA     NA     NA
2  2   1    3   NA   NA   NA   NA   NA    9.0     NA     NA     NA     NA
3  3   1    0   NA   NA   NA   NA   NA    8.5     NA     NA     NA     NA
4  4  NA    0   NA   NA   NA   NA   NA    9.0     NA     NA     NA     NA
5  5  NA    0   NA   NA   NA   NA   NA    8.0     NA     NA     NA     NA
6  6  NA    0   NA   NA   NA   NA   NA    8.5     NA     NA     NA     NA
  sleep6 bed1 bed2 bed3 bed4 bed5 bed6
1     NA   NA   NA   NA   NA   NA   NA
2     NA 21.0   NA   NA   NA   NA   NA
3     NA 21.0   NA   NA   NA   NA   NA
4     NA 21.0   NA   NA   NA   NA   NA
5     NA 21.5   NA   NA   NA   NA   NA
6     NA 21.5   NA   NA   NA   NA   NA
CLPMT4 <- '
  sleep1 ~ muy1*1
  sleep2 ~ inty2*1
  sleep3 ~ inty3*1
  sleep4 ~ inty4*1

  GHQ1 ~ mux1*1
  GHQ2 ~ intx2*1
  GHQ3 ~ intx3*1
  GHQ4 ~ intx4*1

  sleep2 ~ betay2*sleep1 + gammay2*GHQ1
  GHQ2   ~ betax2*GHQ1 + gammax2*sleep1

  sleep3 ~ betay3*sleep2 + gammay3*GHQ2
  GHQ3   ~ betax3*GHQ2 + gammax3*sleep2

  sleep4 ~ betay4*sleep3 + gammay4*GHQ3
  GHQ4   ~ betax4*GHQ3 + gammax4*sleep3

  GHQ2   ~~ Omegax2*GHQ2
  sleep2 ~~ Omegay2*sleep2
  GHQ2   ~~ Omegaxy2*sleep2

  GHQ3   ~~ Omegax3*GHQ3
  sleep3 ~~ Omegay3*sleep3
  GHQ3   ~~ Omegaxy3*sleep3

  GHQ4   ~~ Omegax4*GHQ4
  sleep4 ~~ Omegay4*sleep4
  GHQ4   ~~ Omegaxy4*sleep4
'

fit.clpm <- sem(CLPMT4, data = SLEEP, missing = "fiml") #fimlは完全情報最尤推定法による欠測処理
summary(fit.clpm, fit.measures = TRUE)
lavaan 0.6-21 ended normally after 76 iterations

  Estimator                                         ML
  Optimization method                           NLMINB
  Number of model parameters                        31

                                                  Used       Total
  Number of observations                           826        1067
  Number of missing patterns                        30            

Model Test User Model:
                                                      
  Test statistic                               111.985
  Degrees of freedom                                13
  P-value (Chi-square)                           0.000

Model Test Baseline Model:

  Test statistic                               850.195
  Degrees of freedom                                28
  P-value                                        0.000

User Model versus Baseline Model:

  Comparative Fit Index (CFI)                    0.880
  Tucker-Lewis Index (TLI)                       0.741
                                                      
  Robust Comparative Fit Index (CFI)             0.829
  Robust Tucker-Lewis Index (TLI)                0.631

Loglikelihood and Information Criteria:

  Loglikelihood user model (H0)              -7388.811
  Loglikelihood unrestricted model (H1)      -7332.819
                                                      
  Akaike (AIC)                               14839.622
  Bayesian (BIC)                             14985.837
  Sample-size adjusted Bayesian (SABIC)      14887.392

Root Mean Square Error of Approximation:

  RMSEA                                          0.096
  90 Percent confidence interval - lower         0.080
  90 Percent confidence interval - upper         0.113
  P-value H_0: RMSEA <= 0.050                    0.000
  P-value H_0: RMSEA >= 0.080                    0.951
                                                      
  Robust RMSEA                                   0.174
  90 Percent confidence interval - lower         0.143
  90 Percent confidence interval - upper         0.206
  P-value H_0: Robust RMSEA <= 0.050             0.000
  P-value H_0: Robust RMSEA >= 0.080             1.000

Standardized Root Mean Square Residual:

  SRMR                                           0.095

Parameter Estimates:

  Standard errors                             Standard
  Information                                 Observed
  Observed information based on                Hessian

Regressions:
                   Estimate  Std.Err  z-value  P(>|z|)
  sleep2 ~                                            
    sleep1  (bty2)    0.573    0.045   12.742    0.000
    GHQ1   (gmmy2)   -0.036    0.020   -1.811    0.070
  GHQ2 ~                                              
    GHQ1    (btx2)    0.569    0.045   12.576    0.000
    sleep1 (gmmx2)   -0.032    0.103   -0.310    0.757
  sleep3 ~                                            
    sleep2  (bty3)    0.440    0.041   10.637    0.000
    GHQ2   (gmmy3)   -0.015    0.019   -0.787    0.431
  GHQ3 ~                                              
    GHQ2    (btx3)    0.724    0.046   15.738    0.000
    sleep2 (gmmx3)   -0.172    0.117   -1.469    0.142
  sleep4 ~                                            
    sleep3  (bty4)    0.393    0.043    9.145    0.000
    GHQ3   (gmmy4)   -0.025    0.015   -1.655    0.098
  GHQ4 ~                                              
    GHQ3    (btx4)    0.467    0.041   11.399    0.000
    sleep3 (gmmx4)   -0.208    0.119   -1.752    0.080

Covariances:
                   Estimate  Std.Err  z-value  P(>|z|)
 .sleep2 ~~                                           
   .GHQ2    (Omg2)   -0.150    0.102   -1.463    0.144
 .sleep3 ~~                                           
   .GHQ3    (Omg3)   -0.340    0.114   -2.981    0.003
 .sleep4 ~~                                           
   .GHQ4    (Omg4)   -0.098    0.112   -0.872    0.383

Intercepts:
                   Estimate  Std.Err  z-value  P(>|z|)
    sleep1  (muy1)    7.476    0.043  174.757    0.000
   .sleep2 (inty2)    2.911    0.348    8.364    0.000
   .sleep3 (inty3)    3.907    0.310   12.608    0.000
   .sleep4 (inty4)    4.039    0.315   12.815    0.000
    GHQ1    (mux1)    2.028    0.099   20.557    0.000
   .GHQ2   (intx2)    1.596    0.796    2.005    0.045
   .GHQ3   (intx3)    2.355    0.883    2.666    0.008
   .GHQ4   (intx4)    2.846    0.871    3.267    0.001

Variances:
                   Estimate  Std.Err  z-value  P(>|z|)
   .GHQ2   (Omgx2)    4.768    0.326   14.618    0.000
   .sleep2 (Omgy2)    0.959    0.062   15.344    0.000
   .GHQ3   (Omgx3)    5.262    0.411   12.790    0.000
   .sleep3 (Omgy3)    0.890    0.059   15.189    0.000
   .GHQ4   (Omgx4)    6.250    0.436   14.349    0.000
   .sleep4 (Omgy4)    0.822    0.057   14.405    0.000
    sleep1            1.077    0.063   17.200    0.000
    GHQ1              5.796    0.344   16.870    0.000

9.1 補足

交差遅延パネルモデルからの様々な派生モデルの種類やその概説、lavaanコードについては以下も参照のこと.

宇佐美慧 (2022). 個人内関係の推測と統計モデル―ランダム切片交差遅延パネルモデルを巡って―

Rコード

10 平均構造(Meanstructure)

SEMでは、例えば因子分析モデルを扱う場合のように、共分散(相関)構造のみに関心がもたれる場合も多いが、例えば潜在成長モデル(LGM; 本書の4章)のように平均構造に関心がもたれることも多い。より一般に、平均共分散構造分析(MACS)としてのSEMでは、変数の切片や平均も推定の対象となる。先述した3因子のCFAモデルでは,観測変数の切片を次のようにして加えることができる。

10.1 切片の設定

HS.mean <- '
  visual  =~ x1 + x2 + x3
  textual =~ x4 + x5 + x6
  speed   =~ x7 + x8 + x9

#以下、切片の設定
  x1 ~ 1
  x2 ~ 1
  x3 ~ 1
  x4 ~ 1
  x5 ~ 1
  x6 ~ 1
  x7 ~ 1
  x8 ~ 1
  x9 ~ 1
'

10.2 meanstructure = TRUE

モデル構文で切片を明示しないかわりにcfa()においてmeanstructure=TRUEを加える方法もある:

fit.mean <- cfa(HS.model,
                data = HolzingerSwineford1939,
                meanstructure = TRUE)
summary(fit.mean)
lavaan 0.6-21 ended normally after 35 iterations

  Estimator                                         ML
  Optimization method                           NLMINB
  Number of model parameters                        30

  Number of observations                           301

Model Test User Model:
                                                      
  Test statistic                                85.306
  Degrees of freedom                                24
  P-value (Chi-square)                           0.000

Parameter Estimates:

  Standard errors                             Standard
  Information                                 Expected
  Information saturated (h1) model          Structured

Latent Variables:
                   Estimate  Std.Err  z-value  P(>|z|)
  visual =~                                           
    x1                1.000                           
    x2                0.554    0.100    5.554    0.000
    x3                0.729    0.109    6.685    0.000
  textual =~                                          
    x4                1.000                           
    x5                1.113    0.065   17.014    0.000
    x6                0.926    0.055   16.703    0.000
  speed =~                                            
    x7                1.000                           
    x8                1.180    0.165    7.152    0.000
    x9                1.082    0.151    7.155    0.000

Covariances:
                   Estimate  Std.Err  z-value  P(>|z|)
  visual ~~                                           
    textual           0.408    0.074    5.552    0.000
    speed             0.262    0.056    4.660    0.000
  textual ~~                                          
    speed             0.173    0.049    3.518    0.000

Intercepts:
                   Estimate  Std.Err  z-value  P(>|z|)
   .x1                4.936    0.067   73.473    0.000
   .x2                6.088    0.068   89.855    0.000
   .x3                2.250    0.065   34.579    0.000
   .x4                3.061    0.067   45.694    0.000
   .x5                4.341    0.074   58.452    0.000
   .x6                2.186    0.063   34.667    0.000
   .x7                4.186    0.063   66.766    0.000
   .x8                5.527    0.058   94.854    0.000
   .x9                5.374    0.058   92.546    0.000

Variances:
                   Estimate  Std.Err  z-value  P(>|z|)
   .x1                0.549    0.114    4.833    0.000
   .x2                1.134    0.102   11.146    0.000
   .x3                0.844    0.091    9.317    0.000
   .x4                0.371    0.048    7.779    0.000
   .x5                0.446    0.058    7.642    0.000
   .x6                0.356    0.043    8.277    0.000
   .x7                0.799    0.081    9.823    0.000
   .x8                0.488    0.074    6.573    0.000
   .x9                0.566    0.071    8.003    0.000
    visual            0.809    0.145    5.564    0.000
    textual           0.979    0.112    8.737    0.000
    speed             0.384    0.086    4.451    0.000

11 多群分析(Multiple Groups)

データセット内の群変数(集団=群を表す変数)の名前を引数groupを通して与えることで、多群モデルの分析(本書の6章)ができる。デフォルトでは、各群に対して同じ(つまり、群間の制約がない)モデルが当てはめられる。次の例では、2つの学校(PasteurとGrant-White)のデータセットに対して先述の3因子のCFAモデルを当てはめている。

fit.group <- cfa(HS.model,
                 data = HolzingerSwineford1939,
                 group = "school")
summary(fit.group)
lavaan 0.6-21 ended normally after 57 iterations

  Estimator                                         ML
  Optimization method                           NLMINB
  Number of model parameters                        60

  Number of observations per group:                   
    Pasteur                                        156
    Grant-White                                    145

Model Test User Model:
                                                      
  Test statistic                               115.851
  Degrees of freedom                                48
  P-value (Chi-square)                           0.000
  Test statistic for each group:
    Pasteur                                     64.309
    Grant-White                                 51.542

Parameter Estimates:

  Standard errors                             Standard
  Information                                 Expected
  Information saturated (h1) model          Structured


Group 1 [Pasteur]:

Latent Variables:
                   Estimate  Std.Err  z-value  P(>|z|)
  visual =~                                           
    x1                1.000                           
    x2                0.394    0.122    3.220    0.001
    x3                0.570    0.140    4.076    0.000
  textual =~                                          
    x4                1.000                           
    x5                1.183    0.102   11.613    0.000
    x6                0.875    0.077   11.421    0.000
  speed =~                                            
    x7                1.000                           
    x8                1.125    0.277    4.057    0.000
    x9                0.922    0.225    4.104    0.000

Covariances:
                   Estimate  Std.Err  z-value  P(>|z|)
  visual ~~                                           
    textual           0.479    0.106    4.531    0.000
    speed             0.185    0.077    2.397    0.017
  textual ~~                                          
    speed             0.182    0.069    2.628    0.009

Intercepts:
                   Estimate  Std.Err  z-value  P(>|z|)
   .x1                4.941    0.095   52.249    0.000
   .x2                5.984    0.098   60.949    0.000
   .x3                2.487    0.093   26.778    0.000
   .x4                2.823    0.092   30.689    0.000
   .x5                3.995    0.105   38.183    0.000
   .x6                1.922    0.079   24.321    0.000
   .x7                4.432    0.087   51.181    0.000
   .x8                5.563    0.078   71.214    0.000
   .x9                5.418    0.079   68.440    0.000

Variances:
                   Estimate  Std.Err  z-value  P(>|z|)
   .x1                0.298    0.232    1.286    0.198
   .x2                1.334    0.158    8.464    0.000
   .x3                0.989    0.136    7.271    0.000
   .x4                0.425    0.069    6.138    0.000
   .x5                0.456    0.086    5.292    0.000
   .x6                0.290    0.050    5.780    0.000
   .x7                0.820    0.125    6.580    0.000
   .x8                0.510    0.116    4.406    0.000
   .x9                0.680    0.104    6.516    0.000
    visual            1.097    0.276    3.967    0.000
    textual           0.894    0.150    5.963    0.000
    speed             0.350    0.126    2.778    0.005


Group 2 [Grant-White]:

Latent Variables:
                   Estimate  Std.Err  z-value  P(>|z|)
  visual =~                                           
    x1                1.000                           
    x2                0.736    0.155    4.760    0.000
    x3                0.925    0.166    5.583    0.000
  textual =~                                          
    x4                1.000                           
    x5                0.990    0.087   11.418    0.000
    x6                0.963    0.085   11.377    0.000
  speed =~                                            
    x7                1.000                           
    x8                1.226    0.187    6.569    0.000
    x9                1.058    0.165    6.429    0.000

Covariances:
                   Estimate  Std.Err  z-value  P(>|z|)
  visual ~~                                           
    textual           0.408    0.098    4.153    0.000
    speed             0.276    0.076    3.639    0.000
  textual ~~                                          
    speed             0.222    0.073    3.022    0.003

Intercepts:
                   Estimate  Std.Err  z-value  P(>|z|)
   .x1                4.930    0.095   51.696    0.000
   .x2                6.200    0.092   67.416    0.000
   .x3                1.996    0.086   23.195    0.000
   .x4                3.317    0.093   35.625    0.000
   .x5                4.712    0.096   48.986    0.000
   .x6                2.469    0.094   26.277    0.000
   .x7                3.921    0.086   45.819    0.000
   .x8                5.488    0.087   63.174    0.000
   .x9                5.327    0.085   62.571    0.000

Variances:
                   Estimate  Std.Err  z-value  P(>|z|)
   .x1                0.715    0.126    5.676    0.000
   .x2                0.899    0.123    7.339    0.000
   .x3                0.557    0.103    5.409    0.000
   .x4                0.315    0.065    4.870    0.000
   .x5                0.419    0.072    5.812    0.000
   .x6                0.406    0.069    5.880    0.000
   .x7                0.600    0.091    6.584    0.000
   .x8                0.401    0.094    4.249    0.000
   .x9                0.535    0.089    6.010    0.000
    visual            0.604    0.160    3.762    0.000
    textual           0.942    0.152    6.177    0.000
    speed             0.461    0.118    3.910    0.000

11.1 複数のグループ(群;G=2)に対して母数(因子負荷)のラベリングを行う例

デフォルトでは各集団に対して個別に母数を推定しているので、上の設定は以下のようにラベリングすることと等価である:

HS.model.group <- '
  visual  =~ c(a11,a12)*x1 + c(a21,a22)*x2 + c(a31,a32)*x3
  textual =~ c(a41,a42)*x4 + c(a51,a52)*x5 + c(a61,a62)*x6
  speed   =~ c(a71,a72)*x7 + c(a81,a82)*x8 + c(a91,a92)*x9
'

11.2 グループ間で因子負荷の等値制約を課す例

もしx1,x2,x3の因子負荷に等値制約を課す場合、以下のようにする。 これは、測定不変性の観点から弱不変(本書の6.2節)のモデルを設定していることを意味する。

HS.model.invar <- '
  visual =~ c(a1,a1)*x1 + c(a2,a2)*x2 + c(a3,a3)*x3
'

12 カテゴリカルデータ

lavaan では、外生的なカテゴリカル変数と内生的なカテゴリカル変数を区別して考える必要がある。

12.1 外生的なカテゴリカル変数

  • 2値外生変数:通常のダミー変数(0/1)に変換する
  • 順序外生変数:順序を反映する数値として扱うことが多い
  • (カテゴリ数がKの)名義外生変数:K-1 個のダミー変数に置き換える

12.2 内生的なカテゴリカル変数

12.2.1 方法1:あらかじめ ordered にする

分析を行う前に、データフレーム内の対応する変数を‘順序付き’変数と宣言する方法(ordered()関数を利用)。例えば、データフレーム(’Data’とする)にある4つの変数(item1,item2,item3,item4)の尺度水準が順序尺度と宣言するとき、次のようにする:

Data[, c("item1", "item2", "item3", "item4")] <-
  lapply(Data[, c("item1", "item2", "item3", "item4")], ordered)

12.2.2 方法2:関数内で ordered = を指定する

もう一つは、以下のように、cfa関数において引数ordered=を利用する方法がある:

fit.cat <- cfa(myModel,
               data = myData,
               ordered = c("item1", "item2", "item3", "item4"))

扱う変数が全てカテゴリカル変数の場合、 ordered = TRUE という略記も可能。なお、ordered= を使うと lavaan は自動的に WLSMV 系の推定へ切り替わる(本書の9章)。

13 共分散行列だけを入力する方法

生データがなくても、標本共分散行列とサンプルサイズの情報が入手できればモデルの当てはめが可能である(例えば本書の2.4.4節)。

#共分散行列の読み込み
wheaton.cov <- lav_matrix_lower2full(
  c(11.834,
    6.947, 9.364,
    6.819, 5.091, 12.532,
    4.783, 5.028, 7.495, 9.986,
   -3.839, -3.889, -3.841, -3.625, 9.610,
   -21.899, -18.831, -21.748, -18.775, 35.522, 450.288)
)
#変数名の付与
colnames(wheaton.cov) <- rownames(wheaton.cov) <-
  c("anomia67", "powerless67", "anomia71", "powerless71", "education", "sei")

#分析モデルの設定
wheaton.model <- '
  ses     =~ education + sei
  alien67 =~ anomia67 + powerless67
  alien71 =~ anomia71 + powerless71

  alien71 ~ alien67 + ses
  alien67 ~ ses

  anomia67 ~~ anomia71
  powerless67 ~~ powerless71
'

fit.cov <- sem(wheaton.model,
               sample.cov = wheaton.cov,
               sample.nobs = 932) #サンプルサイズを明示する。
summary(fit.cov, standardized = TRUE)
lavaan 0.6-21 ended normally after 84 iterations

  Estimator                                         ML
  Optimization method                           NLMINB
  Number of model parameters                        17

  Number of observations                           932

Model Test User Model:
                                                      
  Test statistic                                 4.735
  Degrees of freedom                                 4
  P-value (Chi-square)                           0.316

Parameter Estimates:

  Standard errors                             Standard
  Information                                 Expected
  Information saturated (h1) model          Structured

Latent Variables:
                   Estimate  Std.Err  z-value  P(>|z|)   Std.lv  Std.all
  ses =~                                                                
    education         1.000                               2.607    0.842
    sei               5.219    0.422   12.364    0.000   13.609    0.642
  alien67 =~                                                            
    anomia67          1.000                               2.663    0.774
    powerless67       0.979    0.062   15.895    0.000    2.606    0.852
  alien71 =~                                                            
    anomia71          1.000                               2.850    0.805
    powerless71       0.922    0.059   15.498    0.000    2.628    0.832

Regressions:
                   Estimate  Std.Err  z-value  P(>|z|)   Std.lv  Std.all
  alien71 ~                                                             
    alien67           0.607    0.051   11.898    0.000    0.567    0.567
    ses              -0.227    0.052   -4.334    0.000   -0.207   -0.207
  alien67 ~                                                             
    ses              -0.575    0.056  -10.195    0.000   -0.563   -0.563

Covariances:
                   Estimate  Std.Err  z-value  P(>|z|)   Std.lv  Std.all
 .anomia67 ~~                                                           
   .anomia71          1.623    0.314    5.176    0.000    1.623    0.356
 .powerless67 ~~                                                        
   .powerless71       0.339    0.261    1.298    0.194    0.339    0.121

Variances:
                   Estimate  Std.Err  z-value  P(>|z|)   Std.lv  Std.all
   .education         2.801    0.507    5.525    0.000    2.801    0.292
   .sei             264.597   18.126   14.597    0.000  264.597    0.588
   .anomia67          4.731    0.453   10.441    0.000    4.731    0.400
   .powerless67       2.563    0.403    6.359    0.000    2.563    0.274
   .anomia71          4.399    0.515    8.542    0.000    4.399    0.351
   .powerless71       3.070    0.434    7.070    0.000    3.070    0.308
    ses               6.798    0.649   10.475    0.000    1.000    1.000
   .alien67           4.841    0.467   10.359    0.000    0.683    0.683
   .alien71           4.083    0.404   10.104    0.000    0.503    0.503

13.1 sample.cov.rescale

ML 推定では、入力共分散行列が通常 N-1 で割られていることを前提に、内部で (N-1)/N の rescale が行われる。必要に応じて sample.cov.rescale = FALSE を指定する。

14 推定法・欠測値・標準誤差・ブートストラップ

14.1 主な推定法

連続変数のみを扱う場合、デフォルトは ML である。その他、次の推定法が利用可能である(例えば、本書の2.4節や9.2節)。

# estimator = "ML"
# estimator = "GLS"
# estimator = "WLS"
# estimator = "DWLS"
# estimator = "ULS"
# estimator = "DLS"
# estimator = "PML"

ロバスト版としては、例えば次がある(標準誤差、カイ二乗検定統計量、およびそれに基づく適合度評価がロバストに補正され、通常、点推定値自体はMLと同じ)。

# estimator = "MLM"
# estimator = "MLMVS"
# estimator = "MLMV"
# estimator = "MLR"

14.2 補足:カイ二乗統計量(適合度検定統計量)の計算方式

lavaan で最尤推定(estimator = "ML")またはそのロバスト版を用いるとき、デフォルトでは “normal” アプローチが採用される。これは、標本共分散行列を \(N\) で割った形(いわゆる biased sample covariance matrix)に基づいて内部的に計算し、カイ二乗統計量(尤度比検定統計量:-2logλ)も最小化された目的関数値に \(N\) を掛けて求めるという方法である。

これに対して、慣習的な理由で、不偏共分散行列、すなわち各要素を \(N-1\) で割った共分散行列に基づいて計算し、さらにカイ二乗統計量の計算でも \(N-1\)を掛けることも多い(本書の2.4.1節;(N−1)fML)。その場合には、likelihood = "wishart" を指定する。たとえば、先述のCFAモデルでは次のように書く。

fit <- cfa(
  HS.model,
  data = HolzingerSwineford1939,
  likelihood = "wishart"
)
fit
lavaan 0.6-21 ended normally after 35 iterations

  Estimator                                         ML
  Optimization method                           NLMINB
  Number of model parameters                        21

  Number of observations                           301

Model Test User Model:
                                                      
  Test statistic                                85.022
  Degrees of freedom                                24
  P-value (Chi-square)                           0.000

likelihood = "wishart" を指定すると、主として検定統計量のスケーリングが変わる。そのため、同じモデル・同じデータでも、lavaan のデフォルト設定(normal)で得られるカイ二乗値と、Wishart 指定で得られるカイ二乗値はわずかに異なりうる(サンプルサイズが大きくなればこの違いは実用上無視してよい)。本書では便宜的に、3章以降、lavaan のデフォルト設定(normal)の下でのカイ二乗値を示している。

14.3 欠測値

FIML(例えば、本書の11.4節、11.6節)を使った欠測データの処理。

fit.missing <- sem(model.sem,
                   data = PoliticalDemocracy,
                   missing = "fiml")
summary(fit.missing)

14.4 標準誤差

MLR(本書の2.4.5節)によるロバストな標準誤差の出力(ロバスト最尤推定法)。

fit.robust <- sem(model.sem,
                  data = PoliticalDemocracy,
                  estimator = "MLR")
summary(fit.robust, fit.measures = TRUE)

14.5 ブートストラップ

ブートストラップを使った標準誤差の推定。モデルを当てはめるときにse=“bootstrap”または test=“bootstrap”と指定する。

fit.boot <- sem(model.sem,
                data = PoliticalDemocracy,
                se = "bootstrap",
                bootstrap = 1000)
summary(fit.boot)

15 blavaan によるベイズ推定(e.g., 本書の2.4.6節)

blavaan (Merkle & Rosseel, 2018)および、組み込みデータセットPoliticalDemocracyと先述のmodel.sem(共通因子間の回帰モデル)を用いたベイズ推定の例を示す。blavaanは、lavaan に近い構文を用いながらSEMのベイズ推定が実装できるR パッケージであり、sem() に対応する bsem()や、他にもbcfa() などを使ってモデルを推定できる。推定ではJAGSやStanを利用したMCMCが用いられる。blavaanの詳細については上記サイトを参照されたい。

install.packages("blavaan")
パッケージ 'blavaan' は無事に展開され、MD5 サムもチェックされました

ダウンロードされたパッケージは、以下にあります
    C:\Users\Satoshi Usami\AppData\Local\Temp\RtmpIxTDFF\downloaded_packages
library(blavaan)

#事前分布の設定。dpriors()では母数のクラスごとに同じ形式の分布を設定する
mydp <- dpriors(
  nu     = "normal(0,100)",
  alpha  = "normal(0,100)",
  lambda = "normal(0,10)",
  beta   = "normal(0,10)",
  theta  = "gamma(1,0.01)",
  psi    = "gamma(1,0.01)"
)

fit.bayes <- bcfa(model.sem, #共通因子間の回帰モデル(6節)
                  data = PoliticalDemocracy,
                  dp = mydp, #mydpで事前分布を設定する。
                  sample = 2000)

SAMPLING FOR MODEL 'stanmarg' NOW (CHAIN 1).
Chain 1: 
Chain 1: Gradient evaluation took 0.000957 seconds
Chain 1: 1000 transitions using 10 leapfrog steps per transition would take 9.57 seconds.
Chain 1: Adjust your expectations accordingly!
Chain 1: 
Chain 1: 
Chain 1: Iteration:    1 / 2500 [  0%]  (Warmup)
Chain 1: Iteration:  250 / 2500 [ 10%]  (Warmup)
Chain 1: Iteration:  500 / 2500 [ 20%]  (Warmup)
Chain 1: Iteration:  501 / 2500 [ 20%]  (Sampling)
Chain 1: Iteration:  750 / 2500 [ 30%]  (Sampling)
Chain 1: Iteration: 1000 / 2500 [ 40%]  (Sampling)
Chain 1: Iteration: 1250 / 2500 [ 50%]  (Sampling)
Chain 1: Iteration: 1500 / 2500 [ 60%]  (Sampling)
Chain 1: Iteration: 1750 / 2500 [ 70%]  (Sampling)
Chain 1: Iteration: 2000 / 2500 [ 80%]  (Sampling)
Chain 1: Iteration: 2250 / 2500 [ 90%]  (Sampling)
Chain 1: Iteration: 2500 / 2500 [100%]  (Sampling)
Chain 1: 
Chain 1:  Elapsed Time: 3.578 seconds (Warm-up)
Chain 1:                11.07 seconds (Sampling)
Chain 1:                14.648 seconds (Total)
Chain 1: 

SAMPLING FOR MODEL 'stanmarg' NOW (CHAIN 2).
Chain 2: 
Chain 2: Gradient evaluation took 0.00032 seconds
Chain 2: 1000 transitions using 10 leapfrog steps per transition would take 3.2 seconds.
Chain 2: Adjust your expectations accordingly!
Chain 2: 
Chain 2: 
Chain 2: Iteration:    1 / 2500 [  0%]  (Warmup)
Chain 2: Iteration:  250 / 2500 [ 10%]  (Warmup)
Chain 2: Iteration:  500 / 2500 [ 20%]  (Warmup)
Chain 2: Iteration:  501 / 2500 [ 20%]  (Sampling)
Chain 2: Iteration:  750 / 2500 [ 30%]  (Sampling)
Chain 2: Iteration: 1000 / 2500 [ 40%]  (Sampling)
Chain 2: Iteration: 1250 / 2500 [ 50%]  (Sampling)
Chain 2: Iteration: 1500 / 2500 [ 60%]  (Sampling)
Chain 2: Iteration: 1750 / 2500 [ 70%]  (Sampling)
Chain 2: Iteration: 2000 / 2500 [ 80%]  (Sampling)
Chain 2: Iteration: 2250 / 2500 [ 90%]  (Sampling)
Chain 2: Iteration: 2500 / 2500 [100%]  (Sampling)
Chain 2: 
Chain 2:  Elapsed Time: 3.104 seconds (Warm-up)
Chain 2:                11.072 seconds (Sampling)
Chain 2:                14.176 seconds (Total)
Chain 2: 

SAMPLING FOR MODEL 'stanmarg' NOW (CHAIN 3).
Chain 3: 
Chain 3: Gradient evaluation took 0.000441 seconds
Chain 3: 1000 transitions using 10 leapfrog steps per transition would take 4.41 seconds.
Chain 3: Adjust your expectations accordingly!
Chain 3: 
Chain 3: 
Chain 3: Iteration:    1 / 2500 [  0%]  (Warmup)
Chain 3: Iteration:  250 / 2500 [ 10%]  (Warmup)
Chain 3: Iteration:  500 / 2500 [ 20%]  (Warmup)
Chain 3: Iteration:  501 / 2500 [ 20%]  (Sampling)
Chain 3: Iteration:  750 / 2500 [ 30%]  (Sampling)
Chain 3: Iteration: 1000 / 2500 [ 40%]  (Sampling)
Chain 3: Iteration: 1250 / 2500 [ 50%]  (Sampling)
Chain 3: Iteration: 1500 / 2500 [ 60%]  (Sampling)
Chain 3: Iteration: 1750 / 2500 [ 70%]  (Sampling)
Chain 3: Iteration: 2000 / 2500 [ 80%]  (Sampling)
Chain 3: Iteration: 2250 / 2500 [ 90%]  (Sampling)
Chain 3: Iteration: 2500 / 2500 [100%]  (Sampling)
Chain 3: 
Chain 3:  Elapsed Time: 3.596 seconds (Warm-up)
Chain 3:                11.83 seconds (Sampling)
Chain 3:                15.426 seconds (Total)
Chain 3: 
Computing post-estimation metrics (including lvs if requested)...
summary(fit.bayes, fit.measures = TRUE)
blavaan 0.5.10 ended normally after 2000 iterations

  Estimator                                      BAYES
  Optimization method                             MCMC
  Number of model parameters                        31

  Number of observations                            75

  Statistic                                 MargLogLik         PPP
  Value                                             NA       0.505

Parameter Estimates:


Latent Variables:
                   Estimate  Post.SD pi.lower pi.upper     Rhat    Prior     
  ind60 =~                                                                   
    x1                1.000                                                  
    x2                2.303    0.154    2.020    2.625    1.001  normal(0,10)
    x3                1.869    0.169    1.557    2.224    1.000  normal(0,10)
  dem60 =~                                                                   
    y1                1.000                                                  
    y2                1.354    0.215    0.965    1.805    1.003  normal(0,10)
    y3                1.114    0.172    0.804    1.483    1.000  normal(0,10)
    y4                1.364    0.187    1.043    1.781    1.001  normal(0,10)
  dem65 =~                                                                   
    y5                1.000                                                  
    y6                1.248    0.195    0.913    1.678    1.001  normal(0,10)
    y7                1.334    0.185    1.016    1.750    1.001  normal(0,10)
    y8                1.338    0.190    1.007    1.762    1.001  normal(0,10)

Regressions:
                   Estimate  Post.SD pi.lower pi.upper     Rhat    Prior     
  dem60 ~                                                                    
    ind60             1.434    0.395    0.680    2.222    1.000  normal(0,10)
  dem65 ~                                                                    
    ind60             0.511    0.229    0.052    0.964    1.000  normal(0,10)
    dem60             0.882    0.112    0.682    1.119    1.000  normal(0,10)

Covariances:
                   Estimate  Post.SD pi.lower pi.upper     Rhat    Prior     
 .y1 ~~                                                                      
   .y5                0.649    0.385   -0.034    1.489    1.001   lkj_corr(1)
 .y2 ~~                                                                      
   .y4                1.236    0.683   -0.017    2.684    1.000     beta(1,1)
   .y6                2.076    0.735    0.777    3.626    1.001     beta(1,1)
 .y3 ~~                                                                      
   .y7                0.763    0.626   -0.379    2.129    1.000   lkj_corr(1)
 .y4 ~~                                                                      
   .y8                0.240    0.467   -0.626    1.222    1.000     beta(1,1)
 .y6 ~~                                                                      
   .y8                1.341    0.555    0.340    2.542    1.000     beta(1,1)

Variances:
                   Estimate  Post.SD pi.lower pi.upper     Rhat    Prior     
   .x1                0.096    0.021    0.058    0.139    1.000 gamma(1,0.01)
   .x2                0.057    0.059    0.003    0.213    1.001 gamma(1,0.01)
   .x3                0.510    0.096    0.349    0.724    1.000 gamma(1,0.01)
   .y1                2.103    0.509    1.240    3.230    1.001 gamma(1,0.01)
   .y2                7.475    1.358    5.136   10.516    1.000 gamma(1,0.01)
   .y3                5.355    1.024    3.676    7.662    1.000 gamma(1,0.01)
   .y4                3.158    0.808    1.693    4.892    1.001 gamma(1,0.01)
   .y5                2.508    0.521    1.622    3.700    1.000 gamma(1,0.01)
   .y6                5.108    0.934    3.535    7.237    1.000 gamma(1,0.01)
   .y7                3.647    0.766    2.359    5.393    1.000 gamma(1,0.01)
   .y8                3.331    0.739    2.081    4.948    1.000 gamma(1,0.01)
    ind60             0.420    0.086    0.279    0.613    1.000 gamma(1,0.01)
   .dem60             3.604    0.938    2.010    5.686    1.001 gamma(1,0.01)
   .dem65             0.055    0.094    0.003    0.344    1.000 gamma(1,0.01)

16 マルチレベルSEM(Multilevel SEM)

連続データに対する、(random interceptsを導入した、切片[平均部分]に関する集団レベルの変動[集団差]を考慮した)2レベルのSEMが実装できる(本書の8.3節)。以下は、Demo.twolevelのデータを用いた例。

Demo.twolevelデータの構造の確認。変数clusterによって集団が区別されている。

head(Demo.twolevel)
          y1         y2         y3         y4         y5         y6         x1
1  0.2293216  1.3555232 -0.6911702  0.8028079 -0.3011085 -1.7260671  1.1739003
2  0.3085801 -1.8624397 -2.4179783  0.7659289  1.6750617  1.1764210 -1.0039958
3  0.2004934 -1.3400514  0.4376087  1.1974194  1.1951594  1.4988962 -0.4402545
4  1.0447982 -0.9624490 -0.4464898 -0.2027252 -0.4590574  1.1734061 -0.6253657
5  0.6881792 -0.4565633 -0.6422296  0.9900408  1.7662535  0.7944601 -0.8450025
6 -2.0687644 -0.5997856  0.3148418  0.6764432 -0.6519928  1.8405605 -0.7831784
           x2         x3         w1         w2 cluster
1 -0.62315173  0.6470414 -0.2479975 -0.4989800       1
2 -0.56689380  0.0201264 -0.2479975 -0.4989800       1
3 -2.13432572 -0.4591246 -0.2479975 -0.4989800       1
4 -0.33688869  1.2852093 -0.2479975 -0.4989800       1
5 -0.04229954  1.5598970 -0.2479975 -0.4989800       1
6 -0.22441996 -0.3814231 -2.3219338 -0.6910567       2
model.ml <- '
#集団内(個人)レベルの設定
  level: 1
    fw =~ y1 + y2 + y3
    fw ~ x1 + x2 + x3
#集団間(集団)レベルの設定
  level: 2
    fb =~ y1 + y2 + y3
    fb ~ w1 + w2
'

fit.ml <- sem(model = model.ml,
              data = Demo.twolevel,
              cluster = "cluster")
summary(fit.ml)
lavaan 0.6-21 ended normally after 36 iterations

  Estimator                                         ML
  Optimization method                           NLMINB
  Number of model parameters                        20

  Number of observations                          2500
  Number of clusters [cluster]                     200

Model Test User Model:
                                                      
  Test statistic                                 8.092
  Degrees of freedom                                10
  P-value (Chi-square)                           0.620

Parameter Estimates:

  Standard errors                             Standard
  Information                                 Observed
  Observed information based on                Hessian


Level 1 [within]:

Latent Variables:
                   Estimate  Std.Err  z-value  P(>|z|)
  fw =~                                               
    y1                1.000                           
    y2                0.774    0.034   22.671    0.000
    y3                0.734    0.033   22.355    0.000

Regressions:
                   Estimate  Std.Err  z-value  P(>|z|)
  fw ~                                                
    x1                0.510    0.023   22.037    0.000
    x2                0.407    0.022   18.273    0.000
    x3                0.205    0.021    9.740    0.000

Variances:
                   Estimate  Std.Err  z-value  P(>|z|)
   .y1                0.986    0.046   21.591    0.000
   .y2                1.066    0.039   27.271    0.000
   .y3                1.011    0.037   27.662    0.000
   .fw                0.546    0.040   13.539    0.000


Level 2 [cluster]:

Latent Variables:
                   Estimate  Std.Err  z-value  P(>|z|)
  fb =~                                               
    y1                1.000                           
    y2                0.717    0.052   13.824    0.000
    y3                0.587    0.048   12.329    0.000

Regressions:
                   Estimate  Std.Err  z-value  P(>|z|)
  fb ~                                                
    w1                0.165    0.079    2.093    0.036
    w2                0.131    0.076    1.715    0.086

Intercepts:
                   Estimate  Std.Err  z-value  P(>|z|)
   .y1                0.024    0.075    0.327    0.743
   .y2               -0.016    0.060   -0.269    0.788
   .y3               -0.042    0.054   -0.777    0.437

Variances:
                   Estimate  Std.Err  z-value  P(>|z|)
   .y1                0.058    0.047    1.213    0.225
   .y2                0.120    0.031    3.825    0.000
   .y3                0.149    0.028    5.319    0.000
   .fb                0.899    0.118    7.592    0.000

17 媒介分析(横断研究に基づく間接効果の推定;本書の8.1.1節)

従属変数Y、独立変数X、媒介変数Mの3変数がある状況を考える。数値例として、これら3変数を含むデータセットを作り、XがYに直接影響を与え(=直接効果)、XがMを媒介してYに影響を与えること(=間接効果)を仮定した媒介モデルを当てはめる:

set.seed(1234)
X <- rnorm(100)
M <- 0.5 * X + rnorm(100)
Y <- 0.7 * M + rnorm(100)
Data.med <- data.frame(X = X, Y = Y, M = M)

model.med <- '
  Y ~ c*X
  M ~ a*X
  Y ~ b*M

  ab := a*b
  total := c + (a*b)
'

fit.med <- sem(model.med, data = Data.med)
summary(fit.med)
lavaan 0.6-21 ended normally after 1 iteration

  Estimator                                         ML
  Optimization method                           NLMINB
  Number of model parameters                         5

  Number of observations                           100

Model Test User Model:
                                                      
  Test statistic                                 0.000
  Degrees of freedom                                 0

Parameter Estimates:

  Standard errors                             Standard
  Information                                 Expected
  Information saturated (h1) model          Structured

Regressions:
                   Estimate  Std.Err  z-value  P(>|z|)
  Y ~                                                 
    X          (c)    0.036    0.104    0.348    0.728
  M ~                                                 
    X          (a)    0.474    0.103    4.613    0.000
  Y ~                                                 
    M          (b)    0.788    0.092    8.539    0.000

Variances:
                   Estimate  Std.Err  z-value  P(>|z|)
   .Y                 0.898    0.127    7.071    0.000
   .M                 1.054    0.149    7.071    0.000

Defined Parameters:
                   Estimate  Std.Err  z-value  P(>|z|)
    ab                0.374    0.092    4.059    0.000
    total             0.410    0.125    3.287    0.001

ブートストラップ標準誤差を使うなら次のようにする。

fit.med.boot <- sem(model.med,
                    data = Data.med,
                    se = "bootstrap",
                    bootstrap = 1000)
summary(fit.med.boot)

18 修正指標(Modification Indices; 本書の2.5.2節)

分析モデルにおいて(0 などに)固定された全ての母数について計算される、自由度1 のカイ二乗検定に基づく指標。もし固定せずに推定した場合に検定統計量がどの程度低下するか を近似的に求めている。

modindices(fit.cfa, sort = TRUE, maximum.number = 5)
       lhs op rhs     mi    epc sepc.lv sepc.all sepc.nox
30  visual =~  x9 36.411  0.577   0.519    0.515    0.515
76      x7 ~~  x8 34.145  0.536   0.536    0.859    0.859
28  visual =~  x7 18.631 -0.422  -0.380   -0.349   -0.349
78      x8 ~~  x9 14.946 -0.423  -0.423   -0.805   -0.805
33 textual =~  x3  9.151 -0.272  -0.269   -0.238   -0.238

因子負荷(“=~”に相当する母数)だけに絞って確認するには次のようにする。

mi <- modindices(fit.cfa)
mi[mi$op == "=~", ]
       lhs op rhs     mi    epc sepc.lv sepc.all sepc.nox
25  visual =~  x4  1.211  0.077   0.069    0.059    0.059
26  visual =~  x5  7.441 -0.210  -0.189   -0.147   -0.147
27  visual =~  x6  2.843  0.111   0.100    0.092    0.092
28  visual =~  x7 18.631 -0.422  -0.380   -0.349   -0.349
29  visual =~  x8  4.295 -0.210  -0.189   -0.187   -0.187
30  visual =~  x9 36.411  0.577   0.519    0.515    0.515
31 textual =~  x1  8.903  0.350   0.347    0.297    0.297
32 textual =~  x2  0.017 -0.011  -0.011   -0.010   -0.010
33 textual =~  x3  9.151 -0.272  -0.269   -0.238   -0.238
34 textual =~  x7  0.098 -0.021  -0.021   -0.019   -0.019
35 textual =~  x8  3.359 -0.121  -0.120   -0.118   -0.118
36 textual =~  x9  4.796  0.138   0.137    0.136    0.136
37   speed =~  x1  0.014  0.024   0.015    0.013    0.013
38   speed =~  x2  1.580 -0.198  -0.123   -0.105   -0.105
39   speed =~  x3  0.716  0.136   0.084    0.075    0.075
40   speed =~  x4  0.003 -0.005  -0.003   -0.003   -0.003
41   speed =~  x5  0.201 -0.044  -0.027   -0.021   -0.021
42   speed =~  x6  0.273  0.044   0.027    0.025    0.025

等値制約を外した場合の検討には lavTestScore() が有用である。

19 当てはめた分析モデルに関する情報の抽出

19.1 coef()

母数の点推定値の抽出。

coef(fit.sem)
   ind60=~x2    ind60=~x3    dem60=~y2    dem60=~y3    dem60=~y4    dem65=~y6 
       2.180        1.819        1.257        1.058        1.265        1.186 
   dem65=~y7    dem65=~y8  dem60~ind60  dem65~ind60  dem65~dem60       y1~~y5 
       1.280        1.266        1.483        0.572        0.837        0.624 
      y2~~y4       y2~~y6       y3~~y7       y4~~y8       y6~~y8       x1~~x1 
       1.313        2.153        0.795        0.348        1.356        0.082 
      x2~~x2       x3~~x3       y1~~y1       y2~~y2       y3~~y3       y4~~y4 
       0.120        0.467        1.891        7.373        5.067        3.148 
      y5~~y5       y6~~y6       y7~~y7       y8~~y8 ind60~~ind60 dem60~~dem60 
       2.351        4.954        3.431        3.254        0.448        3.956 
dem65~~dem65 
       0.172 

19.2 parameterEstimates()

母数の点推定値やSE, 検定統計量等の抽出。

parameterEstimates(fit.sem)
     lhs op   rhs   est    se      z pvalue ci.lower ci.upper
1  ind60 =~    x1 1.000 0.000     NA     NA    1.000    1.000
2  ind60 =~    x2 2.180 0.139 15.742  0.000    1.909    2.452
3  ind60 =~    x3 1.819 0.152 11.967  0.000    1.521    2.116
4  dem60 =~    y1 1.000 0.000     NA     NA    1.000    1.000
5  dem60 =~    y2 1.257 0.182  6.889  0.000    0.899    1.614
6  dem60 =~    y3 1.058 0.151  6.987  0.000    0.761    1.354
7  dem60 =~    y4 1.265 0.145  8.722  0.000    0.981    1.549
8  dem65 =~    y5 1.000 0.000     NA     NA    1.000    1.000
9  dem65 =~    y6 1.186 0.169  7.024  0.000    0.855    1.517
10 dem65 =~    y7 1.280 0.160  8.002  0.000    0.966    1.593
11 dem65 =~    y8 1.266 0.158  8.007  0.000    0.956    1.576
12 dem60  ~ ind60 1.483 0.399  3.715  0.000    0.701    2.265
13 dem65  ~ ind60 0.572 0.221  2.586  0.010    0.139    1.006
14 dem65  ~ dem60 0.837 0.098  8.514  0.000    0.645    1.030
15    y1 ~~    y5 0.624 0.358  1.741  0.082   -0.079    1.326
16    y2 ~~    y4 1.313 0.702  1.871  0.061   -0.063    2.689
17    y2 ~~    y6 2.153 0.734  2.934  0.003    0.715    3.591
18    y3 ~~    y7 0.795 0.608  1.308  0.191   -0.396    1.986
19    y4 ~~    y8 0.348 0.442  0.787  0.431   -0.519    1.215
20    y6 ~~    y8 1.356 0.568  2.386  0.017    0.242    2.470
21    x1 ~~    x1 0.082 0.019  4.184  0.000    0.043    0.120
22    x2 ~~    x2 0.120 0.070  1.718  0.086   -0.017    0.256
23    x3 ~~    x3 0.467 0.090  5.177  0.000    0.290    0.643
24    y1 ~~    y1 1.891 0.444  4.256  0.000    1.020    2.762
25    y2 ~~    y2 7.373 1.374  5.366  0.000    4.680   10.066
26    y3 ~~    y3 5.067 0.952  5.324  0.000    3.202    6.933
27    y4 ~~    y4 3.148 0.739  4.261  0.000    1.700    4.596
28    y5 ~~    y5 2.351 0.480  4.895  0.000    1.410    3.292
29    y6 ~~    y6 4.954 0.914  5.419  0.000    3.162    6.746
30    y7 ~~    y7 3.431 0.713  4.814  0.000    2.034    4.829
31    y8 ~~    y8 3.254 0.695  4.685  0.000    1.893    4.615
32 ind60 ~~ ind60 0.448 0.087  5.173  0.000    0.279    0.618
33 dem60 ~~ dem60 3.956 0.921  4.295  0.000    2.151    5.762
34 dem65 ~~ dem65 0.172 0.215  0.803  0.422   -0.249    0.593

19.3 standardizedSolution()

標準化解の抽出。

standardizedSolution(fit.sem)
     lhs op   rhs est.std    se      z pvalue ci.lower ci.upper
1  ind60 =~    x1   0.920 0.023 40.005  0.000    0.875    0.965
2  ind60 =~    x2   0.973 0.016 59.110  0.000    0.941    1.005
3  ind60 =~    x3   0.872 0.031 28.081  0.000    0.811    0.933
4  dem60 =~    y1   0.850 0.042 20.200  0.000    0.768    0.933
5  dem60 =~    y2   0.717 0.064 11.157  0.000    0.591    0.843
6  dem60 =~    y3   0.722 0.063 11.439  0.000    0.599    0.846
7  dem60 =~    y4   0.846 0.043 19.462  0.000    0.761    0.931
8  dem65 =~    y5   0.808 0.048 16.847  0.000    0.714    0.902
9  dem65 =~    y6   0.746 0.058 12.924  0.000    0.633    0.859
10 dem65 =~    y7   0.824 0.045 18.392  0.000    0.736    0.911
11 dem65 =~    y8   0.828 0.045 18.477  0.000    0.740    0.916
12 dem60  ~ ind60   0.447 0.103  4.323  0.000    0.244    0.649
13 dem65  ~ ind60   0.182 0.070  2.587  0.010    0.044    0.320
14 dem65  ~ dem60   0.885 0.051 17.405  0.000    0.786    0.985
15    y1 ~~    y5   0.296 0.140  2.117  0.034    0.022    0.570
16    y2 ~~    y4   0.273 0.120  2.265  0.024    0.037    0.508
17    y2 ~~    y6   0.356 0.097  3.690  0.000    0.167    0.545
18    y3 ~~    y7   0.191 0.133  1.432  0.152   -0.070    0.452
19    y4 ~~    y8   0.109 0.131  0.832  0.405   -0.147    0.365
20    y6 ~~    y8   0.338 0.110  3.079  0.002    0.123    0.553
21    x1 ~~    x1   0.154 0.042  3.637  0.000    0.071    0.237
22    x2 ~~    x2   0.053 0.032  1.661  0.097   -0.010    0.116
23    x3 ~~    x3   0.239 0.054  4.419  0.000    0.133    0.346
24    y1 ~~    y1   0.277 0.072  3.865  0.000    0.136    0.417
25    y2 ~~    y2   0.486 0.092  5.269  0.000    0.305    0.666
26    y3 ~~    y3   0.478 0.091  5.242  0.000    0.299    0.657
27    y4 ~~    y4   0.285 0.073  3.875  0.000    0.141    0.429
28    y5 ~~    y5   0.347 0.078  4.478  0.000    0.195    0.499
29    y6 ~~    y6   0.443 0.086  5.149  0.000    0.275    0.612
30    y7 ~~    y7   0.322 0.074  4.359  0.000    0.177    0.466
31    y8 ~~    y8   0.315 0.074  4.242  0.000    0.169    0.460
32 ind60 ~~ ind60   1.000 0.000     NA     NA    1.000    1.000
33 dem60 ~~ dem60   0.800 0.092  8.669  0.000    0.619    0.981
34 dem65 ~~ dem65   0.039 0.048  0.808  0.419   -0.056    0.134

19.4 fitted()

推定されたモデルの共分散構造(および平均構造)の出力。

fitted(fit.cfa)
$cov
      x1    x2    x3    x4    x5    x6    x7    x8    x9
x1 1.358                                                
x2 0.448 1.382                                          
x3 0.590 0.327 1.275                                    
x4 0.408 0.226 0.298 1.351                              
x5 0.454 0.252 0.331 1.090 1.660                        
x6 0.378 0.209 0.276 0.907 1.010 1.196                  
x7 0.262 0.145 0.191 0.173 0.193 0.161 1.183            
x8 0.309 0.171 0.226 0.205 0.228 0.190 0.453 1.022      
x9 0.284 0.157 0.207 0.188 0.209 0.174 0.415 0.490 1.015

19.5 resid() / residuals()

残差平均ベクトル・残差共分散行列・残差相関行列(本書の2.5.2節)

resid(fit.cfa, type = "standardized")
$type
[1] "standardized"

$cov
       x1     x2     x3     x4     x5     x6     x7     x8     x9
x1  0.000                                                        
x2 -1.996  0.000                                                 
x3 -0.997  2.689  0.000                                          
x4  2.679 -0.284 -1.899  0.000                                   
x5 -0.359 -0.591 -4.157  1.545  0.000                            
x6  2.155  0.681 -0.711 -2.588  0.942  0.000                     
x7 -3.773 -3.654 -1.858  0.865 -0.842 -0.326  0.000              
x8 -1.380 -1.119 -0.300 -2.021 -1.099 -0.641  4.823  0.000       
x9  4.077  1.606  3.518  1.225  1.701  1.423 -2.325 -4.132  0.000

19.6 vcov()

推定値の漸近共分散行列(本書の2.4.3節)

vcov(fit.sem)
             i60=~2 i60=~3 d60=~2 d60=~3 d60=~4 d65=~6 d65=~7 d65=~8 d60~60
ind60=~x2     0.019                                                        
ind60=~x3     0.010  0.023                                                 
dem60=~y2     0.000  0.000  0.033                                          
dem60=~y3     0.000  0.000  0.008  0.023                                   
dem60=~y4     0.000  0.000  0.015  0.009  0.021                            
dem65=~y6     0.000  0.000  0.010  0.003  0.004  0.028                     
dem65=~y7     0.000  0.000  0.004  0.006  0.004  0.013  0.026              
dem65=~y8     0.000  0.000  0.004  0.003  0.005  0.017  0.013  0.025       
dem60~ind60   0.009  0.007 -0.012 -0.010 -0.012 -0.004 -0.005 -0.005  0.159
dem65~ind60   0.003  0.003 -0.001 -0.002 -0.001 -0.005 -0.006 -0.005 -0.011
dem65~dem60   0.000  0.000  0.004  0.004  0.004 -0.006 -0.007 -0.007 -0.003
y1~~y5        0.000  0.000  0.009  0.007  0.010  0.009  0.010  0.010 -0.007
y2~~y4        0.000  0.000 -0.025 -0.001 -0.017 -0.002 -0.002 -0.002  0.010
y2~~y6        0.000  0.000 -0.009  0.000  0.002 -0.007  0.000  0.003  0.000
y3~~y7        0.000  0.000  0.000 -0.014  0.000  0.000 -0.011  0.000  0.003
y4~~y8        0.000  0.000 -0.001 -0.001 -0.010 -0.001 -0.001 -0.010  0.005
y6~~y8        0.000  0.000  0.000  0.000  0.000 -0.013 -0.001 -0.011  0.000
x1~~x1        0.001  0.000  0.000  0.000  0.000  0.000  0.000  0.000  0.000
x2~~x2       -0.004  0.000  0.000  0.000  0.000  0.000  0.000  0.000  0.000
x3~~x3        0.001 -0.001  0.000  0.000  0.000  0.000  0.000  0.000  0.000
y1~~y1        0.000  0.000  0.015  0.013  0.017  0.006  0.007  0.007 -0.012
y2~~y2        0.000  0.000 -0.041 -0.001 -0.010 -0.008 -0.002  0.001  0.009
y3~~y3        0.000  0.000 -0.001 -0.023  0.000 -0.001 -0.007  0.000  0.006
y4~~y4        0.000  0.000 -0.015 -0.002 -0.028 -0.003 -0.002 -0.008  0.014
y5~~y5        0.000  0.000  0.005  0.005  0.006  0.012  0.013  0.014 -0.005
y6~~y6        0.000  0.000 -0.004  0.000  0.001 -0.020 -0.001 -0.006  0.000
y7~~y7        0.000  0.000 -0.001 -0.007  0.000 -0.001 -0.019  0.000  0.002
y8~~y8        0.000  0.000 -0.001  0.000 -0.006 -0.009 -0.001 -0.021  0.003
ind60~~ind60 -0.006 -0.004  0.000  0.000  0.000  0.000  0.000  0.000 -0.003
dem60~~dem60  0.001  0.000 -0.063 -0.053 -0.064 -0.021 -0.023 -0.023  0.039
dem65~~dem65  0.000  0.000  0.000 -0.001  0.000 -0.002 -0.003 -0.002  0.001
             dm65~n60 dm65~d60 y1~~y5 y2~~y4 y2~~y6 y3~~y7 y4~~y8 y6~~y8 x1~~x1
ind60=~x2                                                                      
ind60=~x3                                                                      
dem60=~y2                                                                      
dem60=~y3                                                                      
dem60=~y4                                                                      
dem65=~y6                                                                      
dem65=~y7                                                                      
dem65=~y8                                                                      
dem60~ind60                                                                    
dem65~ind60     0.049                                                          
dem65~dem60    -0.006    0.010                                                 
y1~~y5          0.000   -0.003  0.128                                          
y2~~y4         -0.013    0.005 -0.029  0.493                                   
y2~~y6          0.000    0.000  0.000  0.020  0.538                            
y3~~y7          0.006   -0.003 -0.012 -0.028  0.000  0.369                     
y4~~y8          0.003   -0.002 -0.017 -0.010 -0.064 -0.016  0.196              
y6~~y8          0.002    0.003 -0.018 -0.034  0.009 -0.020  0.000  0.323       
x1~~x1          0.000    0.000  0.000  0.000  0.000  0.000  0.000  0.000  0.000
x2~~x2          0.000    0.000  0.000  0.000  0.000  0.000  0.000  0.000 -0.001
x3~~x3          0.000    0.000  0.000  0.000  0.000  0.000  0.000  0.000  0.000
y1~~y1         -0.008    0.007  0.086 -0.025  0.000 -0.023 -0.025 -0.003  0.000
y2~~y2         -0.013    0.005 -0.027  0.486  0.504 -0.025 -0.060 -0.034  0.000
y3~~y3         -0.002    0.001 -0.018 -0.019  0.000  0.230 -0.018 -0.004  0.000
y4~~y4         -0.011    0.003 -0.044  0.273 -0.037 -0.041  0.113 -0.010  0.000
y5~~y5         -0.002   -0.005  0.092 -0.009  0.000 -0.018 -0.018 -0.008  0.000
y6~~y6          0.002    0.003 -0.015 -0.028  0.321 -0.017 -0.040  0.281  0.000
y7~~y7          0.005    0.002 -0.024 -0.001  0.000  0.166 -0.023 -0.005  0.000
y8~~y8          0.004    0.002 -0.030 -0.015 -0.035 -0.031  0.105  0.227  0.000
ind60~~ind60   -0.001    0.000  0.000  0.000  0.000  0.000  0.000  0.000  0.000
dem60~~dem60    0.032   -0.030 -0.042  0.013  0.000  0.017  0.018  0.004  0.000
dem65~~dem65    0.008   -0.004  0.012 -0.028  0.000  0.017  0.008 -0.026  0.000
             x2~~x2 x3~~x3 y1~~y1 y2~~y2 y3~~y3 y4~~y4 y5~~y5 y6~~y6 y7~~y7
ind60=~x2                                                                  
ind60=~x3                                                                  
dem60=~y2                                                                  
dem60=~y3                                                                  
dem60=~y4                                                                  
dem65=~y6                                                                  
dem65=~y7                                                                  
dem65=~y8                                                                  
dem60~ind60                                                                
dem65~ind60                                                                
dem65~dem60                                                                
y1~~y5                                                                     
y2~~y4                                                                     
y2~~y6                                                                     
y3~~y7                                                                     
y4~~y8                                                                     
y6~~y8                                                                     
x1~~x1                                                                     
x2~~x2        0.005                                                        
x3~~x3       -0.002  0.008                                                 
y1~~y1        0.000  0.000  0.198                                          
y2~~y2        0.000  0.000 -0.021  1.888                                   
y3~~y3        0.000  0.000 -0.021 -0.016  0.906                            
y4~~y4        0.000  0.000 -0.046  0.131 -0.034  0.546                     
y5~~y5        0.000  0.000  0.044 -0.007 -0.012 -0.024  0.231              
y6~~y6        0.000  0.000 -0.001  0.140 -0.002 -0.025 -0.005  0.836       
y7~~y7        0.000  0.000 -0.015  0.003  0.067 -0.019 -0.016 -0.003  0.508
y8~~y8        0.000  0.000 -0.019 -0.034 -0.015  0.038 -0.019  0.102 -0.020
ind60~~ind60  0.001  0.000  0.000  0.000  0.000  0.000  0.000  0.000  0.000
dem60~~dem60 -0.001  0.000 -0.086  0.010  0.013  0.027 -0.027  0.002  0.012
dem65~~dem65  0.000  0.000 -0.010 -0.029 -0.004 -0.022 -0.010 -0.025 -0.019
             y8~~y8 i60~~6 d60~~6 d65~~6
ind60=~x2                               
ind60=~x3                               
dem60=~y2                               
dem60=~y3                               
dem60=~y4                               
dem65=~y6                               
dem65=~y7                               
dem65=~y8                               
dem60~ind60                             
dem65~ind60                             
dem65~dem60                             
y1~~y5                                  
y2~~y4                                  
y2~~y6                                  
y3~~y7                                  
y4~~y8                                  
y6~~y8                                  
x1~~x1                                  
x2~~x2                                  
x3~~x3                                  
y1~~y1                                  
y2~~y2                                  
y3~~y3                                  
y4~~y4                                  
y5~~y5                                  
y6~~y6                                  
y7~~y7                                  
y8~~y8        0.482                     
ind60~~ind60  0.000  0.008              
dem60~~dem60  0.015  0.000  0.849       
dem65~~dem65 -0.023  0.000  0.013  0.046

19.7 AIC() / BIC()

分析モデルに関する情報量規準の出力(本書の2.5.1節)

AIC(fit.sem)
[1] 3157.582
BIC(fit.sem)
[1] 3229.424

19.8 fitMeasures()

適合度指標の出力(本書の2.5.1節)

fitMeasures(fit.cfa, c("cfi", "rmsea", "srmr"))
  cfi rmsea  srmr 
0.931 0.092 0.065 

19.9 lavInspect() / inspect()

当てはめたモデルオブジェクトの中身を調べたり、内部情報を取り出したりするための関数。例えば、設定した分析モデルの係数行列や分散共分散行列、母数の設定の具体を確認したいときに利用できる。

lavInspect(fit.cfa)
$lambda
   visual textul speed
x1      0      0     0
x2      1      0     0
x3      2      0     0
x4      0      0     0
x5      0      3     0
x6      0      4     0
x7      0      0     0
x8      0      0     5
x9      0      0     6

$theta
   x1 x2 x3 x4 x5 x6 x7 x8 x9
x1  7                        
x2  0  8                     
x3  0  0  9                  
x4  0  0  0 10               
x5  0  0  0  0 11            
x6  0  0  0  0  0 12         
x7  0  0  0  0  0  0 13      
x8  0  0  0  0  0  0  0 14   
x9  0  0  0  0  0  0  0  0 15

$psi
        visual textul speed
visual      16             
textual     19     17      
speed       20     21    18
inspect(fit.cfa)
$lambda
   visual textul speed
x1      0      0     0
x2      1      0     0
x3      2      0     0
x4      0      0     0
x5      0      3     0
x6      0      4     0
x7      0      0     0
x8      0      0     5
x9      0      0     6

$theta
   x1 x2 x3 x4 x5 x6 x7 x8 x9
x1  7                        
x2  0  8                     
x3  0  0  9                  
x4  0  0  0 10               
x5  0  0  0  0 11            
x6  0  0  0  0  0 12         
x7  0  0  0  0  0  0 13      
x8  0  0  0  0  0  0  0 14   
x9  0  0  0  0  0  0  0  0 15

$psi
        visual textul speed
visual      16             
textual     19     17      
speed       20     21    18

初期値を確認したい場合は次のようにする。

lavInspect(fit.cfa, what = "start")
$lambda
   visual textul speed
x1  1.000  0.000 0.000
x2  0.778  0.000 0.000
x3  1.107  0.000 0.000
x4  0.000  1.000 0.000
x5  0.000  1.133 0.000
x6  0.000  0.924 0.000
x7  0.000  0.000 1.000
x8  0.000  0.000 1.225
x9  0.000  0.000 0.854

$theta
      x1    x2    x3    x4    x5    x6    x7    x8    x9
x1 0.679                                                
x2 0.000 0.691                                          
x3 0.000 0.000 0.637                                    
x4 0.000 0.000 0.000 0.675                              
x5 0.000 0.000 0.000 0.000 0.830                        
x6 0.000 0.000 0.000 0.000 0.000 0.598                  
x7 0.000 0.000 0.000 0.000 0.000 0.000 0.592            
x8 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.511      
x9 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.508

$psi
        visual textul speed
visual    0.05             
textual   0.00   0.05      
speed     0.00   0.00  0.05

20 ESEM と EFA

公式 tutorial では、モデル内に含まれる多重指標(e.g.,心理尺度)に関する因子構造が事前に明確でない場合に、因子の回転を用いて構造を探索しながら(モデル全体を)推定する枠組みとして ESEM / EFA (探索的構造方程式モデリング、探索的因子分析)が紹介されている。具体例は公式 tutorial の ESEM / EFA セクションを参照されたい。

21 補足: read.csv() によるデータの読み込み

本資料のCLPMの例では、別添の semexample.csv (人工データ)を利用する。まず利用者はこのファイルをダウンロードし、lavaan_tutorial_ja.Rmd と同じフォルダに保存する。そうすれば、read.csv() を用いて次のように読み込める。

mydata <- read.csv("semexample.csv", header = TRUE)
head(mydata)

ここで、"semexample.csv" は読み込むCSVファイル名を表しており、header = TRUE は 1 行目を変数名として扱うことを意味する。head(mydata) は、読み込んだデータの冒頭数行を確認するための関数である。

21.1 ファイルパスの考え方

上のような書き方ができるのは、semexample.csvlavaan_tutorial_ja.Rmd と同じフォルダに置いているためである。別の場所に保存した場合には、その保存場所に応じてパスを書き換える必要がある。

  • 相対パス:現在の作業フォルダを基準に書く方法
  • 絶対パス:PC 内の保存場所を最初から最後までフルに書く方法

21.2 絶対パスを用いる例

Windows の場合、たとえば semexample.csv をデスクトップに保存したなら、次のように書ける。

mydata <- read.csv("C:/Users/YourName/Desktop/semexample.csv", header = TRUE)
head(mydata)

宇佐美のPCの場合、以下の通り:

mydata <- read.csv("C:/Users/Satoshi Usami/Desktop/semexample.csv", header = TRUE)
head(mydata)
  ID Sex GHQ1 GHQ2 GHQ3 GHQ4 GHQ5 GHQ6 sleep1 sleep2 sleep3 sleep4 sleep5
1  1  NA    0   NA   NA   NA   NA   NA     NA     NA     NA     NA     NA
2  2   1    3   NA   NA   NA   NA   NA    9.0     NA     NA     NA     NA
3  3   1    0   NA   NA   NA   NA   NA    8.5     NA     NA     NA     NA
4  4  NA    0   NA   NA   NA   NA   NA    9.0     NA     NA     NA     NA
5  5  NA    0   NA   NA   NA   NA   NA    8.0     NA     NA     NA     NA
6  6  NA    0   NA   NA   NA   NA   NA    8.5     NA     NA     NA     NA
  sleep6 bed1 bed2 bed3 bed4 bed5 bed6
1     NA   NA   NA   NA   NA   NA   NA
2     NA 21.0   NA   NA   NA   NA   NA
3     NA 21.0   NA   NA   NA   NA   NA
4     NA 21.0   NA   NA   NA   NA   NA
5     NA 21.5   NA   NA   NA   NA   NA
6     NA 21.5   NA   NA   NA   NA   NA

Mac の場合は、たとえばホームディレクトリの Documents フォルダ内に保存したとき、次のような形になる。

mydata <- read.csv("/Users/YourName/Documents/semexample.csv",
                   header = TRUE)
head(mydata)

21.3 日本語ファイル名や文字コードに関する指定

semexample.csv をそのまま利用する限り通常は特別な指定は不要だが、環境によって文字コードの指定が必要になる場合には、必要に応じて fileEncoding を指定する。

mydata <- read.csv("semexample.csv",
                   header = TRUE,
                   fileEncoding = "UTF-8")
head(mydata)

Windows の環境では、CSV ファイルが CP932Shift-JIS 系で保存されていることもあり、その場合には次のような指定が必要になることがある。

mydata <- read.csv("semexample.csv",
                   header = TRUE,
                   fileEncoding = "CP932")
head(mydata)

21.4 プロジェクトを作成した場合の例

RStudio では、Project を作成して作業する方法が便利である。Project は作業ディレクトリと結び付けられており、分析ごとに 1 つの Project を作成して、その中にデータ・スクリプト・出力ファイルをまとめて管理する方法が勧められている。Project の基本的な説明や作成方法については、Posit 公式ガイドの Projects in RStudio IDE および Get Started を参照のこと。

たとえば、次のようなフォルダ構成を考える。

project_folder/
  lavaan_tutorial_ja.Rmd
  semexample.csv

このとき、lavaan_tutorial_ja.Rmdsemexample.csv が同じ Project 内の同じフォルダにあれば、次のように書ける。

mydata <- read.csv("semexample.csv", header = TRUE)
head(mydata)

この方法の利点は、他の PC にファイル一式を移しても、フォルダ構成が同じであればそのまま動きやすいことである。授業資料や共同研究でファイルを共有する場合にも有用である。

21.5 現在の作業フォルダの確認

ファイルが見つからない場合は、現在どのフォルダを基準にしているかを確認するとよい。これには getwd() を使う。

getwd()

これにより、現在の作業フォルダが表示される。相対パスでうまくいかない場合は、この出力を確認して、どこを基準にファイルを探しているのかを把握するとよい。

22 引用文献

Bollen, K.A. (1989). Structural Equations with Latent Variables. New York: Wiley.

Merkle, E.C., & Rosseel, Y. (2018). blavaan: Bayesian structural equation models via parameter expansion. Journal of Statistical Software, 85, 1–30. https://www.jstatsoft.org/article/view/v085i04.

R Core Team (2024). R: A language and environment for statistical computing. R Foundation for Statistical Computing, Vienna, Austria. https://www.R-project.org/

Rosseel, Y. (2012). lavaan: An R Package for Structural Equation Modeling. Journal of Statistical Software, 48(2), 1–36.

宇佐美慧 (2026). 縦断データ分析:構造方程式モデリングによるアプローチ 東京大学出版会