pacman::p_load(dplyr, rio,ggplot2,psych, nlme,tidyr,sjPlot,osfr,GGally) Dyadic Analysis workshop, Zurich Nov. 24: APIM longitudinal
APIM longitudinal
Libraries
Uploading data from OSF
file <- osf_retrieve_file("6718dd0a2295d821b91b3a1f")
setwd("C:\\Users\\97254\\Downloads")
osf_download(file, progress = T,conflicts = "overwrigt")
temp1 <- import("RSAData.csv")
temp1 <- temp1 %>% select(Couple,Id,DiaryDay,Female,NegMood,PPR)Creating dyadic data for MLM
Male <- temp1 %>% filter(Female==0)
pMale <- temp1 %>% filter(Female==1) %>% rename(pNegMood=NegMood,pPPR=PPR) %>% select(-Female,-Id)
Male <- left_join(Male,pMale,by=c("Couple","DiaryDay"))
Female <- temp1 %>% filter(Female==1)
pFemale <- temp1 %>% filter(Female==0) %>% rename(pNegMood=NegMood,pPPR=PPR) %>% select(-Female,-Id)
Female <- left_join(Female,pFemale,by=c("Couple","DiaryDay"))
temp2 <- bind_rows(Male,Female) %>% arrange(Couple,DiaryDay)
head(temp2) Couple Id DiaryDay Female NegMood PPR pNegMood pPPR
1 200 3000 0 0 32.50 5.333333 27.00 5.333333
2 200 4000 0 1 27.00 5.333333 32.50 5.333333
3 200 4000 1 1 14.50 5.000000 NA NA
4 200 3000 2 0 10.75 5.000000 16.75 5.000000
5 200 4000 2 1 16.75 5.000000 10.75 5.000000
6 200 3000 3 0 26.50 5.666667 15.00 5.000000Plots
temp3 <- temp2 %>% mutate(Male=ifelse(Female==0,1,0),
Role=ifelse(Female==0,"Male","Female"))
#plots
ggpairs(temp3 %>% select(PPR,pPPR), title="Corr Matrix")
ggplot(filter(temp3,Couple==200|Couple==203), aes(x=NegMood, y=PPR, color=Role, shape=Role)) +
geom_point() +
facet_wrap(~ Couple,nrow=2)+
geom_smooth(method=lm, se=FALSE, fullrange=TRUE)+
theme_classic()+scale_color_brewer(palette="Dark2")
ggplot(temp3, aes(x=NegMood, y=PPR, color=Role, shape=Role)) +
geom_point() +
facet_wrap(~ Couple)+
geom_smooth(method=lm, se=FALSE, fullrange=TRUE)+
theme_classic()+scale_color_brewer(palette="Dark2")
ggplot(temp3, aes(x = NegMood, y = PPR, color = Role, group = interaction(Couple, Role))) +
# geom_line() +
geom_smooth(method = lm, se = FALSE, fullrange = TRUE, aes(group = interaction(Couple, Role)), color = "gray", size = 0.5) + # Individual regression lines
geom_smooth(aes(group = 1), method = lm, se = FALSE, fullrange = TRUE, color = "blue", size = 1.2) + # Average regression line
theme_classic() +
scale_color_brewer(palette = "Dark2")
Person-mean centering
#person-mean centering
temp3<-temp3 %>%
group_by(Couple,Male) %>%
mutate(NegMoodmc=NegMood-mean(NegMood, na.rm = T),AvgNegMood=mean(NegMood,na.rm=T),
pNegMoodmc=pNegMood-mean(pNegMood, na.rm = T),pAvgNegMood=mean(pNegMood,na.rm=T)
) %>%
ungroup()
head(temp3)# A tibble: 6 × 14
Couple Id DiaryDay Female NegMood PPR pNegMood pPPR Male Role
<int> <int> <int> <int> <dbl> <dbl> <dbl> <dbl> <dbl> <chr>
1 200 3000 0 0 32.5 5.33 27 5.33 1 Male
2 200 4000 0 1 27 5.33 32.5 5.33 0 Female
3 200 4000 1 1 14.5 5 NA NA 0 Female
4 200 3000 2 0 10.8 5 16.8 5 1 Male
5 200 4000 2 1 16.8 5 10.8 5 0 Female
6 200 3000 3 0 26.5 5.67 15 5 1 Male
# ℹ 4 more variables: NegMoodmc <dbl>, AvgNegMood <dbl>, pNegMoodmc <dbl>,
# pAvgNegMood <dbl>Estimating Multi-Level APIM
model1 <- lme(fixed = PPR ~ -1 +
Male + #intercept for males
Male:NegMoodmc + #Actor effect for males
Male:pNegMoodmc + #Partner effect for males
Female + #intercept for females
Female:NegMoodmc + #Actor effect for females
Female:pNegMoodmc, #Partner effect for females
random = ~ -1 + #random effects between-dyads, G matrix
Male+Male:NegMoodmc+Male:pNegMoodmc+
Female+Female:NegMoodmc+Female:pNegMoodmc| Couple,
weights=varIdent(form = ~1 | Male), # this invokes separate sigma^{2}_{e} for each member
corr=corCompSymm(form = ~1 | Couple/DiaryDay), # this invokes the off-diaginal sigma_{e1e2}
data=temp3,na.action = na.exclude,
control =list(msMaxIter = 1000, msMaxEval = 1000))
summary(model1)Linear mixed-effects model fit by REML
Data: temp3
AIC BIC logLik
5806.811 5984.315 -2873.405
Random effects:
Formula: ~-1 + Male + Male:NegMoodmc + Male:pNegMoodmc + Female + Female:NegMoodmc + Female:pNegMoodmc | Couple
Structure: General positive-definite, Log-Cholesky parametrization
StdDev Corr
Male 0.81741666 Male Female Ml:NgM Ml:pNM NgMd:F
Female 0.77405038 0.540
Male:NegMoodmc 0.01864470 0.392 0.614
Male:pNegMoodmc 0.00902635 0.388 0.212 0.211
NegMoodmc:Female 0.01673060 0.341 0.593 0.640 0.064
pNegMoodmc:Female 0.01400373 0.189 0.393 0.545 -0.515 0.529
Residual 0.63324832
Correlation Structure: Compound symmetry
Formula: ~1 | Couple/DiaryDay
Parameter estimate(s):
Rho
0.1451731
Variance function:
Structure: Different standard deviations per stratum
Formula: ~1 | Male
Parameter estimates:
1 0
1.0000000 0.9124022
Fixed effects: PPR ~ -1 + Male + Male:NegMoodmc + Male:pNegMoodmc + Female + Female:NegMoodmc + Female:pNegMoodmc
Value Std.Error DF t-value p-value
Male 5.870359 0.09196289 2662 63.83400 0.0000
Female 5.983040 0.08696509 2662 68.79818 0.0000
Male:NegMoodmc -0.015298 0.00285172 2662 -5.36442 0.0000
Male:pNegMoodmc -0.006307 0.00192096 2662 -3.28303 0.0010
NegMoodmc:Female -0.014925 0.00243029 2662 -6.14123 0.0000
pNegMoodmc:Female -0.006674 0.00236988 2662 -2.81634 0.0049
Correlation:
Male Female Ml:NgM Ml:pNM NgMd:F
Female 0.526
Male:NegMoodmc 0.279 0.436
Male:pNegMoodmc 0.198 0.108 -0.027
NegMoodmc:Female 0.254 0.445 0.343 0.100
pNegMoodmc:Female 0.122 0.252 0.352 -0.209 0.189
Standardized Within-Group Residuals:
Min Q1 Med Q3 Max
-7.68311135 -0.35994239 0.08015739 0.48011754 4.12884461
Number of Observations: 2749
Number of Groups: 82 tab_model(model1)| PPR | |||
|---|---|---|---|
| Predictors | Estimates | CI | p |
| Male | 5.87 | 5.69 – 6.05 | <0.001 |
| Female | 5.98 | 5.81 – 6.15 | <0.001 |
| Male × NegMoodmc | -0.02 | -0.02 – -0.01 | <0.001 |
| Male × pNegMoodmc | -0.01 | -0.01 – -0.00 | 0.001 |
| NegMoodmc × Female | -0.01 | -0.02 – -0.01 | <0.001 |
| pNegMoodmc × Female | -0.01 | -0.01 – -0.00 | 0.005 |
| Random Effects | |||
| σ2 | 0.40 | ||
| τ00 | |||
| τ00 | |||
| τ11 Couple.Female | 0.60 | ||
| τ11 Couple.Male:NegMoodmc | 0.00 | ||
| τ11 Couple.Male:pNegMoodmc | 0.00 | ||
| τ11 Couple.NegMoodmc:Female | 0.00 | ||
| τ11 Couple.pNegMoodmc:Female | 0.00 | ||
| ρ01 | 0.54 | ||
| 0.39 | |||
| 0.39 | |||
| 0.34 | |||
| 0.19 | |||
| ICC | 0.63 | ||
| N Couple | 82 | ||
| Observations | 2749 | ||
| Marginal R2 / Conditional R2 | 0.036 / 0.646 | ||
The model was too complex and hardly converged, here I ran a simpler G-matrix model (only co/var in intercepts)
model2 <- lme(fixed = PPR ~ -1 +
Male + #intercept for males
Male:NegMoodmc + #Actor effect for males
Male:pNegMoodmc + #Partner effect for males
Female + #intercept for females
Female:NegMoodmc + #Actor effect for females
Female:pNegMoodmc, #Partner effect for females
random = ~ -1 + #random effects between-dyads, G matrix
Male+Female| Couple,
weights=varIdent(form = ~1 | Male), # this invokes separate sigma^{2}_{e} for each member
corr=corCompSymm(form = ~1 | Couple/DiaryDay), # this invokes the off-diaginal sigma_{e1e2}
data=temp3,na.action = na.exclude,
control=list(maxIter=1000))
summary(model2)Linear mixed-effects model fit by REML
Data: temp3
AIC BIC logLik
5961.723 6032.725 -2968.862
Random effects:
Formula: ~-1 + Male + Female | Couple
Structure: General positive-definite, Log-Cholesky parametrization
StdDev Corr
Male 0.8143544 Male
Female 0.7747095 0.542
Residual 0.6725580
Correlation Structure: Compound symmetry
Formula: ~1 | Couple/DiaryDay
Parameter estimate(s):
Rho
0.1734369
Variance function:
Structure: Different standard deviations per stratum
Formula: ~1 | Male
Parameter estimates:
1 0
1.0000000 0.9325946
Fixed effects: PPR ~ -1 + Male + Male:NegMoodmc + Male:pNegMoodmc + Female + Female:NegMoodmc + Female:pNegMoodmc
Value Std.Error DF t-value p-value
Male 5.871260 0.09183373 2662 63.93360 0e+00
Female 5.980983 0.08728935 2662 68.51905 0e+00
Male:NegMoodmc -0.016254 0.00180741 2662 -8.99276 0e+00
Male:pNegMoodmc -0.007850 0.00162768 2662 -4.82274 0e+00
NegMoodmc:Female -0.015908 0.00151637 2662 -10.49100 0e+00
pNegMoodmc:Female -0.006086 0.00168553 2662 -3.61049 3e-04
Correlation:
Male Female Ml:NgM Ml:pNM NgMd:F
Female 0.527
Male:NegMoodmc 0.001 0.000
Male:pNegMoodmc 0.000 0.000 -0.226
NegMoodmc:Female 0.000 0.001 -0.039 0.173
pNegMoodmc:Female 0.000 0.000 0.173 -0.039 -0.225
Standardized Within-Group Residuals:
Min Q1 Med Q3 Max
-7.5365903 -0.3690256 0.0845538 0.5099811 3.6410083
Number of Observations: 2749
Number of Groups: 82 This is a more complex G model - var of all variables, covariance between same predictors
model3 <- lme(fixed = PPR ~ -1 +
Male + #intercept for males
Male:NegMoodmc + #Actor effect for males
Male:pNegMoodmc + #Partner effect for males
Female + #intercept for females
Female:NegMoodmc + #Actor effect for females
Female:pNegMoodmc, #Partner effect for females
#radnom model G matrix var of all model, covariance only btw same predictors client-therapist,actor-actor,partner-partner
random = list(Couple=pdBlocked(list(~0+Male+Female,~0+Male:NegMoodmc+Female:NegMoodmc,
~0+Male:pNegMoodmc+Female:pNegMoodmc))),
weights=varIdent(form = ~1 | Male), # this invokes separate sigma^{2}_{e} for each member
corr=corCompSymm(form = ~1 | Couple/DiaryDay), # this invokes the off-diaginal sigma_{e1e2}
data=temp3,na.action = na.exclude,
control=list(maxIter=1000))
summary(model3)Linear mixed-effects model fit by REML
Data: temp3
AIC BIC logLik
5823.41 5929.912 -2893.705
Random effects:
Composite Structure: Blocked
Block 1: Male, Female
Formula: ~0 + Male + Female | Couple
Structure: General positive-definite
StdDev Corr
Male 0.8176301 Male
Female 0.7744229 0.539
Block 2: Male:NegMoodmc, NegMoodmc:Female
Formula: ~0 + Male:NegMoodmc + Female:NegMoodmc | Couple
Structure: General positive-definite
StdDev Corr
Male:NegMoodmc 0.01804104 Ml:NgM
NegMoodmc:Female 0.01690776 0.623
Block 3: Male:pNegMoodmc, pNegMoodmc:Female
Formula: ~0 + Male:pNegMoodmc + Female:pNegMoodmc | Couple
Structure: General positive-definite
StdDev Corr
Male:pNegMoodmc 0.01041809 Ml:pNM
pNegMoodmc:Female 0.01323298 -0.286
Residual 0.63172769
Correlation Structure: Compound symmetry
Formula: ~1 | Couple/DiaryDay
Parameter estimate(s):
Rho
0.1534077
Variance function:
Structure: Different standard deviations per stratum
Formula: ~1 | Male
Parameter estimates:
1 0
1.0000000 0.9163729
Fixed effects: PPR ~ -1 + Male + Male:NegMoodmc + Male:pNegMoodmc + Female + Female:NegMoodmc + Female:pNegMoodmc
Value Std.Error DF t-value p-value
Male 5.869981 0.09197961 2662 63.81828 0.0000
Female 5.982845 0.08701443 2662 68.75693 0.0000
Male:NegMoodmc -0.014339 0.00283696 2662 -5.05426 0.0000
Male:pNegMoodmc -0.006372 0.00203782 2662 -3.12692 0.0018
NegMoodmc:Female -0.015010 0.00246554 2662 -6.08786 0.0000
pNegMoodmc:Female -0.006678 0.00234633 2662 -2.84621 0.0045
Correlation:
Male Female Ml:NgM Ml:pNM NgMd:F
Female 0.525
Male:NegMoodmc 0.002 0.000
Male:pNegMoodmc 0.000 0.000 -0.106
NegMoodmc:Female 0.000 0.001 0.335 0.075
pNegMoodmc:Female 0.000 -0.001 0.071 -0.117 -0.086
Standardized Within-Group Residuals:
Min Q1 Med Q3 Max
-7.60174269 -0.35495152 0.08945415 0.47199527 3.95920792
Number of Observations: 2749
Number of Groups: 82 Now adding person-level predictors
# Grand mean-centering
temp3 <- temp3 %>% mutate(AvgNegMoodmc=AvgNegMood-mean(AvgNegMood,na.rm=T),
pAvgNegMoodmc=pAvgNegMood-mean(AvgNegMood,na.rm=T)
)
model4 <- lme(fixed = PPR ~ -1 +
Male + #intercept for males
Male:NegMoodmc + #Actor effect for males
Male:pNegMoodmc + #Partner effect for males
Male:AvgNegMoodmc+#Actor effect for males- between-subject
Male:pAvgNegMoodmc+#Partner effect for males- between subject
Female + #intercept for females
Female:NegMoodmc + #Actor effect for females
Female:pNegMoodmc+ #Partner effect for females
Female:AvgNegMoodmc+#Actor effect for females- between-subject
Female:pAvgNegMoodmc,#Partner effect for females- between subject
#radnom model G matrix var of all model, covariance only btw same predictors client-therapist,actor-actor,partner-partner
random = list(Couple=pdBlocked(list(~0+Male+Female,~0+Male:NegMoodmc+Female:NegMoodmc,
~0+Male:pNegMoodmc+Female:pNegMoodmc))),
weights=varIdent(form = ~1 | Male), # this invokes separate sigma^{2}_{e} for each member
corr=corCompSymm(form = ~1 | Couple/DiaryDay), # this invokes the off-diaginal sigma_{e1e2}
data=temp3,na.action = na.exclude,
control=list(maxIter=1000))
summary(model4)Linear mixed-effects model fit by REML
Data: temp3
AIC BIC logLik
5854.045 5984.183 -2905.022
Random effects:
Composite Structure: Blocked
Block 1: Male, Female
Formula: ~0 + Male + Female | Couple
Structure: General positive-definite
StdDev Corr
Male 0.8128935 Male
Female 0.7557917 0.536
Block 2: Male:NegMoodmc, NegMoodmc:Female
Formula: ~0 + Male:NegMoodmc + Female:NegMoodmc | Couple
Structure: General positive-definite
StdDev Corr
Male:NegMoodmc 0.01804175 Ml:NgM
NegMoodmc:Female 0.01689566 0.623
Block 3: Male:pNegMoodmc, pNegMoodmc:Female
Formula: ~0 + Male:pNegMoodmc + Female:pNegMoodmc | Couple
Structure: General positive-definite
StdDev Corr
Male:pNegMoodmc 0.01041759 Ml:pNM
pNegMoodmc:Female 0.01323565 -0.286
Residual 0.63171333
Correlation Structure: Compound symmetry
Formula: ~1 | Couple/DiaryDay
Parameter estimate(s):
Rho
0.1534125
Variance function:
Structure: Different standard deviations per stratum
Formula: ~1 | Male
Parameter estimates:
1 0
1.0000000 0.9164085
Fixed effects: PPR ~ -1 + Male + Male:NegMoodmc + Male:pNegMoodmc + Male:AvgNegMoodmc + Male:pAvgNegMoodmc + Female + Female:NegMoodmc + Female:pNegMoodmc + Female:AvgNegMoodmc + Female:pAvgNegMoodmc
Value Std.Error DF t-value p-value
Male 5.890341 0.09279129 2658 63.47946 0.0000
Female 6.021978 0.08640453 2658 69.69516 0.0000
Male:NegMoodmc -0.014334 0.00283705 2658 -5.05232 0.0000
Male:pNegMoodmc -0.006374 0.00203776 2658 -3.12782 0.0018
Male:AvgNegMoodmc 0.002594 0.00867655 2658 0.29897 0.7650
Male:pAvgNegMoodmc -0.017506 0.01025574 2658 -1.70697 0.0879
NegMoodmc:Female -0.015013 0.00246446 2658 -6.09198 0.0000
pNegMoodmc:Female -0.006678 0.00234656 2658 -2.84582 0.0045
AvgNegMoodmc:Female -0.018762 0.00961748 2658 -1.95087 0.0512
pAvgNegMoodmc:Female 0.016981 0.00809195 2658 2.09855 0.0360
Correlation:
Male Female Ml:NgM Ml:pNM Ml:AvNM Ml:pANM NgMd:F pNgM:F
Female 0.522
Male:NegMoodmc 0.002 0.000
Male:pNegMoodmc 0.000 0.000 -0.106
Male:AvgNegMoodmc 0.130 0.071 0.002 0.000
Male:pAvgNegMoodmc -0.138 -0.077 -0.001 0.000 -0.273
NegMoodmc:Female 0.000 0.001 0.334 0.075 0.000 0.001
pNegMoodmc:Female 0.000 -0.001 0.071 -0.117 0.000 0.000 -0.086
AvgNegMoodmc:Female -0.072 -0.150 -0.001 0.000 -0.143 0.518 0.000 0.000
pAvgNegMoodmc:Female 0.068 0.138 0.001 0.000 0.519 -0.145 0.000 0.000
AvNM:F
Female
Male:NegMoodmc
Male:pNegMoodmc
Male:AvgNegMoodmc
Male:pAvgNegMoodmc
NegMoodmc:Female
pNegMoodmc:Female
AvgNegMoodmc:Female
pAvgNegMoodmc:Female -0.285
Standardized Within-Group Residuals:
Min Q1 Med Q3 Max
-7.63522814 -0.35672458 0.08504122 0.47670683 3.97116322
Number of Observations: 2749
Number of Groups: 82 Testing for gender differences in fixed effects
model4 <- lme(fixed = PPR ~
1 +
Male*NegMoodmc +
Male*pNegMoodmc
,
random = ~ -1 + #random effects between-dyads, G matrix
Male+Female| Couple,
weights=varIdent(form = ~1 | Male), # this invokes separate sigma^{2}_{e} for each member
corr=corCompSymm(form = ~1 | Couple/DiaryDay), # this invokes the off-diaginal sigma_{e1e2}
data=temp3,na.action = na.exclude,
control=list(maxIter=1000))
summary(model4)Linear mixed-effects model fit by REML
Data: temp3
AIC BIC logLik
5961.723 6032.725 -2968.862
Random effects:
Formula: ~-1 + Male + Female | Couple
Structure: General positive-definite, Log-Cholesky parametrization
StdDev Corr
Male 0.8143543 Male
Female 0.7747095 0.542
Residual 0.6725581
Correlation Structure: Compound symmetry
Formula: ~1 | Couple/DiaryDay
Parameter estimate(s):
Rho
0.173437
Variance function:
Structure: Different standard deviations per stratum
Formula: ~1 | Male
Parameter estimates:
1 0
1.0000000 0.9325946
Fixed effects: PPR ~ 1 + Male * NegMoodmc + Male * pNegMoodmc
Value Std.Error DF t-value p-value
(Intercept) 5.980983 0.08728935 2662 68.51905 0.0000
Male -0.109722 0.08718047 2662 -1.25857 0.2083
NegMoodmc -0.015908 0.00151637 2662 -10.49100 0.0000
pNegMoodmc -0.006086 0.00168553 2662 -3.61049 0.0003
Male:NegMoodmc -0.000345 0.00240415 2662 -0.14365 0.8858
Male:pNegMoodmc -0.001764 0.00238836 2662 -0.73869 0.4602
Correlation:
(Intr) Male NgMdmc pNgMdm Ml:NgM
Male -0.446
NegMoodmc 0.001 -0.001
pNegMoodmc 0.000 0.001 -0.225
Male:NegMoodmc -0.001 0.001 -0.660 0.272
Male:pNegMoodmc 0.000 -0.001 0.277 -0.732 -0.382
Standardized Within-Group Residuals:
Min Q1 Med Q3 Max
-7.53659026 -0.36902563 0.08455381 0.50998115 3.64100826
Number of Observations: 2749
Number of Groups: 82 tab_model(model4)| PPR | |||
|---|---|---|---|
| Predictors | Estimates | CI | p |
| (Intercept) | 5.98 | 5.81 – 6.15 | <0.001 |
| Male | -0.11 | -0.28 – 0.06 | 0.208 |
| NegMoodmc | -0.02 | -0.02 – -0.01 | <0.001 |
| pNegMoodmc | -0.01 | -0.01 – -0.00 | <0.001 |
| Male × NegMoodmc | -0.00 | -0.01 – 0.00 | 0.886 |
| Male × pNegMoodmc | -0.00 | -0.01 – 0.00 | 0.460 |
| Random Effects | |||
| σ2 | 0.45 | ||
| τ00 | |||
| τ00 | |||
| τ11 Couple.Female | 0.60 | ||
| ρ01 Couple | 0.54 | ||
| ICC | 0.42 | ||
| N Couple | 82 | ||
| Observations | 2749 | ||
| Marginal R2 / Conditional R2 | 0.055 / 0.455 | ||