Chapter 4

Dyadic MLM

Data load

sex = -1 : 여자 , 1 : 남자
contrib = 가사일과 경제적으로 기여하는 정도가 혼합된 점수로 grand mean centering됨.
future = 커플이 미래를 긍정적으로 보는 정도를 측정한 변수.
culture = 1 : 미국인, -1 : 아시아인
culture는 dyad멤버 모두에게 동일하게 적용되기 때문에 MLM에서 level 2변수로 사용가능함.
MLM에서 level 2 변수는 dyad 멤버들 모두 동일해야함.

setwd("C:/Users/LG/Documents/Dyadic Analysis Book Data")
library(haven)

## Warning: package 'haven' was built under R version 3.5.3

MLMdata <- read_spss("chapter4 table43.sav")
MLMdata

## # A tibble: 20 x 9
##     dyad future person   sex contrib culture intercep    x1    x2
##    <dbl>  <dbl>  <dbl> <dbl>   <dbl>   <dbl>    <dbl> <dbl> <dbl>
##  1     1     75      1    -1     -10       1        1     1     0
##  2     1     90      2     1      -5       1        1     0     1
##  3     2     55      1    -1       0       1        1     1     0
##  4     2     75      2     1      10       1        1     0     1
##  5     3     45      1    -1     -10       1        1     1     0
##  6     3     33      2     1     -15       1        1     0     1
##  7     4     70      1    -1       5       1        1     1     0
##  8     4     75      2     1      15       1        1     0     1
##  9     5     50      1    -1       0       1        1     1     0
## 10     5     40      2     1      -5       1        1     0     1
## 11     6     85      1    -1     -10      -1        1     1     0
## 12     6     90      2     1      20      -1        1     0     1
## 13     7     75      1    -1      -5      -1        1     1     0
## 14     7     80      2     1       0      -1        1     0     1
## 15     8     90      1    -1       5      -1        1     1     0
## 16     8     68      2     1       0      -1        1     0     1
## 17     9     65      1    -1       0      -1        1     1     0
## 18     9     78      2     1      15      -1        1     0     1
## 19    10     88      1    -1     -15      -1        1     1     0
## 20    10     95      2     1       5      -1        1     0     1

Indistinguishable Dyads MLM

heterogeneous variances는 lmer로는 불가능하고, lme 또는 gls로 가능함.
Dyadic data를 MLM으로 분석하기 위해서는 pair wise data형태로 구성되어야 함.

#lme
library(nlme)
imlmmodel <- lme(future ~ contrib + culture + contrib:culture, 
                 random = ~1|dyad, data = MLMdata, 
                 correlation = corCompSymm(form=~1|dyad))
summary(imlmmodel)

## Linear mixed-effects model fit by REML
##  Data: MLMdata 
##        AIC      BIC    logLik
##   156.9258 162.3339 -71.46289
## 
## Random effects:
##  Formula: ~1 | dyad
##         (Intercept) Residual
## StdDev:    13.28266 6.614647
## 
## Correlation Structure: Compound symmetry
##  Formula: ~1 | dyad 
##  Parameter estimate(s):
## Rho 
##   0 
## Fixed effects: future ~ contrib + culture + contrib:culture 
##                    Value Std.Error DF   t-value p-value
## (Intercept)     71.83089  4.469636  8 16.070858  0.0000
## contrib          0.84479  0.255653  8  3.304426  0.0108
## culture         -9.03282  4.469636  8 -2.020929  0.0779
## contrib:culture  0.48726  0.255653  8  1.905943  0.0931
##  Correlation: 
##                 (Intr) contrb cultur
## contrib         0.050               
## culture         0.004  0.086        
## contrib:culture 0.086  0.587  0.050 
## 
## Standardized Within-Group Residuals:
##        Min         Q1        Med         Q3        Max 
## -1.5738389 -0.5010191  0.1034914  0.3899364  1.4818595 
## 
## Number of Observations: 20
## Number of Groups: 10

#variance 확인
#sigma : an optional numeric value used as a multiplier for thestandard deviations
#rdig : an optional integer value specifying the number of digitsused to represent correlation estimates
VarCorr(imlmmodel, sigma = imlmmodel$sigma, rdig = 3)

## dyad = pdLogChol(1) 
##             Variance  StdDev   
## (Intercept) 176.42915 13.282664
## Residual     43.75355  6.614647

Distinguishable Dyads MLM

two intercept model (intercept)

Dyads member들이 indistinguishalbe일 때, MLM을 사용하는 것이 유용함.
Dyads member들이 distinguishable일 때, 사용하는 방법이 two intercept model임.

#dummy변수 만들기
#성별이라는 distinguishable 변수를 사용해서 분석하기 때문에 성별을 더미변수로 만듬. 
twointer <- cbind(MLMdata$x1, MLMdata$x2)
colnames(twointer) <- c("dummy_w", "dummy_m")
twointer

##       dummy_w dummy_m
##  [1,]       1       0
##  [2,]       0       1
##  [3,]       1       0
##  [4,]       0       1
##  [5,]       1       0
##  [6,]       0       1
##  [7,]       1       0
##  [8,]       0       1
##  [9,]       1       0
## [10,]       0       1
## [11,]       1       0
## [12,]       0       1
## [13,]       1       0
## [14,]       0       1
## [15,]       1       0
## [16,]       0       1
## [17,]       1       0
## [18,]       0       1
## [19,]       1       0
## [20,]       0       1

twointerdata <- cbind(MLMdata[1:6],twointer)
twointerdata

##    dyad future person sex contrib culture dummy_w dummy_m
## 1     1     75      1  -1     -10       1       1       0
## 2     1     90      2   1      -5       1       0       1
## 3     2     55      1  -1       0       1       1       0
## 4     2     75      2   1      10       1       0       1
## 5     3     45      1  -1     -10       1       1       0
## 6     3     33      2   1     -15       1       0       1
## 7     4     70      1  -1       5       1       1       0
## 8     4     75      2   1      15       1       0       1
## 9     5     50      1  -1       0       1       1       0
## 10    5     40      2   1      -5       1       0       1
## 11    6     85      1  -1     -10      -1       1       0
## 12    6     90      2   1      20      -1       0       1
## 13    7     75      1  -1      -5      -1       1       0
## 14    7     80      2   1       0      -1       0       1
## 15    8     90      1  -1       5      -1       1       0
## 16    8     68      2   1       0      -1       0       1
## 17    9     65      1  -1       0      -1       1       0
## 18    9     78      2   1      15      -1       0       1
## 19   10     88      1  -1     -15      -1       1       0
## 20   10     95      2   1       5      -1       0       1

#correlation 확인
#더미변수 간 상관은 항상 -1임.
cor(twointerdata$dummy_w,twointerdata$dummy_m)

## [1] -1

#lme로 추정하기
library(nlme)

dmlmmodel <- lme(future ~ dummy_w + dummy_m + contrib + culture + contrib:culture -1, 
                 random = ~ -1 + dummy_w + dummy_m|dyad, 
                 data = twointerdata, 
                 correlation = corCompSymm(form=~1|dyad), 
                 weight = varIdent(form = ~1|sex))
summary(dmlmmodel)

## Linear mixed-effects model fit by REML
##  Data: twointerdata 
##        AIC      BIC    logLik
##   159.4705 167.2591 -68.73527
## 
## Random effects:
##  Formula: ~-1 + dummy_w + dummy_m | dyad
##  Structure: General positive-definite, Log-Cholesky parametrization
##          StdDev    Corr  
## dummy_w  13.160079 dmmy_w
## dummy_m  15.702344 0.872 
## Residual  4.086932       
## 
## Correlation Structure: Compound symmetry
##  Formula: ~1 | dyad 
##  Parameter estimate(s):
## Rho 
##   0 
## Variance function:
##  Structure: Different standard deviations per stratum
##  Formula: ~1 | sex 
##  Parameter estimates:
## -1  1 
##  1  1 
## Fixed effects: future ~ dummy_w + dummy_m + contrib + culture + contrib:culture -      1 
##                    Value Std.Error DF   t-value p-value
## dummy_w         72.88184  4.497691  7 16.204280  0.0000
## dummy_m         70.69451  5.248280  7 13.470035  0.0000
## contrib          0.88516  0.311349  7  2.842971  0.0249
## culture         -9.66655  4.418493  9 -2.187748  0.0565
## contrib:culture  0.45878  0.275408  7  1.665835  0.1397
##  Correlation: 
##                 dmmy_w dmmy_m contrb cultur
## dummy_m          0.726                     
## contrib          0.242 -0.119              
## culture          0.011  0.019  0.078       
## contrib:culture  0.095  0.076  0.564  0.178
## 
## Standardized Within-Group Residuals:
##         Min          Q1         Med          Q3         Max 
## -0.77109796 -0.24395895  0.04312602  0.28546209  0.84622291 
## 
## Number of Observations: 20
## Number of Groups: 10

#더미변수가 random variable이기 때문에 random 함수에 sex를 포함함. random 함수의 group 변수와 correlation의 group 변수가 동일해야함. 

# corCompSymm : representing a compound symmetry structure corresponding to uniform correlation
# -> Defaults to ~ 1, which corresponds to using the order of the observations in the data as a covariate, and no groups, |이후 dyad는 grouping factor가 들어가는 자리임
# varIdent : representing a constant variance function structure
# varIdent에 form에는 v|g 형태이며, v = variance covariate, g = grouping factor임 
# dyad간 상관이 있으니 corCompSymm을 사용하고, 성별간 분산을 다르게 추정할 수 있게 하기 위해서 varIdent를 사용함.

#variance 확인
VarCorr(dmlmmodel, sigma = dmlmmodel$sigma, rdig = 3)

## dyad = pdLogChol(-1 + dummy_w + dummy_m) 
##          Variance  StdDev    Corr  
## dummy_w  173.18767 13.160079 dmmy_w
## dummy_m  246.56360 15.702344 0.872 
## Residual  16.70301  4.086932

#dummy_w, dummy_m의 covariate
#correltaion = covariate/standard deviation
#covariate = correlation * stadnard deviation

#cov(dummy_w, dummy_m) = cor(dummy_w, dummy_m) * sd(dummy_w) * sd(dummy_m)
13.160079*15.702344*0.872

## [1] 180.1936

Two intercept model (predictor)

library(nlme)

pmlmmodel <- lme(future ~ dummy_w + dummy_m + contrib:dummy_w + contrib:dummy_m -1, 
                 random = ~ -1 + dummy_w + dummy_m|dyad, 
                 data = twointerdata, 
                 correlation = corCompSymm(form=~1|dyad), 
                 weight = varIdent(form = ~1|sex))
summary(pmlmmodel)

## Linear mixed-effects model fit by REML
##  Data: twointerdata 
##        AIC      BIC    logLik
##   167.9316 175.6575 -73.96582
## 
## Random effects:
##  Formula: ~-1 + dummy_w + dummy_m | dyad
##  Structure: General positive-definite, Log-Cholesky parametrization
##          StdDev    Corr  
## dummy_w  16.420695 dmmy_w
## dummy_m  16.430275 0.826 
## Residual  4.849573       
## 
## Correlation Structure: Compound symmetry
##  Formula: ~1 | dyad 
##  Parameter estimate(s):
## Rho 
##   0 
## Variance function:
##  Structure: Different standard deviations per stratum
##  Formula: ~1 | sex 
##  Parameter estimates:
## -1  1 
##  1  1 
## Fixed effects: future ~ dummy_w + dummy_m + contrib:dummy_w + contrib:dummy_m -      1 
##                    Value Std.Error DF   t-value p-value
## dummy_w         71.65751  5.827890  7 12.295617  0.0000
## dummy_m         69.23726  5.587998  7 12.390352  0.0000
## dummy_w:contrib  0.46438  0.539012  7  0.861534  0.4175
## dummy_m:contrib  0.79068  0.342645  7  2.307587  0.0544
##  Correlation: 
##                 dmmy_w dmmy_m dmmy_w:
## dummy_m          0.668               
## dummy_w:contrib  0.370 -0.043        
## dummy_m:contrib  0.065 -0.245  0.176 
## 
## Standardized Within-Group Residuals:
##         Min          Q1         Med          Q3         Max 
## -0.58827968 -0.30936029 -0.09717406  0.25534017  0.72901009 
## 
## Number of Observations: 20
## Number of Groups: 10

VarCorr(pmlmmodel, sigma = pmlmmodel$sigma, rdig = 3)

## dyad = pdLogChol(-1 + dummy_w + dummy_m) 
##          Variance  StdDev    Corr  
## dummy_w  269.63922 16.420695 dmmy_w
## dummy_m  269.95394 16.430275 0.826 
## Residual  23.51836  4.849573

#cov
0.826*16.420695*16.430275

## [1] 222.8519

Test Distinguishablity MLM

Model 2개의 fixed effect가 다르기 때문에 추정방법 = ML인 Model1, Model2를 추정후 ANOVA로 비교한다.
Model 1 = Dyads가 indistinguishable인 경우
Model 2 = Dyads가 distinguishable인 경우
ANOVA결과 Model이 유의한 차이가 없으면 Dyads를 indistinguishable로 다룸

library(nlme)

model1 <- lme(future ~ contrib + culture + contrib:culture, random = ~1|dyad, data = MLMdata, correlation = corCompSymm(form=~1|dyad), method = "ML")

summary(model1)

## Linear mixed-effects model fit by maximum likelihood
##  Data: MLMdata 
##        AIC     BIC    logLik
##   163.8859 170.856 -74.94294
## 
## Random effects:
##  Formula: ~1 | dyad
##         (Intercept) Residual
## StdDev:    11.49586 6.065927
## 
## Correlation Structure: Compound symmetry
##  Formula: ~1 | dyad 
##  Parameter estimate(s):
## Rho 
##   0 
## Fixed effects: future ~ contrib + culture + contrib:culture 
##                    Value Std.Error DF   t-value p-value
## (Intercept)     71.82103  4.355355  8 16.490281  0.0000
## contrib          0.83445  0.258247  8  3.231206  0.0120
## culture         -9.04833  4.355355  8 -2.077517  0.0714
## contrib:culture  0.48069  0.258247  8  1.861350  0.0997
##  Correlation: 
##                 (Intr) contrb cultur
## contrib         0.051               
## culture         0.005  0.089        
## contrib:culture 0.089  0.576  0.051 
## 
## Standardized Within-Group Residuals:
##        Min         Q1        Med         Q3        Max 
## -1.7230855 -0.5792601  0.1033656  0.4132005  1.6121334 
## 
## Number of Observations: 20
## Number of Groups: 10

model2 <- lme(future ~ contrib + culture + contrib:culture + sex + contrib:sex + culture:sex + contrib:culture:sex, random = ~1|dyad, data = MLMdata, correlation = corCompSymm(form=~1|dyad), 
              weight = varIdent(form = ~1|sex), method = "ML")

summary(model2)

## Linear mixed-effects model fit by maximum likelihood
##  Data: MLMdata 
##        AIC      BIC    logLik
##   171.8413 183.7901 -73.92063
## 
## Random effects:
##  Formula: ~1 | dyad
##         (Intercept) Residual
## StdDev:    11.56877 2.985059
## 
## Correlation Structure: Compound symmetry
##  Formula: ~1 | dyad 
##  Parameter estimate(s):
##      Rho 
## 0.139964 
## Variance function:
##  Structure: Different standard deviations per stratum
##  Formula: ~1 | sex 
##  Parameter estimates:
##       -1        1 
## 1.000000 2.618752 
## Fixed effects: future ~ contrib + culture + contrib:culture + sex + contrib:sex +      culture:sex + contrib:culture:sex 
##                        Value Std.Error DF   t-value p-value
## (Intercept)         71.22949  5.305882  8 13.424626  0.0000
## contrib              0.92193  0.385782  4  2.389778  0.0752
## culture             -8.55800  5.305882  8 -1.612927  0.1454
## sex                 -1.59781  2.143668  4 -0.745363  0.4975
## contrib:culture      0.34473  0.385782  4  0.893590  0.4220
## contrib:sex          0.06694  0.277351  4  0.241360  0.8211
## culture:sex          1.52632  2.143668  4  0.712015  0.5158
## contrib:culture:sex -0.04794  0.277351  4 -0.172837  0.8712
##  Correlation: 
##                     (Intr) contrb cultur sex    cntrb:c cntrb:s cltr:s
## contrib              0.081                                            
## culture             -0.062  0.133                                     
## sex                  0.224 -0.521 -0.114                              
## contrib:culture      0.133  0.565  0.081 -0.010                       
## contrib:sex         -0.280 -0.274  0.164  0.010 -0.287                
## culture:sex         -0.114 -0.010  0.224 -0.225 -0.521   0.277        
## contrib:culture:sex  0.164 -0.287 -0.280  0.277 -0.274  -0.300   0.010
## 
## Standardized Within-Group Residuals:
##         Min          Q1         Med          Q3         Max 
## -1.53936408 -0.52451330  0.09766917  0.43683003  1.25314878 
## 
## Number of Observations: 20
## Number of Groups: 10

anova(model1, model2)

##        Model df      AIC      BIC    logLik   Test L.Ratio p-value
## model1     1  7 163.8859 170.8560 -74.94294                       
## model2     2 12 171.8413 183.7901 -73.92063 1 vs 2 2.04461  0.8429

Chapter 5

Dyadic CFA

Data load

Table 5.1 데이터를 입력하고 분석함.
남자 일란성 쌍둥이 데이터임.
Table 5.1 은 상관 테이블이며 n = 137, 평균, 표준편차가 나타나 있음.

A <- matrix(c(1,0.380,0.351,0.531,0.411,0.313,0.118,0.214,0,
              1,0.483,0.386,0.161,0.453,0.080,0.148,0,0,
              1,0.385,0.142,0.266,0.092,0.129,0,0,0,
              1,0.228,0.245,0.099,0.357,0,0,0,0,
              1,0.348,0.323,0.403,0,0,0,0,0,
              1,0.381,0.431,0,0,0,0,0,0,
              1,0.256,0,0,0,0,0,0,0,
              1), 8,8)
A

##       [,1]  [,2]  [,3]  [,4]  [,5]  [,6]  [,7] [,8]
## [1,] 1.000 0.000 0.000 0.000 0.000 0.000 0.000    0
## [2,] 0.380 1.000 0.000 0.000 0.000 0.000 0.000    0
## [3,] 0.351 0.483 1.000 0.000 0.000 0.000 0.000    0
## [4,] 0.531 0.386 0.385 1.000 0.000 0.000 0.000    0
## [5,] 0.411 0.161 0.142 0.228 1.000 0.000 0.000    0
## [6,] 0.313 0.453 0.266 0.245 0.348 1.000 0.000    0
## [7,] 0.118 0.080 0.092 0.099 0.323 0.381 1.000    0
## [8,] 0.214 0.148 0.129 0.357 0.403 0.431 0.256    1

library(gdata) #상관행렬 upper 또는 lower 채우는 packages

A[upper.tri(A)]<-t(lowerTriangle(A, diag=F, byrow=TRUE))
data <- round(A,3)
data

##       [,1]  [,2]  [,3]  [,4]  [,5]  [,6]  [,7]  [,8]
## [1,] 1.000 0.380 0.351 0.531 0.411 0.313 0.118 0.214
## [2,] 0.380 1.000 0.483 0.386 0.161 0.453 0.080 0.148
## [3,] 0.351 0.483 1.000 0.385 0.142 0.266 0.092 0.129
## [4,] 0.531 0.386 0.385 1.000 0.228 0.245 0.099 0.357
## [5,] 0.411 0.161 0.142 0.228 1.000 0.348 0.323 0.403
## [6,] 0.313 0.453 0.266 0.245 0.348 1.000 0.381 0.431
## [7,] 0.118 0.080 0.092 0.099 0.323 0.381 1.000 0.256
## [8,] 0.214 0.148 0.129 0.357 0.403 0.431 0.256 1.000

colnames(data) <- c("x1","x2","x3","x4","x5","x6","x7","x8")
rownames(data) <- colnames(data)
data

##       x1    x2    x3    x4    x5    x6    x7    x8
## x1 1.000 0.380 0.351 0.531 0.411 0.313 0.118 0.214
## x2 0.380 1.000 0.483 0.386 0.161 0.453 0.080 0.148
## x3 0.351 0.483 1.000 0.385 0.142 0.266 0.092 0.129
## x4 0.531 0.386 0.385 1.000 0.228 0.245 0.099 0.357
## x5 0.411 0.161 0.142 0.228 1.000 0.348 0.323 0.403
## x6 0.313 0.453 0.266 0.245 0.348 1.000 0.381 0.431
## x7 0.118 0.080 0.092 0.099 0.323 0.381 1.000 0.256
## x8 0.214 0.148 0.129 0.357 0.403 0.431 0.256 1.000

#x1 : 멤버 1 재평가, x2 : 멤버 1 예지력, x3 : 멤버 1 통찰력, x4 : 멤버 1 자기주도성
#x5 : 멤버 2 재평가, x6 : 멤버 2 예지력, x7 : 멤버 2 통찰력, x8 : 멤버 2 자기주도성

mean <- c(3.042,2.571,2.903,3.095,3.074,2.474,2.913,3.144)
mean

## [1] 3.042 2.571 2.903 3.095 3.074 2.474 2.913 3.144

sd <- c(0.579,0.626,0.607,0.653,0.662,0.605,0.615,0.694)
sd

## [1] 0.579 0.626 0.607 0.653 0.662 0.605 0.615 0.694

library(lavaan) #cor2cov는 lavaan에 있음.

covdata <- cor2cov(data, sd)
cordata <- cov2cor(covdata)

round(covdata, 3)

##       x1    x2    x3    x4    x5    x6    x7    x8
## x1 0.335 0.138 0.123 0.201 0.158 0.110 0.042 0.086
## x2 0.138 0.392 0.184 0.158 0.067 0.172 0.031 0.064
## x3 0.123 0.184 0.368 0.153 0.057 0.098 0.034 0.054
## x4 0.201 0.158 0.153 0.426 0.099 0.097 0.040 0.162
## x5 0.158 0.067 0.057 0.099 0.438 0.139 0.132 0.185
## x6 0.110 0.172 0.098 0.097 0.139 0.366 0.142 0.181
## x7 0.042 0.031 0.034 0.040 0.132 0.142 0.378 0.109
## x8 0.086 0.064 0.054 0.162 0.185 0.181 0.109 0.482

round(cordata, 3)

##       x1    x2    x3    x4    x5    x6    x7    x8
## x1 1.000 0.380 0.351 0.531 0.411 0.313 0.118 0.214
## x2 0.380 1.000 0.483 0.386 0.161 0.453 0.080 0.148
## x3 0.351 0.483 1.000 0.385 0.142 0.266 0.092 0.129
## x4 0.531 0.386 0.385 1.000 0.228 0.245 0.099 0.357
## x5 0.411 0.161 0.142 0.228 1.000 0.348 0.323 0.403
## x6 0.313 0.453 0.266 0.245 0.348 1.000 0.381 0.431
## x7 0.118 0.080 0.092 0.099 0.323 0.381 1.000 0.256
## x8 0.214 0.148 0.129 0.357 0.403 0.431 0.256 1.000

Indistinguishable Dyads CFA

Null Model

Null Model은 평균, 분산만 추정하는 Model임.
Dyadic Null Model은 covariate을 모두 0으로 고정함.

library(lavaan)

nullmodel <- '
x1 ~ m1*1 # Means 
x5 ~ m1*1
x2 ~ m2*1
x6 ~ m2*1
x3 ~ m3*1
x7 ~ m3*1
x4 ~ m4*1
x8 ~ m4*1

x1 ~~ r1*x1 #Variance
x2 ~~ r2*x2
x3 ~~ r3*x3
x4 ~~ r4*x4
x5 ~~ r1*x5
x6 ~~ r2*x6
x7 ~~ r3*x7
x8 ~~ r4*x8

x1 ~~ 0*x2 #Covariate
x1 ~~ 0*x3 
x1 ~~ 0*x4
x1 ~~ 0*x5
x1 ~~ 0*x6
x1 ~~ 0*x7
x1 ~~ 0*x8
x2 ~~ 0*x3
x2 ~~ 0*x4
x2 ~~ 0*x5
x2 ~~ 0*x6
x2 ~~ 0*x7
x2 ~~ 0*x8
x3 ~~ 0*x4
x3 ~~ 0*x5
x3 ~~ 0*x6 
x3 ~~ 0*x7 
x3 ~~ 0*x8
x4 ~~ 0*x5
x4 ~~ 0*x6
x4 ~~ 0*x7
x4 ~~ 0*x8
x5 ~~ 0*x6
x5 ~~ 0*x7
x5 ~~ 0*x8
x6 ~~ 0*x7 
x6 ~~ 0*x8
x7 ~~ 0*x8
'
nullfit <- cfa(nullmodel, sample.cov = covdata, sample.mean = mean, sample.nobs = 137, 
               std.lv = F, mimic = "Mplus") # std.lv = TRUE 는 요인 분산 1로 고정하는 명령어.
summary(nullfit, fit.measures = TRUE)

## lavaan 0.6-5 ended normally after 23 iterations
## 
##   Estimator                                         ML
##   Optimization method                           NLMINB
##   Number of free parameters                         16
##   Number of equality constraints                     8
##   Row rank of the constraints matrix                 8
##                                                       
##   Number of observations                           137
##                                                       
## Model Test User Model:
##                                                       
##   Test statistic                               285.154
##   Degrees of freedom                                36
##   P-value (Chi-square)                           0.000
## 
## Model Test Baseline Model:
## 
##   Test statistic                               279.739
##   Degrees of freedom                                28
##   P-value                                        0.000
## 
## User Model versus Baseline Model:
## 
##   Comparative Fit Index (CFI)                    0.010
##   Tucker-Lewis Index (TLI)                       0.230
## 
## Loglikelihood and Information Criteria:
## 
##   Loglikelihood user model (H0)              -1046.033
##   Loglikelihood unrestricted model (H1)       -903.455
##                                                       
##   Akaike (AIC)                                2108.065
##   Bayesian (BIC)                              2131.425
##   Sample-size adjusted Bayesian (BIC)         2106.117
## 
## Root Mean Square Error of Approximation:
## 
##   RMSEA                                          0.225
##   90 Percent confidence interval - lower         0.201
##   90 Percent confidence interval - upper         0.249
##   P-value RMSEA <= 0.05                          0.000
## 
## Standardized Root Mean Square Residual:
## 
##   SRMR                                           0.264
## 
## Parameter Estimates:
## 
##   Information                                 Observed
##   Observed information based on                Hessian
##   Standard errors                             Standard
## 
## Covariances:
##                    Estimate  Std.Err  z-value  P(>|z|)
##   x1 ~~                                               
##     x2                0.000                           
##     x3                0.000                           
##     x4                0.000                           
##     x5                0.000                           
##     x6                0.000                           
##     x7                0.000                           
##     x8                0.000                           
##   x2 ~~                                               
##     x3                0.000                           
##     x4                0.000                           
##     x5                0.000                           
##     x6                0.000                           
##     x7                0.000                           
##     x8                0.000                           
##   x3 ~~                                               
##     x4                0.000                           
##     x5                0.000                           
##     x6                0.000                           
##     x7                0.000                           
##     x8                0.000                           
##   x4 ~~                                               
##     x5                0.000                           
##     x6                0.000                           
##     x7                0.000                           
##     x8                0.000                           
##   x5 ~~                                               
##     x6                0.000                           
##     x7                0.000                           
##     x8                0.000                           
##   x6 ~~                                               
##     x7                0.000                           
##     x8                0.000                           
##   x7 ~~                                               
##     x8                0.000                           
## 
## Intercepts:
##                    Estimate  Std.Err  z-value  P(>|z|)
##     x1        (m1)    3.058    0.037   81.667    0.000
##     x5        (m1)    3.058    0.037   81.667    0.000
##     x2        (m2)    2.522    0.037   67.866    0.000
##     x6        (m2)    2.522    0.037   67.866    0.000
##     x3        (m3)    2.908    0.037   79.067    0.000
##     x7        (m3)    2.908    0.037   79.067    0.000
##     x4        (m4)    3.120    0.041   76.864    0.000
##     x8        (m4)    3.120    0.041   76.864    0.000
## 
## Variances:
##                    Estimate  Std.Err  z-value  P(>|z|)
##     x1        (r1)    0.384    0.033   11.705    0.000
##     x2        (r2)    0.379    0.032   11.705    0.000
##     x3        (r3)    0.371    0.032   11.705    0.000
##     x4        (r4)    0.451    0.039   11.705    0.000
##     x5        (r1)    0.384    0.033   11.705    0.000
##     x6        (r2)    0.379    0.032   11.705    0.000
##     x7        (r3)    0.371    0.032   11.705    0.000
##     x8        (r4)    0.451    0.039   11.705    0.000

I-SAT Model

I-SAT : saturated model for indistinguishable dyads
I-SAT Model은 분산-공분산 행렬을 input으로 사용함.
잠재변수가 Model에 포함되지 않음.
평균, 분산, 공분산에 제약을 가함.
intrapersonal 공분산, interpersonal 공분산을 동일하게 제약함.
intrapersonal : 멤버내에 측정된 변수들 간 상관이 서로 동일하게 제약된 것.
interpersonal : 멤버간 측정 변수들 간 상관이 서로 동일하게 제약된 것.
일반적으로 측정된 변수가 k개 라면, I-SAT Model은 총 k(k+1)개의 제약을 가짐
- k개 평균.
- k개 분산.
- k(k-1)개 공분산(intrapersonal 공분산 + interpersonal 공분산).
- k(k+1)의 자유도.
I-SAT의 chi-square는 해석하지 않고 수정된 chi-square, 수정된 df를 위해서 사용함.

library(lavaan)

isatmodel <- ' 
x1 ~ m1*1 #Means
x5 ~ m1*1
x2 ~ m2*1
x6 ~ m2*1
x3 ~ m3*1
x7 ~ m3*1
x4 ~ m4*1
x8 ~ m4*1

x1 ~~ v1*x1 #Variance
x2 ~~ v2*x2
x3 ~~ v3*x3
x4 ~~ v4*x4
x5 ~~ v1*x5
x6 ~~ v2*x6
x7 ~~ v3*x7
x8 ~~ v4*x8

x1 ~~ c12*x2 #Within
x1 ~~ c13*x3
x1 ~~ c14*x4
x2 ~~ c23*x3
x2 ~~ c24*x4
x3 ~~ c34*x4

x5 ~~ c12*x6
x5 ~~ c13*x7
x5 ~~ c14*x8
x6 ~~ c23*x7
x6 ~~ c24*x8
x7 ~~ c34*x8

x1 ~~ x5 #residual correlation
x2 ~~ x6
x3 ~~ x7
x4 ~~ x8

x1 ~~ cc12*x6 #Between
x1 ~~ cc13*x7 
x1 ~~ cc14*x8
x2 ~~ cc23*x7
x2 ~~ cc24*x8
x3 ~~ cc34*x8

x5 ~~ cc12*x2
x5 ~~ cc13*x3
x5 ~~ cc14*x4
x6 ~~ cc23*x3
x6 ~~ cc24*x4
x7 ~~ cc34*x4
'  
isatfit <- cfa(isatmodel, sample.cov = covdata, sample.mean = mean, sample.nobs = 137, 
               std.lv = F, mimic = "Mplus")
summary(isatfit, fit.measures = TRUE)

## lavaan 0.6-5 ended normally after 59 iterations
## 
##   Estimator                                         ML
##   Optimization method                           NLMINB
##   Number of free parameters                         44
##   Number of equality constraints                    20
##   Row rank of the constraints matrix                20
##                                                       
##   Number of observations                           137
##                                                       
## Model Test User Model:
##                                                       
##   Test statistic                                21.839
##   Degrees of freedom                                20
##   P-value (Chi-square)                           0.349
## 
## Model Test Baseline Model:
## 
##   Test statistic                               279.739
##   Degrees of freedom                                28
##   P-value                                        0.000
## 
## User Model versus Baseline Model:
## 
##   Comparative Fit Index (CFI)                    0.993
##   Tucker-Lewis Index (TLI)                       0.990
## 
## Loglikelihood and Information Criteria:
## 
##   Loglikelihood user model (H0)               -914.375
##   Loglikelihood unrestricted model (H1)       -903.455
##                                                       
##   Akaike (AIC)                                1876.750
##   Bayesian (BIC)                              1946.829
##   Sample-size adjusted Bayesian (BIC)         1870.904
## 
## Root Mean Square Error of Approximation:
## 
##   RMSEA                                          0.026
##   90 Percent confidence interval - lower         0.000
##   90 Percent confidence interval - upper         0.080
##   P-value RMSEA <= 0.05                          0.709
## 
## Standardized Root Mean Square Residual:
## 
##   SRMR                                           0.085
## 
## Parameter Estimates:
## 
##   Information                                 Observed
##   Observed information based on                Hessian
##   Standard errors                             Standard
## 
## Covariances:
##                    Estimate  Std.Err  z-value  P(>|z|)
##   x1 ~~                                               
##     x2       (c12)    0.137    0.027    5.086    0.000
##     x3       (c13)    0.127    0.025    5.140    0.000
##     x4       (c14)    0.192    0.030    6.438    0.000
##   x2 ~~                                               
##     x3       (c23)    0.161    0.025    6.359    0.000
##     x4       (c24)    0.167    0.029    5.735    0.000
##   x3 ~~                                               
##     x4       (c34)    0.130    0.026    4.917    0.000
##   x5 ~~                                               
##     x6       (c12)    0.137    0.027    5.086    0.000
##     x7       (c13)    0.127    0.025    5.140    0.000
##     x8       (c14)    0.192    0.030    6.438    0.000
##   x6 ~~                                               
##     x7       (c23)    0.161    0.025    6.359    0.000
##     x8       (c24)    0.167    0.029    5.735    0.000
##   x7 ~~                                               
##     x8       (c34)    0.130    0.026    4.917    0.000
##   x1 ~~                                               
##     x5                0.156    0.035    4.407    0.000
##   x2 ~~                                               
##     x6                0.168    0.035    4.747    0.000
##   x3 ~~                                               
##     x7                0.034    0.032    1.071    0.284
##   x4 ~~                                               
##     x8                0.160    0.041    3.911    0.000
##   x1 ~~                                               
##     x6      (cc12)    0.088    0.027    3.284    0.001
##     x7      (cc13)    0.049    0.025    1.994    0.046
##     x8      (cc14)    0.091    0.030    3.059    0.002
##   x2 ~~                                               
##     x7      (cc23)    0.064    0.025    2.525    0.012
##     x8      (cc24)    0.081    0.029    2.788    0.005
##   x3 ~~                                               
##     x8      (cc34)    0.047    0.026    1.761    0.078
##   x2 ~~                                               
##     x5      (cc12)    0.088    0.027    3.284    0.001
##   x3 ~~                                               
##     x5      (cc13)    0.049    0.025    1.994    0.046
##   x4 ~~                                               
##     x5      (cc14)    0.091    0.030    3.059    0.002
##   x3 ~~                                               
##     x6      (cc23)    0.064    0.025    2.525    0.012
##   x4 ~~                                               
##     x6      (cc24)    0.081    0.029    2.788    0.005
##     x7      (cc34)    0.047    0.026    1.761    0.078
## 
## Intercepts:
##                    Estimate  Std.Err  z-value  P(>|z|)
##     x1        (m1)    3.058    0.044   68.864    0.000
##     x5        (m1)    3.058    0.044   68.864    0.000
##     x2        (m2)    2.523    0.045   56.482    0.000
##     x6        (m2)    2.523    0.045   56.482    0.000
##     x3        (m3)    2.908    0.038   75.666    0.000
##     x7        (m3)    2.908    0.038   75.666    0.000
##     x4        (m4)    3.120    0.047   66.043    0.000
##     x8        (m4)    3.120    0.047   66.043    0.000
## 
## Variances:
##                    Estimate  Std.Err  z-value  P(>|z|)
##     x1        (v1)    0.384    0.035   10.843    0.000
##     x2        (v2)    0.379    0.035   10.699    0.000
##     x3        (v3)    0.371    0.032   11.656    0.000
##     x4        (v4)    0.451    0.041   11.032    0.000
##     x5        (v1)    0.384    0.035   10.843    0.000
##     x6        (v2)    0.379    0.035   10.699    0.000
##     x7        (v3)    0.371    0.032   11.656    0.000
##     x8        (v4)    0.451    0.041   11.032    0.000

Indistinguishable dyads CFA Model

NA*는 첫번째 indicator를 free로 추정하는 방법.
marker indicator가 4번째 indicator이고 요인분산까지 추정하는 모델임.
free로 추정하는 첫번째 indicator의 factor loading을 고정하고 싶다면, 다시 indicator변수를 입력하고 factor loadging을 고정하면 됨.
indistinguishable CFA는 figure 5.3에 나타난 그림과 동일하게 제약을 가함.
- 요인분산, 요인부하량, 절편, 오차분산을 멤버간 동일하게 제약
- 잠재변수 평균, 오차 = 0
- 멤버간 동일한 indicator의 오차공분산을 추정함.
- k = 4인 모델에서 총 17개의 unique parameter를 가짐.

library(lavaan)

icfamodel <-'
f1 =~ NA*x1 + a*x1 + b*x2 + c*x3 + 1*x4  #factor loading
f2 =~ NA*x5 + a*x5 + b*x6 + c*x7 + 1*x8


f1 ~ 0*1 #factor means
f2 ~ 0*1

f1 ~~ i*f2 #factor covariate

f1 ~~ h*f1 #factor variance
f2 ~~ h*f2


x1 ~ m1*1 #intercept
x2 ~ m2*1
x3 ~ m3*1
x4 ~ m4*1
x5 ~ m1*1
x6 ~ m2*1
x7 ~ m3*1
x8 ~ m4*1

x1 ~~ v1*x1 #residual variance
x2 ~~ v2*x2
x3 ~~ v3*x3
x4 ~~ v4*x4
x5 ~~ v1*x5
x6 ~~ v2*x6
x7 ~~ v3*x7
x8 ~~ v4*x8

x1 ~~ c1*x5 #residual correlation
x2 ~~ c2*x6
x3 ~~ c3*x7
x4 ~~ c4*x8
'
icfafit <- cfa(icfamodel, sample.cov = covdata, sample.mean = mean, sample.nobs = 137, 
               std.lv = F, mimic = "Mplus")
summary(icfafit, fit.measures = TRUE)

## lavaan 0.6-5 ended normally after 39 iterations
## 
##   Estimator                                         ML
##   Optimization method                           NLMINB
##   Number of free parameters                         29
##   Number of equality constraints                    12
##   Row rank of the constraints matrix                12
##                                                       
##   Number of observations                           137
##                                                       
## Model Test User Model:
##                                                       
##   Test statistic                                34.588
##   Degrees of freedom                                27
##   P-value (Chi-square)                           0.150
## 
## Model Test Baseline Model:
## 
##   Test statistic                               279.739
##   Degrees of freedom                                28
##   P-value                                        0.000
## 
## User Model versus Baseline Model:
## 
##   Comparative Fit Index (CFI)                    0.970
##   Tucker-Lewis Index (TLI)                       0.969
## 
## Loglikelihood and Information Criteria:
## 
##   Loglikelihood user model (H0)               -920.749
##   Loglikelihood unrestricted model (H1)       -903.455
##                                                       
##   Akaike (AIC)                                1875.499
##   Bayesian (BIC)                              1925.138
##   Sample-size adjusted Bayesian (BIC)         1871.358
## 
## Root Mean Square Error of Approximation:
## 
##   RMSEA                                          0.045
##   90 Percent confidence interval - lower         0.000
##   90 Percent confidence interval - upper         0.085
##   P-value RMSEA <= 0.05                          0.537
## 
## Standardized Root Mean Square Residual:
## 
##   SRMR                                           0.091
## 
## Parameter Estimates:
## 
##   Information                                 Observed
##   Observed information based on                Hessian
##   Standard errors                             Standard
## 
## Latent Variables:
##                    Estimate  Std.Err  z-value  P(>|z|)
##   f1 =~                                               
##     x1         (a)    0.864    0.121    7.154    0.000
##     x2         (b)    0.880    0.137    6.414    0.000
##     x3         (c)    0.775    0.127    6.084    0.000
##     x4                1.000                           
##   f2 =~                                               
##     x5         (a)    0.864    0.121    7.154    0.000
##     x6         (b)    0.880    0.137    6.414    0.000
##     x7         (c)    0.775    0.127    6.084    0.000
##     x8                1.000                           
## 
## Covariances:
##                    Estimate  Std.Err  z-value  P(>|z|)
##   f1 ~~                                               
##     f2         (i)    0.093    0.030    3.084    0.002
##  .x1 ~~                                               
##    .x5        (c1)    0.076    0.027    2.875    0.004
##  .x2 ~~                                               
##    .x6        (c2)    0.082    0.026    3.130    0.002
##  .x3 ~~                                               
##    .x7        (c3)   -0.005    0.025   -0.195    0.845
##  .x4 ~~                                               
##    .x8        (c4)    0.078    0.030    2.572    0.010
## 
## Intercepts:
##                    Estimate  Std.Err  z-value  P(>|z|)
##     f1                0.000                           
##     f2                0.000                           
##    .x1        (m1)    3.058    0.044   69.765    0.000
##    .x2        (m2)    2.523    0.044   57.533    0.000
##    .x3        (m3)    2.908    0.039   74.175    0.000
##    .x4        (m4)    3.120    0.048   65.264    0.000
##    .x5        (m1)    3.058    0.044   69.765    0.000
##    .x6        (m2)    2.523    0.044   57.533    0.000
##    .x7        (m3)    2.908    0.039   74.175    0.000
##    .x8        (m4)    3.120    0.048   65.264    0.000
## 
## Variances:
##                    Estimate  Std.Err  z-value  P(>|z|)
##     f1         (h)    0.194    0.042    4.590    0.000
##     f2         (h)    0.194    0.042    4.590    0.000
##    .x1        (v1)    0.236    0.029    8.081    0.000
##    .x2        (v2)    0.223    0.029    7.627    0.000
##    .x3        (v3)    0.254    0.028    9.209    0.000
##    .x4        (v4)    0.261    0.034    7.593    0.000
##    .x5        (v1)    0.236    0.029    8.081    0.000
##    .x6        (v2)    0.223    0.029    7.627    0.000
##    .x7        (v3)    0.254    0.028    9.209    0.000
##    .x8        (v4)    0.261    0.034    7.593    0.000

standardizedSolution(icfafit)

##    lhs op rhs est.std    se      z pvalue ci.lower ci.upper
## 1   f1 =~  x1   0.617 0.059 10.497  0.000    0.502    0.732
## 2   f1 =~  x2   0.635 0.059 10.686  0.000    0.518    0.751
## 3   f1 =~  x3   0.561 0.058  9.628  0.000    0.447    0.675
## 4   f1 =~  x4   0.653 0.056 11.675  0.000    0.544    0.763
## 5   f2 =~  x5   0.617 0.059 10.497  0.000    0.502    0.732
## 6   f2 =~  x6   0.635 0.059 10.686  0.000    0.518    0.751
## 7   f2 =~  x7   0.561 0.058  9.628  0.000    0.447    0.675
## 8   f2 =~  x8   0.653 0.056 11.675  0.000    0.544    0.763
## 9   f1 ~1       0.000 0.000     NA     NA    0.000    0.000
## 10  f2 ~1       0.000 0.000     NA     NA    0.000    0.000
## 11  f1 ~~  f2   0.478 0.099  4.838  0.000    0.285    0.672
## 12  f1 ~~  f1   1.000 0.000     NA     NA    1.000    1.000
## 13  f2 ~~  f2   1.000 0.000     NA     NA    1.000    1.000
## 14  x1 ~1       4.956 0.235 21.126  0.000    4.496    5.415
## 15  x2 ~1       4.129 0.201 20.587  0.000    3.736    4.522
## 16  x3 ~1       4.779 0.216 22.123  0.000    4.355    5.202
## 17  x4 ~1       4.626 0.223 20.762  0.000    4.189    5.063
## 18  x5 ~1       4.956 0.235 21.126  0.000    4.496    5.415
## 19  x6 ~1       4.129 0.201 20.587  0.000    3.736    4.522
## 20  x7 ~1       4.779 0.216 22.123  0.000    4.355    5.202
## 21  x8 ~1       4.626 0.223 20.762  0.000    4.189    5.063
## 22  x1 ~~  x1   0.619 0.073  8.540  0.000    0.477    0.761
## 23  x2 ~~  x2   0.597 0.075  7.923  0.000    0.449    0.745
## 24  x3 ~~  x3   0.685 0.065 10.484  0.000    0.557    0.813
## 25  x4 ~~  x4   0.573 0.073  7.836  0.000    0.430    0.716
## 26  x5 ~~  x5   0.619 0.073  8.540  0.000    0.477    0.761
## 27  x6 ~~  x6   0.597 0.075  7.923  0.000    0.449    0.745
## 28  x7 ~~  x7   0.685 0.065 10.484  0.000    0.557    0.813
## 29  x8 ~~  x8   0.573 0.073  7.836  0.000    0.430    0.716
## 30  x1 ~~  x5   0.324 0.094  3.456  0.001    0.140    0.507
## 31  x2 ~~  x6   0.366 0.093  3.943  0.000    0.184    0.548
## 32  x3 ~~  x7  -0.019 0.100 -0.195  0.846   -0.215    0.176
## 33  x4 ~~  x8   0.301 0.100  3.021  0.003    0.106    0.496

anova(nullfit, isatfit)

## Chi-Squared Difference Test
## 
##         Df    AIC    BIC   Chisq Chisq diff Df diff Pr(>Chisq)    
## isatfit 20 1876.8 1946.8  21.839                                  
## nullfit 36 2108.1 2131.4 285.154     263.31      16  < 2.2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

anova(icfafit, isatfit)

## Chi-Squared Difference Test
## 
##         Df    AIC    BIC  Chisq Chisq diff Df diff Pr(>Chisq)  
## isatfit 20 1876.8 1946.8 21.839                                
## icfafit 27 1875.5 1925.1 34.588     12.749       7    0.07846 .
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Chapter 6

Dyadic CFA

Data load

Chapter 6에서 사용하는 데이터는 dyads member들이 성별로 구별이 가능함. 즉, distinguishable by gender임.
하지만 설명을 위해서 indistinguishable로 가정함.

setwd("C:/Users/LG/Documents/Dyadic Analysis Book Data")

library(haven)
indis <- read_spss("chapter6  Acitelli.sav")

names(indis) <- c("x1","x2","x3","x4","x5","x6","LM") 
#x1=아내 친밀성, x2=아내 헌신, x3=아내 만족도, x4=남편 친밀성, x5=남편 헌신, x6=남편 만족도, LM=결혼기간
#각 indicator는 4점 척도임.
#LM = the length of the marriage in years로 between-dyads variable임. 
indis

## # A tibble: 148 x 7
##             x1         x2         x3         x4        x5          x6    LM
##      <dbl+lbl>  <dbl+lbl>  <dbl+lbl>  <dbl+lbl> <dbl+lbl>   <dbl+lbl> <dbl>
##  1 4 [vclose]  4 [very]   4 [vsat]   4 [vclose] 4 [very]  4 [vsat]     2   
##  2 3 [swtclos~ 4 [very]   4 [vsat]   3 [swtclo~ 4 [very]  3 [smwtsat]  3.83
##  3 4 [vclose]  3 [fairly] 3 [smwtsa~ 2 [nvclos~ 3 [fairl~ 3 [smwtsat]  7.17
##  4 4 [vclose]  4 [very]   4 [vsat]   3 [swtclo~ 3 [fairl~ 3 [smwtsat] 22.8 
##  5 3 [swtclos~ 3 [fairly] 3 [smwtsa~ 3 [swtclo~ 4 [very]  2 [smwtdis~ 15.5 
##  6 3 [swtclos~ 3 [fairly] 3 [smwtsa~ 3 [swtclo~ 3 [fairl~ 3 [smwtsat]  7.92
##  7 4 [vclose]  4 [very]   4 [vsat]   4 [vclose] 4 [very]  3 [smwtsat]  9.67
##  8 4 [vclose]  4 [very]   4 [vsat]   4 [vclose] 4 [very]  3 [smwtsat] 11.2 
##  9 4 [vclose]  4 [very]   4 [vsat]   3 [swtclo~ 3 [fairl~ 3 [smwtsat]  1.92
## 10 3 [swtclos~ 4 [very]   4 [vsat]   3 [swtclo~ 4 [very]  4 [vsat]     2.17
## # ... with 138 more rows

library(psych)
describe(indis)

##    vars   n  mean   sd median trimmed   mad min   max range  skew kurtosis
## x1    1 148  3.43 0.73   4.00    3.56  0.00   1  4.00  3.00 -1.17     0.95
## x2    2 148  3.70 0.58   4.00    3.82  0.00   2  4.00  2.00 -1.73     1.91
## x3    3 148  3.61 0.66   4.00    3.73  0.00   1  4.00  3.00 -1.69     2.59
## x4    4 148  3.43 0.64   4.00    3.51  0.00   2  4.00  2.00 -0.65    -0.59
## x5    5 148  3.74 0.52   4.00    3.84  0.00   1  4.00  3.00 -2.19     5.45
## x6    6 148  3.59 0.64   4.00    3.69  0.00   1  4.00  3.00 -1.58     2.53
## LM    7 148 11.21 7.72  10.12   10.93 10.13   0 26.25 26.25  0.23    -1.22
##      se
## x1 0.06
## x2 0.05
## x3 0.05
## x4 0.05
## x5 0.04
## x6 0.05
## LM 0.63

result1 <- describe(indis)

indis1 <- indis[,1:3]
indis2 <- indis[,4:6]
names(indis2) <- c("x1","x2","x3")
doubleindis <- rbind(indis1, indis2)

describe(doubleindis)

##    vars   n mean   sd median trimmed mad min max range  skew kurtosis   se
## x1    1 296 3.43 0.69      4    3.54   0   1   4     3 -0.97     0.45 0.04
## x2    2 296 3.72 0.55      4    3.83   0   1   4     3 -1.95     3.43 0.03
## x3    3 296 3.60 0.65      4    3.71   0   1   4     3 -1.64     2.60 0.04

result2 <- describe(doubleindis)

DS <- rbind(result1[,3:4], result2[,3:4])
print(DS)

##      mean   sd
## x1   3.43 0.73
## x2   3.70 0.58
## x3   3.61 0.66
## x4   3.43 0.64
## x5   3.74 0.52
## x6   3.59 0.64
## LM  11.21 7.72
## x11  3.43 0.69
## x21  3.72 0.55
## x31  3.60 0.65

Correlation

corindis <- indis[,1:6]
corindis

## # A tibble: 148 x 6
##              x1         x2          x3           x4        x5            x6
##       <dbl+lbl>  <dbl+lbl>   <dbl+lbl>    <dbl+lbl> <dbl+lbl>     <dbl+lbl>
##  1 4 [vclose]   4 [very]   4 [vsat]    4 [vclose]   4 [very]  4 [vsat]     
##  2 3 [swtclose] 4 [very]   4 [vsat]    3 [swtclose] 4 [very]  3 [smwtsat]  
##  3 4 [vclose]   3 [fairly] 3 [smwtsat] 2 [nvclose]  3 [fairl~ 3 [smwtsat]  
##  4 4 [vclose]   4 [very]   4 [vsat]    3 [swtclose] 3 [fairl~ 3 [smwtsat]  
##  5 3 [swtclose] 3 [fairly] 3 [smwtsat] 3 [swtclose] 4 [very]  2 [smwtdissa~
##  6 3 [swtclose] 3 [fairly] 3 [smwtsat] 3 [swtclose] 3 [fairl~ 3 [smwtsat]  
##  7 4 [vclose]   4 [very]   4 [vsat]    4 [vclose]   4 [very]  3 [smwtsat]  
##  8 4 [vclose]   4 [very]   4 [vsat]    4 [vclose]   4 [very]  3 [smwtsat]  
##  9 4 [vclose]   4 [very]   4 [vsat]    3 [swtclose] 3 [fairl~ 3 [smwtsat]  
## 10 3 [swtclose] 4 [very]   4 [vsat]    3 [swtclose] 4 [very]  4 [vsat]     
## # ... with 138 more rows

corta <- cor(corindis, corindis)
ccor <- corta

corta[upper.tri(corta, diag = F)] <- NA #uppertri 지우는 방법
corta

##           x1        x2        x3        x4        x5 x6
## x1 1.0000000        NA        NA        NA        NA NA
## x2 0.5383563 1.0000000        NA        NA        NA NA
## x3 0.4983389 0.6868043 1.0000000        NA        NA NA
## x4 0.4479787 0.4433826 0.3838977 1.0000000        NA NA
## x5 0.3283128 0.5262843 0.4381323 0.5318835 1.0000000 NA
## x6 0.4582692 0.5420793 0.5540184 0.5825756 0.6181255  1

round(corta,3)

##       x1    x2    x3    x4    x5 x6
## x1 1.000    NA    NA    NA    NA NA
## x2 0.538 1.000    NA    NA    NA NA
## x3 0.498 0.687 1.000    NA    NA NA
## x4 0.448 0.443 0.384 1.000    NA NA
## x5 0.328 0.526 0.438 0.532 1.000 NA
## x6 0.458 0.542 0.554 0.583 0.618  1

Covariate

library(lavaan)

ccormean <- c(3.432432, 3.695946, 3.608108, 3.425676, 3.743243, 3.587838)
ccorsd <- c(0.7299442, 0.5792181, 0.6560564, 0.6398471, 0.5232242, 0.6381207)
covv <- cor2cov(ccor, ccorsd)
covv

##           x1        x2        x3        x4        x5        x6
## x1 0.5328185 0.2276154 0.2386468 0.2092296 0.1253907 0.2134584
## x2 0.2276154 0.3354936 0.2609855 0.1643225 0.1594962 0.2003585
## x3 0.2386468 0.2609855 0.4304100 0.1611509 0.1503953 0.2319360
## x4 0.2092296 0.1643225 0.1611509 0.4094043 0.1780658 0.2378654
## x5 0.1253907 0.1594962 0.1503953 0.1780658 0.2737636 0.2063799
## x6 0.2134584 0.2003585 0.2319360 0.2378654 0.2063799 0.4071980

round(covv, 3)

##       x1    x2    x3    x4    x5    x6
## x1 0.533 0.228 0.239 0.209 0.125 0.213
## x2 0.228 0.335 0.261 0.164 0.159 0.200
## x3 0.239 0.261 0.430 0.161 0.150 0.232
## x4 0.209 0.164 0.161 0.409 0.178 0.238
## x5 0.125 0.159 0.150 0.178 0.274 0.206
## x6 0.213 0.200 0.232 0.238 0.206 0.407

Distinguishable Dyads CFA

Distinguishable Model

lavaan syntax : =~ : indicator와 factor연결, ~ 1 : 절편, ~~ : 상관, 공분산을 나타냄.
figure5.1, figure6.3을 추정하는 distinguishable dyads CFA Model임.
남편과 아내 모두 동일하게 측정한 indicator에 대해서 factor loading을 고정하는 것은 타당한 제약임.
x1=아내 친밀성, x2=아내 헌신, x3=아내 만족도, x4=남편 친밀성, x5=남편 헌신, x6=남편 만족도

library(lavaan)

ddyadmodel <- '
wife =~ 1*x1 + a*x2 + b*x3
husband =~ 1*x4 + a*x5 + b*x6

wife ~~ husband

x1 ~ 1
x2 ~ 1
x3 ~ 1
x4 ~ 1
x5 ~ 1
x6 ~ 1

x1 ~~ x4
x2 ~~ x5
x3 ~~ x6
'
ddyadfit <- cfa(ddyadmodel, sample.cov = covv, sample.nobs = 148, sample.mean = ccormean, 
                std.lv = FALSE, mimic = "Mplus") 
summary(ddyadfit, fit.measures = TRUE)

## lavaan 0.6-5 ended normally after 35 iterations
## 
##   Estimator                                         ML
##   Optimization method                           NLMINB
##   Number of free parameters                         22
##   Number of equality constraints                     2
##   Row rank of the constraints matrix                 2
##                                                       
##   Number of observations                           148
##                                                       
## Model Test User Model:
##                                                       
##   Test statistic                                 7.124
##   Degrees of freedom                                 7
##   P-value (Chi-square)                           0.416
## 
## Model Test Baseline Model:
## 
##   Test statistic                               388.633
##   Degrees of freedom                                15
##   P-value                                        0.000
## 
## User Model versus Baseline Model:
## 
##   Comparative Fit Index (CFI)                    1.000
##   Tucker-Lewis Index (TLI)                       0.999
## 
## Loglikelihood and Information Criteria:
## 
##   Loglikelihood user model (H0)               -648.024
##   Loglikelihood unrestricted model (H1)       -644.462
##                                                       
##   Akaike (AIC)                                1336.048
##   Bayesian (BIC)                              1395.992
##   Sample-size adjusted Bayesian (BIC)         1332.699
## 
## Root Mean Square Error of Approximation:
## 
##   RMSEA                                          0.011
##   90 Percent confidence interval - lower         0.000
##   90 Percent confidence interval - upper         0.102
##   P-value RMSEA <= 0.05                          0.646
## 
## Standardized Root Mean Square Residual:
## 
##   SRMR                                           0.077
## 
## Parameter Estimates:
## 
##   Information                                 Observed
##   Observed information based on                Hessian
##   Standard errors                             Standard
## 
## Latent Variables:
##                    Estimate  Std.Err  z-value  P(>|z|)
##   wife =~                                             
##     x1                1.000                           
##     x2         (a)    0.974    0.097   10.089    0.000
##     x3         (b)    1.155    0.115   10.002    0.000
##   husband =~                                          
##     x4                1.000                           
##     x5         (a)    0.974    0.097   10.089    0.000
##     x6         (b)    1.155    0.115   10.002    0.000
## 
## Covariances:
##                    Estimate  Std.Err  z-value  P(>|z|)
##   wife ~~                                             
##     husband           0.156    0.034    4.628    0.000
##  .x1 ~~                                               
##    .x4                0.053    0.025    2.086    0.037
##  .x2 ~~                                               
##    .x5                0.016    0.015    1.093    0.274
##  .x3 ~~                                               
##    .x6                0.018    0.020    0.883    0.377
## 
## Intercepts:
##                    Estimate  Std.Err  z-value  P(>|z|)
##    .x1                3.432    0.060   56.755    0.000
##    .x2                3.696    0.046   79.908    0.000
##    .x3                3.608    0.055   65.664    0.000
##    .x4                3.426    0.052   65.985    0.000
##    .x5                3.743    0.044   84.237    0.000
##    .x6                3.588    0.051   70.317    0.000
##     wife              0.000                           
##     husband           0.000                           
## 
## Variances:
##                    Estimate  Std.Err  z-value  P(>|z|)
##    .x1                0.314    0.042    7.504    0.000
##    .x2                0.101    0.022    4.681    0.000
##    .x3                0.144    0.029    4.991    0.000
##    .x4                0.210    0.030    7.014    0.000
##    .x5                0.113    0.020    5.679    0.000
##    .x6                0.134    0.028    4.861    0.000
##     wife              0.227    0.047    4.883    0.000
##     husband           0.189    0.038    4.932    0.000

standardizedSolution(ddyadfit)

##        lhs op     rhs est.std    se      z pvalue ci.lower ci.upper
## 1     wife =~      x1   0.648 0.048 13.420  0.000    0.553    0.742
## 2     wife =~      x2   0.825 0.043 19.261  0.000    0.741    0.909
## 3     wife =~      x3   0.823 0.039 21.147  0.000    0.747    0.900
## 4  husband =~      x4   0.687 0.049 13.995  0.000    0.591    0.784
## 5  husband =~      x5   0.782 0.041 19.075  0.000    0.702    0.862
## 6  husband =~      x6   0.808 0.045 18.094  0.000    0.720    0.895
## 7     wife ~~ husband   0.754 0.055 13.824  0.000    0.647    0.861
## 8       x1 ~1           4.665 0.267 17.481  0.000    4.142    5.188
## 9       x2 ~1           6.568 0.375 17.511  0.000    5.833    7.304
## 10      x3 ~1           5.398 0.320 16.853  0.000    4.770    6.025
## 11      x4 ~1           5.424 0.306 17.700  0.000    4.823    6.025
## 12      x5 ~1           6.924 0.403 17.177  0.000    6.134    7.714
## 13      x6 ~1           5.780 0.329 17.577  0.000    5.135    6.425
## 14      x1 ~~      x4   0.206 0.091  2.254  0.024    0.027    0.384
## 15      x2 ~~      x5   0.148 0.125  1.190  0.234   -0.096    0.393
## 16      x3 ~~      x6   0.128 0.134  0.954  0.340   -0.135    0.390
## 17      x1 ~~      x1   0.580 0.063  9.278  0.000    0.458    0.703
## 18      x2 ~~      x2   0.320 0.071  4.525  0.000    0.181    0.458
## 19      x3 ~~      x3   0.322 0.064  5.027  0.000    0.197    0.448
## 20      x4 ~~      x4   0.527 0.068  7.807  0.000    0.395    0.660
## 21      x5 ~~      x5   0.388 0.064  6.054  0.000    0.263    0.514
## 22      x6 ~~      x6   0.348 0.072  4.822  0.000    0.206    0.489
## 23    wife ~~    wife   1.000 0.000     NA     NA    1.000    1.000
## 24 husband ~~ husband   1.000 0.000     NA     NA    1.000    1.000
## 25    wife ~1           0.000 0.000     NA     NA    0.000    0.000
## 26 husband ~1           0.000 0.000     NA     NA    0.000    0.000

CFA Model

library(lavaan) 

cfamodel <- '
wife =~ 1*x1 + x2 + x3
husband =~ 1*x4 + x5 + x6

wife ~~ husband

x1 ~ 1
x2 ~ 1
x3 ~ 1
x4 ~ 1
x5 ~ 1
x6 ~ 1

x1 ~~ x4
x2 ~~ x5
x3 ~~ x6
'
cfafit <- cfa(cfamodel, sample.cov = covv, sample.nobs = 148, sample.mean = ccormean, 
              std.lv = FALSE, mimic = "Mplus") 
summary(cfafit, fit.measures = TRUE)

## lavaan 0.6-5 ended normally after 35 iterations
## 
##   Estimator                                         ML
##   Optimization method                           NLMINB
##   Number of free parameters                         22
##                                                       
##   Number of observations                           148
##                                                       
## Model Test User Model:
##                                                       
##   Test statistic                                 3.057
##   Degrees of freedom                                 5
##   P-value (Chi-square)                           0.691
## 
## Model Test Baseline Model:
## 
##   Test statistic                               388.633
##   Degrees of freedom                                15
##   P-value                                        0.000
## 
## User Model versus Baseline Model:
## 
##   Comparative Fit Index (CFI)                    1.000
##   Tucker-Lewis Index (TLI)                       1.016
## 
## Loglikelihood and Information Criteria:
## 
##   Loglikelihood user model (H0)               -645.990
##   Loglikelihood unrestricted model (H1)       -644.462
##                                                       
##   Akaike (AIC)                                1335.980
##   Bayesian (BIC)                              1401.919
##   Sample-size adjusted Bayesian (BIC)         1332.297
## 
## Root Mean Square Error of Approximation:
## 
##   RMSEA                                          0.000
##   90 Percent confidence interval - lower         0.000
##   90 Percent confidence interval - upper         0.087
##   P-value RMSEA <= 0.05                          0.826
## 
## Standardized Root Mean Square Residual:
## 
##   SRMR                                           0.015
## 
## Parameter Estimates:
## 
##   Information                                 Observed
##   Observed information based on                Hessian
##   Standard errors                             Standard
## 
## Latent Variables:
##                    Estimate  Std.Err  z-value  P(>|z|)
##   wife =~                                             
##     x1                1.000                           
##     x2                1.070    0.142    7.537    0.000
##     x3                1.117    0.149    7.480    0.000
##   husband =~                                          
##     x4                1.000                           
##     x5                0.874    0.113    7.733    0.000
##     x6                1.209    0.153    7.900    0.000
## 
## Covariances:
##                    Estimate  Std.Err  z-value  P(>|z|)
##   wife ~~                                             
##     husband           0.154    0.033    4.594    0.000
##  .x1 ~~                                               
##    .x4                0.050    0.025    2.005    0.045
##  .x2 ~~                                               
##    .x5                0.018    0.015    1.231    0.218
##  .x3 ~~                                               
##    .x6                0.021    0.020    1.031    0.302
## 
## Intercepts:
##                    Estimate  Std.Err  z-value  P(>|z|)
##    .x1                3.432    0.060   57.335    0.000
##    .x2                3.696    0.048   77.803    0.000
##    .x3                3.608    0.054   67.223    0.000
##    .x4                3.426    0.052   65.514    0.000
##    .x5                3.743    0.043   87.275    0.000
##    .x6                3.588    0.052   68.602    0.000
##     wife              0.000                           
##     husband           0.000                           
## 
## Variances:
##                    Estimate  Std.Err  z-value  P(>|z|)
##    .x1                0.315    0.042    7.454    0.000
##    .x2                0.088    0.023    3.783    0.000
##    .x3                0.158    0.029    5.356    0.000
##    .x4                0.209    0.031    6.819    0.000
##    .x5                0.123    0.020    6.054    0.000
##    .x6                0.119    0.029    4.042    0.000
##     wife              0.215    0.054    4.010    0.000
##     husband           0.196    0.044    4.419    0.000

standardizedSolution(cfafit)

##        lhs op     rhs est.std    se      z pvalue ci.lower ci.upper
## 1     wife =~      x1   0.637 0.058 11.037  0.000    0.524    0.750
## 2     wife =~      x2   0.859 0.042 20.464  0.000    0.777    0.941
## 3     wife =~      x3   0.793 0.046 17.298  0.000    0.704    0.883
## 4  husband =~      x4   0.695 0.054 12.906  0.000    0.590    0.801
## 5  husband =~      x5   0.741 0.051 14.475  0.000    0.641    0.841
## 6  husband =~      x6   0.840 0.045 18.648  0.000    0.752    0.929
## 7     wife ~~ husband   0.750 0.054 13.816  0.000    0.643    0.856
## 8       x1 ~1           4.713 0.286 16.481  0.000    4.152    5.273
## 9       x2 ~1           6.395 0.381 16.794  0.000    5.649    7.142
## 10      x3 ~1           5.526 0.331 16.686  0.000    4.877    6.175
## 11      x4 ~1           5.385 0.322 16.708  0.000    4.754    6.017
## 12      x5 ~1           7.174 0.425 16.876  0.000    6.341    8.007
## 13      x6 ~1           5.639 0.338 16.689  0.000    4.977    6.301
## 14      x1 ~~      x4   0.196 0.091  2.152  0.031    0.018    0.375
## 15      x2 ~~      x5   0.173 0.127  1.365  0.172   -0.075    0.422
## 16      x3 ~~      x6   0.152 0.134  1.136  0.256   -0.111    0.415
## 17      x1 ~~      x1   0.594 0.074  8.086  0.000    0.450    0.738
## 18      x2 ~~      x2   0.262 0.072  3.639  0.000    0.121    0.404
## 19      x3 ~~      x3   0.371 0.073  5.091  0.000    0.228    0.513
## 20      x4 ~~      x4   0.517 0.075  6.899  0.000    0.370    0.663
## 21      x5 ~~      x5   0.451 0.076  5.942  0.000    0.302    0.600
## 22      x6 ~~      x6   0.294 0.076  3.881  0.000    0.145    0.442
## 23    wife ~~    wife   1.000 0.000     NA     NA    1.000    1.000
## 24 husband ~~ husband   1.000 0.000     NA     NA    1.000    1.000
## 25    wife ~1           0.000 0.000     NA     NA    0.000    0.000
## 26 husband ~1           0.000 0.000     NA     NA    0.000    0.000

Model Comparistion

기본 CFA Model과 distinguishable dyads CFA Model을 비교함.
ANOVA결과 남편, 아내의 factor loading을 고정함.
남편과 아내는 동일한 indicator의 결합으로 구성되어 있고, 각 indicator가 남편과 아내에게 가지는 의미 또한 동일하다고 할 수 있음.

anova(cfafit, ddyadfit)

## Chi-Squared Difference Test
## 
##          Df  AIC    BIC  Chisq Chisq diff Df diff Pr(>Chisq)
## cfafit    5 1336 1401.9 3.0571                              
## ddyadfit  7 1336 1396.0 7.1244     4.0674       2     0.1309

Indistinguishable Dyads CFA

I-SAT Model

I-SAT Model결과 dyads member들은 성별로 구분할 수 없고, dyads member를 indistinguishable dyads로 다뤄야함.
Dyads분석을 위해서 dyads가 indistinguishable인지 distinguishable인지 확인하는 것이 필요함.

library(lavaan)

ISATmodel <- '
x1 ~ m1*1 #Means
x4 ~ m1*1

x2 ~ m2*1
x5 ~ m2*1

x3 ~ m3*1
x6 ~ m3*1

x1 ~~ v1*x1 #Variances
x4 ~~ v1*x4

x2 ~~ v2*x2
x5 ~~ v2*x5

x3 ~~ v3*x3
x6 ~~ v3*x6

x1 ~~ x4 #Error correlation
x2 ~~ x5
x3 ~~ x6

x1 ~~ c12*x2 #Within member correlation
x1 ~~ c13*x3
x2 ~~ c23*x3

x4 ~~ c12*x5
x4 ~~ c13*x6
x5 ~~ c23*x6

x1 ~~ cc12*x5 #Between member correlation 
x1 ~~ cc13*x6
x2 ~~ cc23*x6

x4 ~~ cc12*x2
x4 ~~ cc13*x3
x5 ~~ cc23*x3
'
ISATfit <- cfa(ISATmodel, sample.cov = covv, sample.nobs = 148, sample.mean = ccormean, 
               std.lv = F, mimic = "Mplus") 
summary(ISATfit, fit.measures = TRUE)

## lavaan 0.6-5 ended normally after 57 iterations
## 
##   Estimator                                         ML
##   Optimization method                           NLMINB
##   Number of free parameters                         27
##   Number of equality constraints                    12
##   Row rank of the constraints matrix                12
##                                                       
##   Number of observations                           148
##                                                       
## Model Test User Model:
##                                                       
##   Test statistic                                14.919
##   Degrees of freedom                                12
##   P-value (Chi-square)                           0.246
## 
## Model Test Baseline Model:
## 
##   Test statistic                               388.633
##   Degrees of freedom                                15
##   P-value                                        0.000
## 
## User Model versus Baseline Model:
## 
##   Comparative Fit Index (CFI)                    0.992
##   Tucker-Lewis Index (TLI)                       0.990
## 
## Loglikelihood and Information Criteria:
## 
##   Loglikelihood user model (H0)               -651.921
##   Loglikelihood unrestricted model (H1)       -644.462
##                                                       
##   Akaike (AIC)                                1333.842
##   Bayesian (BIC)                              1378.800
##   Sample-size adjusted Bayesian (BIC)         1331.331
## 
## Root Mean Square Error of Approximation:
## 
##   RMSEA                                          0.041
##   90 Percent confidence interval - lower         0.000
##   90 Percent confidence interval - upper         0.098
##   P-value RMSEA <= 0.05                          0.546
## 
## Standardized Root Mean Square Residual:
## 
##   SRMR                                           0.150
## 
## Parameter Estimates:
## 
##   Information                                 Observed
##   Observed information based on                Hessian
##   Standard errors                             Standard
## 
## Covariances:
##                    Estimate  Std.Err  z-value  P(>|z|)
##   x1 ~~                                               
##     x4                0.208    0.042    4.938    0.000
##   x2 ~~                                               
##     x5                0.158    0.028    5.619    0.000
##   x3 ~~                                               
##     x6                0.230    0.039    5.891    0.000
##   x1 ~~                                               
##     x2       (c12)    0.201    0.028    7.133    0.000
##     x3       (c13)    0.237    0.034    7.055    0.000
##   x2 ~~                                               
##     x3       (c23)    0.232    0.029    8.032    0.000
##   x4 ~~                                               
##     x5       (c12)    0.201    0.028    7.133    0.000
##     x6       (c13)    0.237    0.034    7.055    0.000
##   x5 ~~                                               
##     x6       (c23)    0.232    0.029    8.032    0.000
##   x1 ~~                                               
##     x5      (cc12)    0.144    0.028    5.099    0.000
##     x6      (cc13)    0.186    0.034    5.544    0.000
##   x2 ~~                                               
##     x6      (cc23)    0.174    0.029    6.043    0.000
##   x4 ~~                                               
##     x2      (cc12)    0.144    0.028    5.099    0.000
##     x3      (cc13)    0.186    0.034    5.544    0.000
##   x5 ~~                                               
##     x3      (cc23)    0.174    0.029    6.043    0.000
## 
## Intercepts:
##                    Estimate  Std.Err  z-value  P(>|z|)
##     x1        (m1)    3.429    0.048   71.768    0.000
##     x4        (m1)    3.429    0.048   71.768    0.000
##     x2        (m2)    3.720    0.039   94.253    0.000
##     x5        (m2)    3.720    0.039   94.253    0.000
##     x3        (m3)    3.598    0.047   76.997    0.000
##     x6        (m3)    3.598    0.047   76.997    0.000
## 
## Variances:
##                    Estimate  Std.Err  z-value  P(>|z|)
##     x1        (v1)    0.468    0.042   11.118    0.000
##     x4        (v1)    0.468    0.042   11.118    0.000
##     x2        (v2)    0.303    0.028   10.790    0.000
##     x5        (v2)    0.303    0.028   10.790    0.000
##     x3        (v3)    0.416    0.039   10.644    0.000
##     x6        (v3)    0.416    0.039   10.644    0.000

NULL Model

library(lavaan)

NULLmodel <- '
x1 ~ m1*1 #Means 
x4 ~ m1*1
x2 ~ m2*1
x5 ~ m2*1
x3 ~ m3*1
x6 ~ m3*1

x1 ~~ r1*x1 #Variance
x2 ~~ r2*x2
x3 ~~ r3*x3
x4 ~~ r1*x4
x5 ~~ r2*x5
x6 ~~ r3*x6

x1 ~~ 0*x2 #Covariate
x1 ~~ 0*x3 
x1 ~~ 0*x4
x1 ~~ 0*x5
x1 ~~ 0*x6
x2 ~~ 0*x3
x2 ~~ 0*x4
x2 ~~ 0*x5
x2 ~~ 0*x6
x3 ~~ 0*x4
x3 ~~ 0*x5
x3 ~~ 0*x6 
x4 ~~ 0*x5
x4 ~~ 0*x6
x5 ~~ 0*x6
'
NULLfit <- cfa(NULLmodel, sample.cov = covv, sample.nobs = 148, sample.mean = ccormean, 
               std.lv = F, mimic = "Mplus")
summary(NULLfit, fit.measures = TRUE)

## lavaan 0.6-5 ended normally after 22 iterations
## 
##   Estimator                                         ML
##   Optimization method                           NLMINB
##   Number of free parameters                         12
##   Number of equality constraints                     6
##   Row rank of the constraints matrix                 6
##                                                       
##   Number of observations                           148
##                                                       
## Model Test User Model:
##                                                       
##   Test statistic                               393.462
##   Degrees of freedom                                21
##   P-value (Chi-square)                           0.000
## 
## Model Test Baseline Model:
## 
##   Test statistic                               388.633
##   Degrees of freedom                                15
##   P-value                                        0.000
## 
## User Model versus Baseline Model:
## 
##   Comparative Fit Index (CFI)                    0.003
##   Tucker-Lewis Index (TLI)                       0.288
## 
## Loglikelihood and Information Criteria:
## 
##   Loglikelihood user model (H0)               -841.193
##   Loglikelihood unrestricted model (H1)       -644.462
##                                                       
##   Akaike (AIC)                                1694.386
##   Bayesian (BIC)                              1712.369
##   Sample-size adjusted Bayesian (BIC)         1693.381
## 
## Root Mean Square Error of Approximation:
## 
##   RMSEA                                          0.346
##   90 Percent confidence interval - lower         0.317
##   90 Percent confidence interval - upper         0.377
##   P-value RMSEA <= 0.05                          0.000
## 
## Standardized Root Mean Square Residual:
## 
##   SRMR                                           0.410
## 
## Parameter Estimates:
## 
##   Information                                 Observed
##   Observed information based on                Hessian
##   Standard errors                             Standard
## 
## Covariances:
##                    Estimate  Std.Err  z-value  P(>|z|)
##   x1 ~~                                               
##     x2                0.000                           
##     x3                0.000                           
##     x4                0.000                           
##     x5                0.000                           
##     x6                0.000                           
##   x2 ~~                                               
##     x3                0.000                           
##     x4                0.000                           
##     x5                0.000                           
##     x6                0.000                           
##   x3 ~~                                               
##     x4                0.000                           
##     x5                0.000                           
##     x6                0.000                           
##   x4 ~~                                               
##     x5                0.000                           
##     x6                0.000                           
##   x5 ~~                                               
##     x6                0.000                           
## 
## Intercepts:
##                    Estimate  Std.Err  z-value  P(>|z|)
##     x1        (m1)    3.429    0.040   86.243    0.000
##     x4        (m1)    3.429    0.040   86.243    0.000
##     x2        (m2)    3.720    0.032  116.233    0.000
##     x5        (m2)    3.720    0.032  116.233    0.000
##     x3        (m3)    3.598    0.037   95.966    0.000
##     x6        (m3)    3.598    0.037   95.966    0.000
## 
## Variances:
##                    Estimate  Std.Err  z-value  P(>|z|)
##     x1        (r1)    0.468    0.038   12.166    0.000
##     x2        (r2)    0.303    0.025   12.166    0.000
##     x3        (r3)    0.416    0.034   12.166    0.000
##     x4        (r1)    0.468    0.038   12.166    0.000
##     x5        (r2)    0.303    0.025   12.166    0.000
##     x6        (r3)    0.416    0.034   12.166    0.000

Indistinguishable Model

library(lavaan)

idyadmodel <-'
wife =~ 1*x1 + a*x2 + b*x3 #factor loading
husband =~ 1*x4 + a*x5 + b*x6

wife ~ 0*1 #factor means
husband ~ 0*1

wife ~~ i*husband #factor covariate

wife ~~ h*wife #factor variance
husband ~~ h*husband

x1 ~ m1*1 #intercept
x2 ~ m2*1
x3 ~ m3*1
x4 ~ m1*1
x5 ~ m2*1
x6 ~ m3*1

x1 ~~ v1*x1 #residual variance
x2 ~~ v2*x2
x3 ~~ v3*x3
x4 ~~ v1*x4
x5 ~~ v2*x5
x6 ~~ v3*x6

x1 ~~ c1*x4 #residual correlation
x2 ~~ c2*x5
x3 ~~ c3*x6
'
idyadfit <- cfa(idyadmodel, sample.cov = covv, sample.nobs = 148, sample.mean = ccormean, 
                std.lv = F, mimic = "Mplus")
summary(idyadfit, fit.measures = TRUE)

## lavaan 0.6-5 ended normally after 40 iterations
## 
##   Estimator                                         ML
##   Optimization method                           NLMINB
##   Number of free parameters                         22
##   Number of equality constraints                     9
##   Row rank of the constraints matrix                 9
##                                                       
##   Number of observations                           148
##                                                       
## Model Test User Model:
##                                                       
##   Test statistic                                15.702
##   Degrees of freedom                                14
##   P-value (Chi-square)                           0.332
## 
## Model Test Baseline Model:
## 
##   Test statistic                               388.633
##   Degrees of freedom                                15
##   P-value                                        0.000
## 
## User Model versus Baseline Model:
## 
##   Comparative Fit Index (CFI)                    0.995
##   Tucker-Lewis Index (TLI)                       0.995
## 
## Loglikelihood and Information Criteria:
## 
##   Loglikelihood user model (H0)               -652.313
##   Loglikelihood unrestricted model (H1)       -644.462
##                                                       
##   Akaike (AIC)                                1330.625
##   Bayesian (BIC)                              1369.589
##   Sample-size adjusted Bayesian (BIC)         1328.449
## 
## Root Mean Square Error of Approximation:
## 
##   RMSEA                                          0.029
##   90 Percent confidence interval - lower         0.000
##   90 Percent confidence interval - upper         0.087
##   P-value RMSEA <= 0.05                          0.659
## 
## Standardized Root Mean Square Residual:
## 
##   SRMR                                           0.150
## 
## Parameter Estimates:
## 
##   Information                                 Observed
##   Observed information based on                Hessian
##   Standard errors                             Standard
## 
## Latent Variables:
##                    Estimate  Std.Err  z-value  P(>|z|)
##   wife =~                                             
##     x1                1.000                           
##     x2         (a)    0.971    0.097    9.996    0.000
##     x3         (b)    1.164    0.117    9.906    0.000
##   husband =~                                          
##     x4                1.000                           
##     x5         (a)    0.971    0.097    9.996    0.000
##     x6         (b)    1.164    0.117    9.906    0.000
## 
## Covariances:
##                    Estimate  Std.Err  z-value  P(>|z|)
##   wife ~~                                             
##     husband    (i)    0.155    0.034    4.579    0.000
##  .x1 ~~                                               
##    .x4        (c1)    0.053    0.026    2.079    0.038
##  .x2 ~~                                               
##    .x5        (c2)    0.016    0.015    1.068    0.285
##  .x3 ~~                                               
##    .x6        (c3)    0.016    0.020    0.812    0.417
## 
## Intercepts:
##                    Estimate  Std.Err  z-value  P(>|z|)
##     wife              0.000                           
##     husband           0.000                           
##    .x1        (m1)    3.429    0.048   71.755    0.000
##    .x2        (m2)    3.720    0.040   93.847    0.000
##    .x3        (m3)    3.598    0.047   77.297    0.000
##    .x4        (m1)    3.429    0.048   71.755    0.000
##    .x5        (m2)    3.720    0.040   93.847    0.000
##    .x6        (m3)    3.598    0.047   77.297    0.000
## 
## Variances:
##                    Estimate  Std.Err  z-value  P(>|z|)
##     wife       (h)    0.205    0.039    5.264    0.000
##     husband    (h)    0.205    0.039    5.264    0.000
##    .x1        (v1)    0.263    0.027    9.899    0.000
##    .x2        (v2)    0.110    0.016    6.676    0.000
##    .x3        (v3)    0.137    0.023    6.042    0.000
##    .x4        (v1)    0.263    0.027    9.899    0.000
##    .x5        (v2)    0.110    0.016    6.676    0.000
##    .x6        (v3)    0.137    0.023    6.042    0.000

standardizedSolution(idyadfit)

##        lhs op     rhs est.std    se      z pvalue ci.lower ci.upper
## 1     wife =~      x1   0.662 0.044 15.221  0.000    0.577    0.748
## 2     wife =~      x2   0.798 0.037 21.763  0.000    0.726    0.870
## 3     wife =~      x3   0.818 0.036 22.593  0.000    0.747    0.889
## 4  husband =~      x4   0.662 0.044 15.221  0.000    0.577    0.748
## 5  husband =~      x5   0.798 0.037 21.763  0.000    0.726    0.870
## 6  husband =~      x6   0.818 0.036 22.593  0.000    0.747    0.889
## 7     wife ~1           0.000 0.000     NA     NA    0.000    0.000
## 8  husband ~1           0.000 0.000     NA     NA    0.000    0.000
## 9     wife ~~ husband   0.753 0.055 13.673  0.000    0.645    0.861
## 10    wife ~~    wife   1.000 0.000     NA     NA    1.000    1.000
## 11 husband ~~ husband   1.000 0.000     NA     NA    1.000    1.000
## 12      x1 ~1           5.013 0.236 21.267  0.000    4.551    5.475
## 13      x2 ~1           6.750 0.322 20.933  0.000    6.118    7.382
## 14      x3 ~1           5.582 0.271 20.615  0.000    5.051    6.112
## 15      x4 ~1           5.013 0.236 21.267  0.000    4.551    5.475
## 16      x5 ~1           6.750 0.322 20.933  0.000    6.118    7.382
## 17      x6 ~1           5.582 0.271 20.615  0.000    5.051    6.112
## 18      x1 ~~      x1   0.561 0.058  9.735  0.000    0.448    0.674
## 19      x2 ~~      x2   0.363 0.059  6.193  0.000    0.248    0.477
## 20      x3 ~~      x3   0.331 0.059  5.584  0.000    0.215    0.447
## 21      x4 ~~      x4   0.561 0.058  9.735  0.000    0.448    0.674
## 22      x5 ~~      x5   0.363 0.059  6.193  0.000    0.248    0.477
## 23      x6 ~~      x6   0.331 0.059  5.584  0.000    0.215    0.447
## 24      x1 ~~      x4   0.203 0.091  2.243  0.025    0.026    0.381
## 25      x2 ~~      x5   0.143 0.123  1.157  0.247   -0.099    0.384
## 26      x3 ~~      x6   0.120 0.137  0.875  0.382   -0.149    0.388

Model Comparistion

anova(ISATfit, NULLfit)

## Chi-Squared Difference Test
## 
##         Df    AIC    BIC   Chisq Chisq diff Df diff Pr(>Chisq)    
## ISATfit 12 1333.8 1378.8  14.919                                  
## NULLfit 21 1694.4 1712.4 393.462     378.54       9  < 2.2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

anova(ISATfit, idyadfit)

## Chi-Squared Difference Test
## 
##          Df    AIC    BIC  Chisq Chisq diff Df diff Pr(>Chisq)
## ISATfit  12 1333.8 1378.8 14.919                              
## idyadfit 14 1330.6 1369.6 15.702    0.78301       2      0.676

Dyadic Data Analysis

Chapter 4

Dyadic MLM

Data load

Indistinguishable Dyads MLM

Distinguishable Dyads MLM

two intercept model (intercept)

Two intercept model (predictor)

Test Distinguishablity MLM

Chapter 5

Dyadic CFA

Data load

Indistinguishable Dyads CFA

Null Model

I-SAT Model

Indistinguishable dyads CFA Model

Chapter 6

Dyadic CFA

Data load

Correlation

Covariate

Distinguishable Dyads CFA

Distinguishable Model

CFA Model

Model Comparistion

Indistinguishable Dyads CFA

I-SAT Model

NULL Model

Indistinguishable Model

Model Comparistion