APIM with MLM
DATA
setwd("C:/Users/LG/Documents/Dyadic Analysis Book Data")
#load data
library(haven)
## Warning: package 'haven' was built under R version 3.5.3
two <- read_spss("chapter7 campbell_pairwise.sav")
two
## # A tibble: 196 x 6
## gender dyad A_Distress P_Distress Act_Neuro Part_Neuro
## <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 1 1 2.52 4.52 -0.928 0.787
## 2 1 2 2.68 3.8 -0.499 -0.642
## 3 1 3 3.12 4.32 0.215 -1.07
## 4 1 4 3.76 5 -0.642 -1.07
## 5 1 5 3.08 1.76 1.07 -0.928
## 6 1 8 3.56 2.92 1.79 -0.356
## 7 1 10 4.52 5.28 0.215 -0.928
## 8 1 11 3.8 4.8 0.0724 0.0724
## 9 1 12 4.32 4.04 1.93 -1.07
## 10 1 14 5 5.32 0.691 0.215
## # ... with 186 more rows
library(dplyr)
## Warning: package 'dplyr' was built under R version 3.5.3
##
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
##
## filter, lag
## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, union
rtwo <- arrange(two, dyad)
newtwo <- rtwo[c(2,1,3,4,5,6)]
newtwo
## # A tibble: 196 x 6
## dyad gender A_Distress P_Distress Act_Neuro Part_Neuro
## <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 1 1 2.52 4.52 -0.928 0.787
## 2 1 -1 3.64 2.52 0.787 -0.928
## 3 2 1 2.68 3.8 -0.499 -0.642
## 4 2 -1 3.2 2.68 -0.642 -0.499
## 5 3 1 3.12 4.32 0.215 -1.07
## 6 3 -1 3.2 3.12 -1.07 0.215
## 7 4 1 3.76 5 -0.642 -1.07
## 8 4 -1 3.12 3.76 -1.07 -0.642
## 9 5 1 3.08 1.76 1.07 -0.928
## 10 5 -1 3.52 3.08 -0.928 1.07
## # ... with 186 more rows
#gender coding을 -1,1을 바꿔야해서 새롭게 만든 gender (책과 실제 data에서 gender가 서로 다름)
gender <- rep(-1:1, 98)
gender <- gender[!(gender %in% c(0))]
gender
## [1] -1 1 -1 1 -1 1 -1 1 -1 1 -1 1 -1 1 -1 1 -1 1 -1 1 -1 1 -1
## [24] 1 -1 1 -1 1 -1 1 -1 1 -1 1 -1 1 -1 1 -1 1 -1 1 -1 1 -1 1
## [47] -1 1 -1 1 -1 1 -1 1 -1 1 -1 1 -1 1 -1 1 -1 1 -1 1 -1 1 -1
## [70] 1 -1 1 -1 1 -1 1 -1 1 -1 1 -1 1 -1 1 -1 1 -1 1 -1 1 -1 1
## [93] -1 1 -1 1 -1 1 -1 1 -1 1 -1 1 -1 1 -1 1 -1 1 -1 1 -1 1 -1
## [116] 1 -1 1 -1 1 -1 1 -1 1 -1 1 -1 1 -1 1 -1 1 -1 1 -1 1 -1 1
## [139] -1 1 -1 1 -1 1 -1 1 -1 1 -1 1 -1 1 -1 1 -1 1 -1 1 -1 1 -1
## [162] 1 -1 1 -1 1 -1 1 -1 1 -1 1 -1 1 -1 1 -1 1 -1 1 -1 1 -1 1
## [185] -1 1 -1 1 -1 1 -1 1 -1 1 -1 1
#변수 순서 바꾸기
cgender <- cbind(newtwo, gender)
ccdata <- cgender[c(1,7,3,4,5,6,7)]
DATA <- ccdata[1:6]
head(DATA) #1:female, -1:male
## dyad gender A_Distress P_Distress Act_Neuro Part_Neuro
## 1 1 -1 2.52 4.52 -0.92756 0.78674
## 2 1 1 3.64 2.52 0.78674 -0.92756
## 3 2 -1 2.68 3.80 -0.49896 -0.64186
## 4 2 1 3.20 2.68 -0.64186 -0.49896
## 5 3 -1 3.12 4.32 0.21524 -1.07046
## 6 3 1 3.20 3.12 -1.07046 0.21524
#dummay variable
m <- rep(1:0, 98)
f <- rep(0:1, 98)
dgender <- cbind(m,f)
colnames(dgender) <- c("M","F")
DDATA <- cbind(DATA, dgender)
head(DDATA)
## dyad gender A_Distress P_Distress Act_Neuro Part_Neuro M F
## 1 1 -1 2.52 4.52 -0.92756 0.78674 1 0
## 2 1 1 3.64 2.52 0.78674 -0.92756 0 1
## 3 2 -1 2.68 3.80 -0.49896 -0.64186 1 0
## 4 2 1 3.20 2.68 -0.64186 -0.49896 0 1
## 5 3 -1 3.12 4.32 0.21524 -1.07046 1 0
## 6 3 1 3.20 3.12 -1.07046 0.21524 0 1
#dummy variable correlation
cor(DDATA$F,DDATA$M)
## [1] -1
Distinguishable dyad
#interaction model
library(nlme)
##
## Attaching package: 'nlme'
## The following object is masked from 'package:dplyr':
##
## collapse
interact <- gls(A_Distress ~ gender + Act_Neuro + Part_Neuro + Act_Neuro:gender + Part_Neuro:gender,
data = DDATA, correlation = corCompSymm(form = ~1|dyad),
weight = varIdent(form = ~1|gender))
summary(interact)
## Generalized least squares fit by REML
## Model: A_Distress ~ gender + Act_Neuro + Part_Neuro + Act_Neuro:gender + Part_Neuro:gender
## Data: DDATA
## AIC BIC logLik
## 405.5982 434.8215 -193.7991
##
## Correlation Structure: Compound symmetry
## Formula: ~1 | dyad
## Parameter estimate(s):
## Rho
## 0.6734969
## Variance function:
## Structure: Different standard deviations per stratum
## Formula: ~1 | gender
## Parameter estimates:
## -1 1
## 1.0000000 0.8046827
##
## Coefficients:
## Value Std.Error t-value p-value
## (Intercept) 3.659444 0.08644876 42.33078 0.0000
## gender -0.032524 0.03926944 -0.82822 0.4086
## Act_Neuro 0.407937 0.07437800 5.48465 0.0000
## Part_Neuro 0.343205 0.07692513 4.46154 0.0000
## gender:Act_Neuro -0.023402 0.07240984 -0.32319 0.7469
## gender:Part_Neuro 0.203649 0.07502383 2.71445 0.0072
##
## Correlation:
## (Intr) gender Act_Nr Prt_Nr gn:A_N
## gender -0.284
## Act_Neuro 0.022 0.229
## Part_Neuro 0.162 -0.273 0.680
## gender:Act_Neuro 0.565 -0.120 -0.054 0.098
## gender:Part_Neuro -0.570 0.201 -0.115 -0.363 -0.635
##
## Standardized residuals:
## Min Q1 Med Q3 Max
## -2.23577604 -0.78278150 -0.08646031 0.64305392 2.37930761
##
## Residual standard error: 0.8064111
## Degrees of freedom: 196 total; 190 residual
#two intercept model
library(nlme)
genderfit <- gls(A_Distress ~ -1 + M + F + Act_Neuro:M + Act_Neuro:F + Part_Neuro:M + Part_Neuro:F,
data = DDATA, correlation = corCompSymm(form = ~1|dyad), weight = varIdent(form = ~1|gender))
summary(genderfit)
## Generalized least squares fit by REML
## Model: A_Distress ~ -1 + M + F + Act_Neuro:M + Act_Neuro:F + Part_Neuro:M + Part_Neuro:F
## Data: DDATA
## AIC BIC logLik
## 401.4394 430.6626 -191.7197
##
## Correlation Structure: Compound symmetry
## Formula: ~1 | dyad
## Parameter estimate(s):
## Rho
## 0.6734969
## Variance function:
## Structure: Different standard deviations per stratum
## Formula: ~1 | gender
## Parameter estimates:
## -1 1
## 1.0000000 0.8046827
##
## Coefficients:
## Value Std.Error t-value p-value
## M 3.691967 0.10461520 35.29093 0.0000
## F 3.626920 0.08418203 43.08425 0.0000
## M:Act_Neuro 0.431339 0.10659170 4.04665 0.0001
## F:Act_Neuro 0.384535 0.10093935 3.80957 0.0002
## M:Part_Neuro 0.139556 0.12543995 1.11253 0.2673
## F:Part_Neuro 0.546853 0.08577249 6.37563 0.0000
##
## Correlation:
## M F M:Ac_N F:Ac_N M:Pr_N
## F 0.673
## M:Act_Neuro -0.395 -0.266
## F:Act_Neuro 0.317 0.471 0.027
## M:Part_Neuro 0.471 0.317 0.040 0.673
## F:Part_Neuro -0.266 -0.395 0.673 0.040 0.027
##
## Standardized residuals:
## Min Q1 Med Q3 Max
## -2.23577604 -0.78278150 -0.08646031 0.64305392 2.37930761
##
## Residual standard error: 0.8064111
## Degrees of freedom: 196 total; 190 residual
# 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임. ~1 : 비교집단이 2개면 2-1로 1을 사용함.
# dyad간 상관이 있으니 corCompSymm을 사용하고, 성별간 분산을 다르게 추정할 수 있게 하기 위해서 varIdent를 사용함.
#lmer로 two intercept model 추정하기
library(lme4)
## Loading required package: Matrix
##
## Attaching package: 'lme4'
## The following object is masked from 'package:nlme':
##
## lmList
fit <- lmer(A_Distress ~ -1 + M + F + Act_Neuro:M + Act_Neuro:F + Part_Neuro:M + Part_Neuro:F + (1|dyad),
data = DDATA)
summary(fit)
## Linear mixed model fit by REML ['lmerMod']
## Formula:
## A_Distress ~ -1 + M + F + Act_Neuro:M + Act_Neuro:F + Part_Neuro:M +
## Part_Neuro:F + (1 | dyad)
## Data: DDATA
##
## REML criterion at convergence: 391.4
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.45432 -0.47763 0.04179 0.46508 2.85081
##
## Random effects:
## Groups Name Variance Std.Dev.
## dyad (Intercept) 0.3524 0.5937
## Residual 0.1833 0.4281
## Number of obs: 196, groups: dyad, 98
##
## Fixed effects:
## Estimate Std. Error t value
## M 3.69197 0.09495 38.883
## F 3.62692 0.09495 38.198
## M:Act_Neuro 0.43134 0.09674 4.459
## F:Act_Neuro 0.38454 0.11385 3.378
## M:Part_Neuro 0.13956 0.11385 1.226
## F:Part_Neuro 0.54685 0.09674 5.653
##
## Correlation of Fixed Effects:
## M F M:Ac_N F:Ac_N M:Pr_N
## F 0.658
## M:Act_Neuro -0.395 -0.260
## F:Act_Neuro 0.310 0.471 0.026
## M:Part_Neur 0.471 0.310 0.040 0.658
## F:Part_Neur -0.260 -0.395 0.658 0.040 0.026
Indistinguishable dyad
#Data
head(DDATA)
## dyad gender A_Distress P_Distress Act_Neuro Part_Neuro M F
## 1 1 -1 2.52 4.52 -0.92756 0.78674 1 0
## 2 1 1 3.64 2.52 0.78674 -0.92756 0 1
## 3 2 -1 2.68 3.80 -0.49896 -0.64186 1 0
## 4 2 1 3.20 2.68 -0.64186 -0.49896 0 1
## 5 3 -1 3.12 4.32 0.21524 -1.07046 1 0
## 6 3 1 3.20 3.12 -1.07046 0.21524 0 1
#lmer
library(lme4)
MLM_I <- lmer(A_Distress ~ Act_Neuro + Part_Neuro + (1|dyad), data = DDATA)
summary(MLM_I)
## Linear mixed model fit by REML ['lmerMod']
## Formula: A_Distress ~ Act_Neuro + Part_Neuro + (1 | dyad)
## Data: DDATA
##
## REML criterion at convergence: 395.7
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.41003 -0.49328 -0.00054 0.45183 2.68299
##
## Random effects:
## Groups Name Variance Std.Dev.
## dyad (Intercept) 0.3542 0.5951
## Residual 0.1973 0.4442
## Number of obs: 196, groups: dyad, 98
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 3.75182 0.06798 55.192
## Act_Neuro 0.43393 0.07324 5.925
## Part_Neuro 0.35536 0.07324 4.852
##
## Correlation of Fixed Effects:
## (Intr) Act_Nr
## Act_Neuro 0.000
## Part_Neuro 0.000 0.780
#gls
library(nlme)
MLM_Igls <- gls(A_Distress ~ Act_Neuro + Part_Neuro, data = DDATA, correlation = corCompSymm(form = ~1|dyad))
summary(MLM_Igls)
## Generalized least squares fit by REML
## Model: A_Distress ~ Act_Neuro + Part_Neuro
## Data: DDATA
## AIC BIC logLik
## 405.7365 422.05 -197.8683
##
## Correlation Structure: Compound symmetry
## Formula: ~1 | dyad
## Parameter estimate(s):
## Rho
## 0.642194
##
## Coefficients:
## Value Std.Error t-value p-value
## (Intercept) 3.751817 0.06797726 55.19224 0
## Act_Neuro 0.433930 0.07324068 5.92472 0
## Part_Neuro 0.355358 0.07324068 4.85192 0
##
## Correlation:
## (Intr) Act_Nr
## Act_Neuro 0.00
## Part_Neuro 0.00 0.78
##
## Standardized residuals:
## Min Q1 Med Q3 Max
## -2.40084373 -0.80152859 -0.09535917 0.65075976 2.55129375
##
## Residual standard error: 0.7426419
## Degrees of freedom: 196 total; 193 residual
APIM with SEM
DATA
#SEM data 만들기
write.csv(DATA,'SEM data.csv')
library(dbplyr)
## Warning: package 'dbplyr' was built under R version 3.5.3
##
## Attaching package: 'dbplyr'
## The following objects are masked from 'package:dplyr':
##
## ident, sql
mdata <- subset(DATA, gender==-1)
fdata <- subset(DATA, gender==1)
sem <- cbind(mdata, fdata)
a <- sem[1:6]
b <- sem[8:12]
semdata <- cbind(a,b)
colnames(semdata) <- c("dyad", "male", "MA_distress", "MP_distress", "MA_neuro", "MP_neuro", "female", "FA_distress", "FP_distress", "FA_neuro", "FP_neuro")
SEM <- semdata[,c(1,2,3,5,7,8,10)]
#남녀 neuro, distress를 사용하면 distinguishable dyads 분석이 가능함. 데이터가 남녀 neuro, distress로 구분되어 있지 않음.
#남녀 gender 변수를 기반으로 Actor에 해당하는 neuro, distress만 사용하는 것이 남녀 neuro, distress변수를 사용하는 것이라 생각해서 사용함.
#partner neuro, distress는 사용하지 않음.
#남녀 neuro가 grand center 되어 있음.
mx <- SEM[,4]
fx <- SEM[,7]
X <- rbind(mx, fx)
cm <- mx - mean(X) #남녀 neuro 변수 grand center
cf <- fx - mean(X)
cd <- SEM[, c(1,3,6)]
C_data <- cbind(cd, cm, cf)
colnames(C_data) <- c("dyad","male_distress", "female_distress", "gcenter_maleneuro", "gcenter_femaleneuro")
head(C_data)
## dyad male_distress female_distress gcenter_maleneuro
## 1 1 2.52 3.64 -0.9274554
## 3 2 2.68 3.20 -0.4988554
## 5 3 3.12 3.20 0.2153446
## 7 4 3.76 3.12 -0.6417554
## 9 5 3.08 3.52 1.0725046
## 11 8 3.56 3.96 1.7868046
## gcenter_femaleneuro
## 1 0.7868446
## 3 -0.6417554
## 5 -1.0703554
## 7 -1.0703554
## 9 -0.9274554
## 11 -0.3560554
Distinguishable dyad
library(lavaan)
## Warning: package 'lavaan' was built under R version 3.5.3
## This is lavaan 0.6-5
## lavaan is BETA software! Please report any bugs.
model1 <- '
male_distress ~ a1*gcenter_maleneuro # ACT EFFECT
female_distress ~ a2*gcenter_femaleneuro
male_distress ~ p1*gcenter_femaleneuro # PARTNER EFFECT
female_distress ~ p2*gcenter_maleneuro
gcenter_maleneuro ~ mx1*1 # MEAN X
gcenter_femaleneuro ~ mx2*1
male_distress ~ my1*1 # MEAN Y
female_distress ~ my2*1
gcenter_maleneuro ~~ vx1*gcenter_maleneuro # X VARIANCE
gcenter_femaleneuro ~~ vx2*gcenter_femaleneuro
male_distress ~~ e1*male_distress # Y ERROR
female_distress ~~ e2*female_distress
gcenter_maleneuro ~~ cx*gcenter_femaleneuro # COVARIANCE
male_distress ~~ cy*female_distress
'
APIM_D <- sem(model1, fixed.x = F, data = C_data, missing = "fiml")
summary(APIM_D, fit.measures = T)
## lavaan 0.6-5 ended normally after 29 iterations
##
## Estimator ML
## Optimization method NLMINB
## Number of free parameters 14
##
## Number of observations 98
## Number of missing patterns 1
##
## Model Test User Model:
##
## Test statistic 0.000
## Degrees of freedom 0
##
## Model Test Baseline Model:
##
## Test statistic 119.418
## Degrees of freedom 6
## P-value 0.000
##
## User Model versus Baseline Model:
##
## Comparative Fit Index (CFI) 1.000
## Tucker-Lewis Index (TLI) 1.000
##
## Loglikelihood and Information Criteria:
##
## Loglikelihood user model (H0) -391.512
## Loglikelihood unrestricted model (H1) -391.512
##
## Akaike (AIC) 811.023
## Bayesian (BIC) 847.213
## Sample-size adjusted Bayesian (BIC) 803.003
##
## Root Mean Square Error of Approximation:
##
## RMSEA 0.000
## 90 Percent confidence interval - lower 0.000
## 90 Percent confidence interval - upper 0.000
## P-value RMSEA <= 0.05 NA
##
## Standardized Root Mean Square Residual:
##
## SRMR 0.000
##
## Parameter Estimates:
##
## Information Observed
## Observed information based on Hessian
## Standard errors Standard
##
## Regressions:
## Estimate Std.Err z-value P(>|z|)
## male_distress ~
## gcntr_mln (a1) 0.431 0.105 4.110 0.000
## female_distress ~
## gcntr_fml (a2) 0.385 0.099 3.869 0.000
## male_distress ~
## gcntr_fml (p1) 0.140 0.124 1.130 0.258
## female_distress ~
## gcntr_mln (p2) 0.547 0.084 6.476 0.000
##
## Covariances:
## Estimate Std.Err z-value P(>|z|)
## gcenter_maleneuro ~~
## gcntr_fml (cx) -0.020 0.050 -0.394 0.693
## .male_distress ~~
## .fml_dstrs (cy) 0.342 0.062 5.530 0.000
##
## Intercepts:
## Estimate Std.Err z-value P(>|z|)
## gcntr_ml (mx1) 0.407 0.077 5.266 0.000
## gcntr_fm (mx2) -0.407 0.066 -6.197 0.000
## .ml_dstrs (my1) 3.692 0.103 35.844 0.000
## .fml_dstr (my2) 3.627 0.083 43.759 0.000
##
## Variances:
## Estimate Std.Err z-value P(>|z|)
## gcntr_ml (vx1) 0.585 0.084 7.000 0.000
## gcntr_fm (vx2) 0.422 0.060 7.000 0.000
## .ml_dstrs (e1) 0.630 0.090 7.000 0.000
## .fml_dstr (e2) 0.408 0.058 7.000 0.000
#SEM paths
library(semPlot)
## Warning: package 'semPlot' was built under R version 3.5.3
semPaths(APIM_D,
# general lay-out for a nice APIM:
fade = F, "est", layout='tree2', rotation = 2, style = "ram",
intercepts = F, residuals = F,
optimizeLatRes = T, curve = 3,
# labels and their sizes:
nodeLabels=c("Male_Distress", "Female_Distress", "Male_Neuro", "Female_Neuro"), sizeMan=20, sizeMan2=10,
# position and size of parameter estimates:
edge.label.position = 0.45, edge.label.cex=1.5)
#equal actor effect
model2 <- '
male_distress ~ a*gcenter_maleneuro
female_distress ~ a*gcenter_femaleneuro
male_distress ~ p1*gcenter_femaleneuro
female_distress ~ p2*gcenter_maleneuro
gcenter_maleneuro ~ mx1*1
gcenter_femaleneuro ~ mx2*1
male_distress ~ my1*1
female_distress ~ my2*1
gcenter_maleneuro ~~ vx1*gcenter_maleneuro
gcenter_femaleneuro ~~ vx2*gcenter_femaleneuro
male_distress ~~ e1*male_distress
female_distress ~~ e2*female_distress
gcenter_maleneuro ~~ cx*gcenter_femaleneuro
male_distress ~~ cy*female_distress
'
APIM_EA <- sem(model2, fixed.x = F, data = C_data, missing = "fiml")
summary(APIM_EA, fit.measures = T)
## lavaan 0.6-5 ended normally after 28 iterations
##
## Estimator ML
## Optimization method NLMINB
## Number of free parameters 14
## Number of equality constraints 1
## Row rank of the constraints matrix 1
##
## Number of observations 98
## Number of missing patterns 1
##
## Model Test User Model:
##
## Test statistic 0.108
## Degrees of freedom 1
## P-value (Chi-square) 0.743
##
## Model Test Baseline Model:
##
## Test statistic 119.418
## Degrees of freedom 6
## P-value 0.000
##
## User Model versus Baseline Model:
##
## Comparative Fit Index (CFI) 1.000
## Tucker-Lewis Index (TLI) 1.047
##
## Loglikelihood and Information Criteria:
##
## Loglikelihood user model (H0) -391.565
## Loglikelihood unrestricted model (H1) -391.512
##
## Akaike (AIC) 809.131
## Bayesian (BIC) 842.735
## Sample-size adjusted Bayesian (BIC) 801.683
##
## Root Mean Square Error of Approximation:
##
## RMSEA 0.000
## 90 Percent confidence interval - lower 0.000
## 90 Percent confidence interval - upper 0.186
## P-value RMSEA <= 0.05 0.771
##
## Standardized Root Mean Square Residual:
##
## SRMR 0.010
##
## Parameter Estimates:
##
## Information Observed
## Observed information based on Hessian
## Standard errors Standard
##
## Regressions:
## Estimate Std.Err z-value P(>|z|)
## male_distress ~
## gcntr_mln (a) 0.407 0.073 5.557 0.000
## female_distress ~
## gcntr_fml (a) 0.407 0.073 5.557 0.000
## male_distress ~
## gcntr_fml (p1) 0.157 0.111 1.419 0.156
## female_distress ~
## gcntr_mln (p2) 0.534 0.075 7.150 0.000
##
## Covariances:
## Estimate Std.Err z-value P(>|z|)
## gcenter_maleneuro ~~
## gcntr_fml (cx) -0.020 0.050 -0.394 0.693
## .male_distress ~~
## .fml_dstrs (cy) 0.342 0.062 5.530 0.000
##
## Intercepts:
## Estimate Std.Err z-value P(>|z|)
## gcntr_ml (mx1) 0.407 0.077 5.266 0.000
## gcntr_fm (mx2) -0.407 0.066 -6.197 0.000
## .ml_dstrs (my1) 3.709 0.089 41.907 0.000
## .fml_dstr (my2) 3.641 0.071 51.570 0.000
##
## Variances:
## Estimate Std.Err z-value P(>|z|)
## gcntr_ml (vx1) 0.585 0.084 7.000 0.000
## gcntr_fm (vx2) 0.422 0.060 7.000 0.000
## .ml_dstrs (e1) 0.631 0.090 6.999 0.000
## .fml_dstr (e2) 0.409 0.058 6.998 0.000
anova(APIM_D,APIM_EA)
## Chi-Squared Difference Test
##
## Df AIC BIC Chisq Chisq diff Df diff Pr(>Chisq)
## APIM_D 0 811.02 847.21 0.0000
## APIM_EA 1 809.13 842.74 0.1077 0.10771 1 0.7428
#equal partner effect
model3 <- '
male_distress ~ a1*gcenter_maleneuro
female_distress ~ a2*gcenter_femaleneuro
male_distress ~ p*gcenter_femaleneuro
female_distress ~ p*gcenter_maleneuro
gcenter_maleneuro ~ mx1*1
gcenter_femaleneuro ~ mx2*1
male_distress ~ my1*1
female_distress ~ my2*1
gcenter_maleneuro ~~ vx1*gcenter_maleneuro
gcenter_femaleneuro ~~ vx2*gcenter_femaleneuro
male_distress ~~ e1*male_distress
female_distress ~~ e2*female_distress
gcenter_maleneuro ~~ cx*gcenter_femaleneuro
male_distress ~~ cy*female_distress
'
APIM_EP <- sem(model3, fixed.x = F, data = C_data, missing = "fiml")
summary(APIM_EP, fit.measures = T)
## lavaan 0.6-5 ended normally after 26 iterations
##
## Estimator ML
## Optimization method NLMINB
## Number of free parameters 14
## Number of equality constraints 1
## Row rank of the constraints matrix 1
##
## Number of observations 98
## Number of missing patterns 1
##
## Model Test User Model:
##
## Test statistic 7.386
## Degrees of freedom 1
## P-value (Chi-square) 0.007
##
## Model Test Baseline Model:
##
## Test statistic 119.418
## Degrees of freedom 6
## P-value 0.000
##
## User Model versus Baseline Model:
##
## Comparative Fit Index (CFI) 0.944
## Tucker-Lewis Index (TLI) 0.662
##
## Loglikelihood and Information Criteria:
##
## Loglikelihood user model (H0) -395.205
## Loglikelihood unrestricted model (H1) -391.512
##
## Akaike (AIC) 816.409
## Bayesian (BIC) 850.014
## Sample-size adjusted Bayesian (BIC) 808.962
##
## Root Mean Square Error of Approximation:
##
## RMSEA 0.255
## 90 Percent confidence interval - lower 0.108
## 90 Percent confidence interval - upper 0.441
## P-value RMSEA <= 0.05 0.014
##
## Standardized Root Mean Square Residual:
##
## SRMR 0.084
##
## Parameter Estimates:
##
## Information Observed
## Observed information based on Hessian
## Standard errors Standard
##
## Regressions:
## Estimate Std.Err z-value P(>|z|)
## male_distress ~
## gcntr_mln (a1) 0.328 0.102 3.224 0.001
## female_distress ~
## gcntr_fml (a2) 0.535 0.086 6.239 0.000
## male_distress ~
## gcntr_fml (p) 0.420 0.074 5.702 0.000
## female_distress ~
## gcntr_mln (p) 0.420 0.074 5.702 0.000
##
## Covariances:
## Estimate Std.Err z-value P(>|z|)
## gcenter_maleneuro ~~
## gcntr_fml (cx) -0.020 0.050 -0.394 0.693
## .male_distress ~~
## .fml_dstrs (cy) 0.368 0.067 5.527 0.000
##
## Intercepts:
## Estimate Std.Err z-value P(>|z|)
## gcntr_ml (mx1) 0.407 0.077 5.266 0.000
## gcntr_fm (mx2) -0.407 0.066 -6.197 0.000
## .ml_dstrs (my1) 3.848 0.089 43.255 0.000
## .fml_dstr (my2) 3.740 0.074 50.394 0.000
##
## Variances:
## Estimate Std.Err z-value P(>|z|)
## gcntr_ml (vx1) 0.585 0.084 7.000 0.000
## gcntr_fm (vx2) 0.422 0.060 7.000 0.000
## .ml_dstrs (e1) 0.671 0.097 6.925 0.000
## .fml_dstr (e2) 0.428 0.062 6.875 0.000
anova(APIM_D, APIM_EP)
## Chi-Squared Difference Test
##
## Df AIC BIC Chisq Chisq diff Df diff Pr(>Chisq)
## APIM_D 0 811.02 847.21 0.0000
## APIM_EP 1 816.41 850.01 7.3861 7.3861 1 0.006573 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Indistinguishable dyad
#Data
head(C_data)
## dyad male_distress female_distress gcenter_maleneuro
## 1 1 2.52 3.64 -0.9274554
## 3 2 2.68 3.20 -0.4988554
## 5 3 3.12 3.20 0.2153446
## 7 4 3.76 3.12 -0.6417554
## 9 5 3.08 3.52 1.0725046
## 11 8 3.56 3.96 1.7868046
## gcenter_femaleneuro
## 1 0.7868446
## 3 -0.6417554
## 5 -1.0703554
## 7 -1.0703554
## 9 -0.9274554
## 11 -0.3560554
#indistinguishable dyad
#indistinguishalbe dyad는 6가지 제약을 가함.
#1. equal actor effect,
#2. equal partner effect,
#3. equal X means,
#4. equal X variance,
#5. equal Y intercept,
#6. equal Y error variances
library(lavaan)
model4 <- '
male_distress ~ a*gcenter_maleneuro + p*gcenter_femaleneuro # ACT, PART EFFECT
female_distress ~ a*gcenter_femaleneuro + p*gcenter_maleneuro
gcenter_maleneuro ~ mx*1 # MEAN X
gcenter_femaleneuro ~ mx*1
male_distress ~ my*1 # MEAN Y
female_distress ~ my*1
gcenter_maleneuro ~~ vx*gcenter_maleneuro # X VARIANCE
gcenter_femaleneuro ~~ vx*gcenter_femaleneuro
male_distress ~~ ve*male_distress # Y ERROR
female_distress ~~ ve*female_distress
gcenter_maleneuro ~~ cx*gcenter_femaleneuro # COVARIANCE
male_distress ~~ cy*female_distress
k := p/a # RELATIVE EFFECT SIZE
'
APIM_I <- sem(model4, fixed.x = F, data = C_data, missing = "fiml")
summary(APIM_I, fit.measures = T)
## lavaan 0.6-5 ended normally after 22 iterations
##
## Estimator ML
## Optimization method NLMINB
## Number of free parameters 14
## Number of equality constraints 6
## Row rank of the constraints matrix 6
##
## Number of observations 98
## Number of missing patterns 1
##
## Model Test User Model:
##
## Test statistic 71.100
## Degrees of freedom 6
## P-value (Chi-square) 0.000
##
## Model Test Baseline Model:
##
## Test statistic 119.418
## Degrees of freedom 6
## P-value 0.000
##
## User Model versus Baseline Model:
##
## Comparative Fit Index (CFI) 0.426
## Tucker-Lewis Index (TLI) 0.426
##
## Loglikelihood and Information Criteria:
##
## Loglikelihood user model (H0) -427.061
## Loglikelihood unrestricted model (H1) -391.512
##
## Akaike (AIC) 870.123
## Bayesian (BIC) 890.802
## Sample-size adjusted Bayesian (BIC) 865.540
##
## Root Mean Square Error of Approximation:
##
## RMSEA 0.333
## 90 Percent confidence interval - lower 0.266
## 90 Percent confidence interval - upper 0.404
## P-value RMSEA <= 0.05 0.000
##
## Standardized Root Mean Square Residual:
##
## SRMR 0.304
##
## Parameter Estimates:
##
## Information Observed
## Observed information based on Hessian
## Standard errors Standard
##
## Regressions:
## Estimate Std.Err z-value P(>|z|)
## male_distress ~
## gcntr_mlnr (a) 0.434 0.073 5.983 0.000
## gcntr_fmln (p) 0.355 0.073 4.899 0.000
## female_distress ~
## gcntr_fmln (a) 0.434 0.073 5.983 0.000
## gcntr_mlnr (p) 0.355 0.073 4.899 0.000
##
## Covariances:
## Estimate Std.Err z-value P(>|z|)
## gcenter_maleneuro ~~
## gcntr_fml (cx) -0.185 0.070 -2.642 0.008
## .male_distress ~~
## .fml_dstrs (cy) 0.346 0.065 5.331 0.000
##
## Intercepts:
## Estimate Std.Err z-value P(>|z|)
## gcntr_mln (mx) -0.000 0.050 -0.000 1.000
## gcntr_fml (mx) -0.000 0.050 -0.000 1.000
## .ml_dstrss (my) 3.752 0.067 55.763 0.000
## .fml_dstrs (my) 3.752 0.067 55.763 0.000
##
## Variances:
## Estimate Std.Err z-value P(>|z|)
## gcntr_mln (vx) 0.669 0.070 9.540 0.000
## gcntr_fml (vx) 0.669 0.070 9.540 0.000
## .ml_dstrss (ve) 0.541 0.065 8.341 0.000
## .fml_dstrs (ve) 0.541 0.065 8.341 0.000
##
## Defined Parameters:
## Estimate Std.Err z-value P(>|z|)
## k 0.819 0.105 7.787 0.000
library(semPlot)
semPaths(APIM_I,
# general lay-out for a nice APIM:
fade = F, "est", layout='tree2', rotation = 2, style = "ram",
intercepts = F, residuals = F,
optimizeLatRes = T, curve = 3,
# labels and their sizes:
nodeLabels=c("Male_Distress", "Female_Distress", "Male_Neuro", "Female_Neuro"), sizeMan=20, sizeMan2=10,
# position and size of parameter estimates:
edge.label.position = 0.45, edge.label.cex=1.5)
Housework data
setwd("C:/Users/LG/Documents/Dyadic Analysis Book Data")
library(haven)
house <- read_spss("chapter7 housework.sav")
house[c(1,8,7,4,2,3,5,6)]
## # A tibble: 40 x 8
## Dyad person P_Gender SATISFACTION ACT_HOUSE PART_HOUSE PSATIS
## <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 1 1 1 6 1.2 0.6 7
## 2 1 2 1 7 0.6 1.2 6
## 3 2 1 -1 5 4.3 2.3 9
## 4 2 2 -1 9 2.3 4.3 5
## 5 3 1 1 3 0.4 0.6 4
## 6 3 2 1 4 0.6 0.4 3
## 7 4 1 1 6 0.3 0.5 8
## 8 4 2 1 8 0.5 0.3 6
## 9 5 1 -1 2 3.2 1 6
## 10 5 2 -1 6 1 3.2 2
## # ... with 30 more rows, and 1 more variable: ACT_Gender <dbl>
library(nlme)
housefit <- lme(SATISFACTION ~ ACT_HOUSE + PART_HOUSE ,
random = ~1|Dyad, data = house,
correlation = corCompSymm(form=~1|Dyad))
summary(housefit)
## Linear mixed-effects model fit by REML
## Data: house
## AIC BIC logLik
## 160.3613 170.0269 -74.18067
##
## Random effects:
## Formula: ~1 | Dyad
## (Intercept) Residual
## StdDev: 0.9753118 1.284307
##
## Correlation Structure: Compound symmetry
## Formula: ~1 | Dyad
## Parameter estimate(s):
## Rho
## -0.002142783
## Fixed effects: SATISFACTION ~ ACT_HOUSE + PART_HOUSE
## Value Std.Error DF t-value p-value
## (Intercept) 4.790978 0.6386893 19 7.501266 0.0000
## ACT_HOUSE -0.590827 0.2493903 18 -2.369085 0.0292
## PART_HOUSE 0.887773 0.2493903 18 3.559772 0.0022
## Correlation:
## (Intr) ACT_HO
## ACT_HOUSE -0.615
## PART_HOUSE -0.615 -0.034
##
## Standardized Within-Group Residuals:
## Min Q1 Med Q3 Max
## -1.61969201 -0.59300609 -0.07710746 0.60362287 1.92807281
##
## Number of Observations: 40
## Number of Groups: 20
housefit2 <- lme(SATISFACTION ~ ACT_HOUSE + P_Gender + ACT_HOUSE:P_Gender ,
random = ~1|Dyad, data = house,
correlation = corCompSymm(form=~1|Dyad))
summary(housefit2)
## Linear mixed-effects model fit by REML
## Data: house
## AIC BIC logLik
## 170.9676 182.0522 -78.48381
##
## Random effects:
## Formula: ~1 | Dyad
## (Intercept) Residual
## StdDev: 1.181317 1.429074
##
## Correlation Structure: Compound symmetry
## Formula: ~1 | Dyad
## Parameter estimate(s):
## Rho
## -0.009234356
## Fixed effects: SATISFACTION ~ ACT_HOUSE + P_Gender + ACT_HOUSE:P_Gender
## Value Std.Error DF t-value p-value
## (Intercept) 6.184736 0.6034627 18 10.248746 0.0000
## ACT_HOUSE -0.438573 0.3208041 18 -1.367106 0.1884
## P_Gender -0.736880 0.6034627 18 -1.221086 0.2378
## ACT_HOUSE:P_Gender 0.481140 0.3208041 18 1.499793 0.1510
## Correlation:
## (Intr) ACT_HOUSE P_Gndr
## ACT_HOUSE -0.793
## P_Gender -0.100 -0.081
## ACT_HOUSE:P_Gender -0.081 0.342 -0.793
##
## Standardized Within-Group Residuals:
## Min Q1 Med Q3 Max
## -1.7504096 -0.6419040 -0.1819624 0.6182831 1.7572948
##
## Number of Observations: 40
## Number of Groups: 20
housefit3 <- lme(SATISFACTION ~ ACT_HOUSE + PART_HOUSE + ACT_HOUSE:ACT_Gender + PART_HOUSE:P_Gender ,
random = ~1|Dyad, data = house,
correlation = corCompSymm(form=~1|Dyad))
summary(housefit3)
## Linear mixed-effects model fit by REML
## Data: house
## AIC BIC logLik
## 164.3323 176.7751 -74.16617
##
## Random effects:
## Formula: ~1 | Dyad
## (Intercept) Residual
## StdDev: 0.8993966 1.297676
##
## Correlation Structure: Compound symmetry
## Formula: ~1 | Dyad
## Parameter estimate(s):
## Rho
## -0.00506752
## Fixed effects: SATISFACTION ~ ACT_HOUSE + PART_HOUSE + ACT_HOUSE:ACT_Gender + PART_HOUSE:P_Gender
## Value Std.Error DF t-value p-value
## (Intercept) 4.442690 0.6551773 19 6.780898 0.0000
## ACT_HOUSE -0.372207 0.2824928 16 -1.317582 0.2062
## PART_HOUSE 0.952839 0.2824928 16 3.372969 0.0039
## ACT_HOUSE:ACT_Gender 0.296500 0.2276240 16 1.302586 0.2112
## PART_HOUSE:P_Gender -0.014722 0.2276240 16 -0.064677 0.9492
## Correlation:
## (Intr) ACT_HOUSE PART_HOUSE ACT_HOUSE:
## ACT_HOUSE -0.602
## PART_HOUSE -0.602 -0.091
## ACT_HOUSE:ACT_Gender -0.137 0.464 -0.204
## PART_HOUSE:P_Gender -0.137 -0.204 0.464 -0.680
##
## Standardized Within-Group Residuals:
## Min Q1 Med Q3 Max
## -1.71601669 -0.55916876 -0.02829742 0.59927155 1.90264555
##
## Number of Observations: 40
## Number of Groups: 20