

資料集
lavaan包的內建資料集PoliticalDemocracy,記錄了發展中國家政治民主和工業化的多種指標,共涉及75個物件(國家)和11個變數(政治民主或工業化指標)。
library(lavaan)
#資料集,詳情 ?PoliticalDemocracy
dat <- PoliticalDemocracy
head(dat)

對於11個變數,分別為:
y1,1960年新聞自由的專家評級;
y2,1960年政治反對派自由;
y3,1960年選舉的公正性;
y4,1960年當選立法機關的效力;
y5,1965年新聞自由的專家評級;
y6,1965年的政治反對自由;
y7,1965年選舉的公正性;
y8,1965年當選立法機關的效力;
x1,1960年人均國民生產總值(GNP);
x2,1960年人均能源能耗;
x3,1960年工業勞動人口的百分比。
這些變數之間具有如下的結構關係:
y1、y2、y3、y4代表了1960年的政治民主指標(dem60 ~ y1 + y2 + y3 + y4);
y5、y6、y7、y8代表了1965年的政治民主指標(dem65 ~ y5 + y6 + y7 + y8);
x1、x2、x3代表了1960年的工業化指標(ind60 ~ x1 + x2 + x3)。
那麼可知,dem60、dem65和ind60即代表了資料集的潛在變數,反映了這些發展中國家的政治民主和工業化發展程度。不可否認,一個國家的政治民主發展水平和工業化發展水平之間存在密不可分的關係,因此我們可以推斷,潛在變數之間還存在這樣的結構關係:
dem60 ~ ind60;
dem65 ~ ind60 + dem60
綜合考慮(潛)變數和(潛)變數之間的因果關係,整個資料集可以建立如下模型:

那麼,我們的推斷是否合理呢?此時潛變數結構方程模型是個不錯的方法。
注:作為示例演示,以下操作忽略資料集是否滿足多元正態性。
lavaan的驗證性因子分析評估潛在變數
首先就是需要確定資料集中的潛在變數是否合適,我們使用驗證性因子分析來實現驗證過程。
##驗證性因子分析(CFA)
#假定潛變數
cfa_model <- '
#latent variables
ind60 =~ x1 + x2 + x3
dem60 =~ y1 + y2 + y3 + y4
dem65 =~ y5 + y6 + y7 + y8
#note that lavaan automatically includes latent covariances
#but we can add here anyway to be explicit
ind60 ~~ dem60
ind60 ~~ dem65
dem60 ~~ dem65
'
#執行 CFA,詳情 ?cfa
cfa_fit <- cfa(model = cfa_model, data = dat)
summary(cfa_fit)
#模型擬合度,詳情 ?fitmeasures
fitmeasures(cfa_fit, c('chisq', 'rmsea', 'cfi', 'aic'))

透過在驗證性因子分析中,預指定的潛在變數(潛在因子)和觀測變數間的關係,並比較預測矩陣和觀測矩陣的差異量化模型擬合度。幾種常見的評估擬合度的指標顯示,CFI值反映出模型擬合良好,但RMSEA值表明模型擬合度一般。如果期望改進模型(尋找更合適的潛在變數和觀測變數間關係),可參考前文驗證性因子分析中的方法(該文也有對CFI、RMSEA值等概念的描述,對於AIC的理解可參考該文)。
這裡作為示例,暫且認為潛在變數是合適的,即原始定義的模型是合理的。然後下一步進行潛變數結構方程建模。
lavaan的結構方程建模
我們將本篇一開始時假定的代表政治民主和工業化變數間因果關係結構輸入,進行SEM建模。
SEM建模與模型評估
lavaan中,結構方程模型透過sem()函式實現。
##結構方程模型(SEM)
#假定變數結構
sem_model <- '
#latent variables
ind60 =~ x1 + x2 + x3
dem60 =~ y1 + y2 + y3 + y4
dem65 =~ y5 + y6 + y7 + y8
#regressions
dem60 ~ ind60
dem65 ~ ind60 + dem60
#residual covariances
y1 ~~ y5
y2 ~~ y4 + y6
y3 ~~ y7
y4 ~~ y8
y6 ~~ y8
'
#執行 SEM,詳情 ?sem
#這裡同時包含了觀測變數以及潛在變數的建模,即一個潛變數結構方程建模
#如果 model 中不存在潛在變數,僅為觀測變數,則該模型即為一種普通的路徑分析模型
sem_fit <- sem(model = sem_model, data = dat, se = 'bootstrap', bootstrap = 100)
summary(sem_fit, standardized = TRUE, fit.measures = TRUE, rsquare = TRUE)
#模型擬合度,詳情 ?fitmeasures
fitmeasures(sem_fit, c('chisq', 'rmsea', 'cfi', 'aic'))
#路徑圖展示模型中的因果關係,詳情 ?semPaths
#例如,連線中的數值用於反映標準化的迴歸係數(標準化的引數估計值)
library(semPlot)
semPaths(sem_fit, what = 'est', layout = 'tree', residuals = FALSE,
edge.label.cex = 1, edge.color = 'red')
根據模型擬合檢驗統計量(卡方統計量,Chi-square)、RMSEA等的p值(p>0.05)可知,模型擬合優度的提升空間幾乎很低了,模型中的有效路徑基本識別完全。

且模型擬合優度評估顯示,所構建的SEM是很可觀的。

(潛)變數間的路徑關係也都很顯著,表明政治民主發展水平確實受工業化發展水平的推動。

路徑圖展示圖,我們可根據引數估計值的大小,並結合路徑的顯著性(檢視summary(sem_fit)中的概要,包括引數估計值,以及顯著性p值),評估關鍵的變數關係。對於示例資料,就可獲知哪種工業化或政治民主指標對於國家發展是更重要的。

SEM模型改進
#透過 MI 值評估是否需要考慮一些遺漏的變數間關係
mf <- modificationindices(sem_fit)
mf <- mf[order(mf$mi, decreasing = TRUE), ]
head(mf)

lhs、op、rhs指示了建議新增的(潛)變數路徑。MI值越高,代表這種路徑的新增更有利於改善現有模型。(新增時,將這種關係在SEM計算的第一步,加入至假定變數結構公式中)
#例如將 ind60 與 y4 的關係考慮至模型中
sem_model2 <- '
#latent variables
ind60 =~ x1 + x2 + x3 + y4 #在原式中繼續新增即可
dem60 =~ y1 + y2 + y3 + y4
dem65 =~ y5 + y6 + y7 + y8
#regressions
dem60 ~ ind60
dem65 ~ ind60 + dem60
#residual covariances
y1 ~~ y5
y2 ~~ y4 + y6
y3 ~~ y7
y4 ~~ y8
y6 ~~ y8
'
#執行 SEM
sem_fit2 <- sem(model = sem_model2, data = dat, se = 'bootstrap', bootstrap = 100)
summary(sem_fit2)
#擬合度評估
fitmeasures(sem_fit2)
不過切記這種“建議”只是基於某種數學上的統計指標給定,不可盲目接受,還需仔細考慮它們的存在意義是否是合適的。就上述示例而言,對於ind60(1960年的工業化指標)與y4(1960年當選立法機關的效力)的關係,真的是種合理的推斷嗎?此時應當謹慎對待。
參考資料
http://lavaan.ugent.be
https://benwhalley.github.io/just-enough-r/cfa.html



