library(dplyr)
library(haven)
library(psych)
library(PerformanceAnalytics)
library(corrplot)
library(Hmisc)
library(lavaan)
library(semTools)
library(mvnTest)
library(mvnormalTest)

###Loading matrix###

CRS_dec <- read_sav("CRS_dec.sav") %>%
  filter(COUNTRY == 3) %>%
  select(starts_with("CRS_"))

  ###Reliability analysis###
  #Tarkkonen's rho syntax
  rho <- function(scale, nfactors)  {
    scale = na.omit(scale)
    res = psych::fa(scale, nfactors, rotate="varimax")
    ### Tarkkonen Rho ###
    C = cor(scale)
    B = res$loading
    A = res$weights
    phi = 1/(res$n.obs-1) * (t(res$scores) %*% res$scores)
    d = dim(res$loading)[2]
    
    I1 = diag(d)
    diag(I1) = diag(t(A) %*% B %*% phi %*% t(B) %*% A)
    I2 = diag(d)
    diag(I2) = diag(t(A) %*% C %*% A)
    Rho = diag(I1 %*% solve(I2))
    print(Rho) }
  
  #Coparenting Agreement
  agree <- CRS_dec %>% dplyr::select(CRS_6, CRS_9r, CRS_11r, CRS_15r)
  ###3: Realiability measures
  psych::alpha(agree) #Cronbach’s alpha coefficient
  rho(agree, 1) #Tarkonnen' rho
 
  #Coparenting Closeness
  close <- CRS_dec %>% dplyr::select(CRS_2, CRS_17, CRS_24, CRS_28r, CRS_30)
  ###3: Realiability measures
  psych::alpha(close) #Cronbach’s alpha coefficient
  rho(close, 1) #Tarkonnen' rho
  
  #Exposure to Conflict
  exp <- CRS_dec %>% dplyr::select(CRS_31, CRS_32, CRS_33, CRS_34, CRS_35)
  ###3: Realiability measures
  psych::alpha(exp) #Cronbach’s alpha coefficient
  rho(exp, 1) #Tarkonnen' rho
  
  #Coparenting Support
  supp <- CRS_dec %>% dplyr::select(CRS_3, CRS_10, CRS_19, CRS_25, CRS_26, CRS_27)
  ###3: Realiability measures
  psych::alpha(supp) #Cronbach’s alpha coefficient
  rho(supp, 1) #Tarkonnen' rho
  
  #Coparenting Undermining
  und <- CRS_dec %>% dplyr::select(CRS_8, CRS_12, CRS_13, CRS_16, CRS_21, CRS_22)
  ###3: Realiability measures
  psych::alpha(und) #Cronbach’s alpha coefficient
  rho(und, 1) #Tarkonnen' rho
  
  #Endorsement Partner’s Parenting
  par <- CRS_dec %>% dplyr::select(CRS_1, CRS_4, CRS_7r, CRS_14, CRS_18, CRS_23, CRS_29r)
  ###3: Realiability measures
  psych::alpha(par) #Cronbach’s alpha coefficient
  rho(par, 1) #Tarkonnen' rho
  
  #Division of Labor
  div <- CRS_dec %>% dplyr::select(CRS_20r, CRS_5r)
  ###3: Realiability measures
  psych::alpha(div) #Cronbach’s alpha coefficient
  rho(div, 1) #Tarkonnen' rho


  ###For fast omega calculation use this model as an input:
  Reliability_CRS <- 'agree =~ CRS_6 + CRS_9r + CRS_11r + CRS_15r
          close =~ CRS_2 + CRS_17 + CRS_24 + CRS_28r + CRS_30
          exp =~ CRS_31 + CRS_32 + CRS_33 + CRS_34 + CRS_35
          supp =~ CRS_3 + CRS_10 + CRS_19 + CRS_25 + CRS_26 + CRS_27
          und =~ CRS_8 + CRS_12 + CRS_13 + CRS_16 + CRS_21 + CRS_22
          par =~ CRS_1 + CRS_4 + CRS_7r + CRS_14 + CRS_18 + CRS_23 + CRS_29r
          div =~ CRS_5r + CRS_20r'
  
  res.Reliability_CRS <- lavaan::cfa(Reliability_CRS, CRS_dec, estimator = "MLR", mimic="MPLUS")
  summary(res.Reliability_CRS, fit = TRUE, standardized = TRUE)
  
  reliability(res.Reliability_CRS)
  
  #CRS short version
  CRS_Brief <- CRS_dec %>% dplyr::select(CRS_1, CRS_2, CRS_4, CRS_5r, CRS_6, CRS_9r, CRS_16r, CRS_20r, CRS_22r, CRS_24, CRS_25, CRS_27, CRS_33r, CRS_34r)
  ###3: Realiability measures
  psych::alpha(CRS_Brief) #Cronbach’s alpha coefficient
  psych::omega(CRS_Brief) #McDonald’s omega (ω)
  rho(CRS_Brief, 1) #Tarkkonen’s rho (ρ)
  
  #CRS full version
  CRS_Full <- CRS_dec %>% dplyr::select(CRS_6, CRS_9r, CRS_11r, CRS_15r, CRS_2, CRS_17, CRS_24, CRS_28r, CRS_30, CRS_31r, CRS_32r, CRS_33r, CRS_34r, CRS_35r, CRS_3, CRS_10, CRS_19, CRS_25, CRS_26, CRS_27, CRS_8r, CRS_12r, CRS_13r, CRS_16r, CRS_21r, CRS_22r, CRS_1, CRS_4, CRS_7r, CRS_14, CRS_18, CRS_23, CRS_29r, CRS_5r, CRS_20r)
  ###3: Realiability measures
  psych::alpha(CRS_Full) #Cronbach’s alpha coefficient
  psych::omega(CRS_Full) #McDonald’s omega (ω)
  rho(CRS_Full, 1) #Tarkonnen' rho
  
  #CRS short version without CRS_5 and CRS_20
  CRS_Brief <- CRS_dec %>% dplyr::select(CRS_1, CRS_2, CRS_4, CRS_6, CRS_9r, CRS_16r, CRS_22r, CRS_24, CRS_25, CRS_27, CRS_33r, CRS_34r)
  ###3: Realiability measures
  psych::alpha(CRS_Brief) #Cronbach’s alpha coefficient
  psych::omega(CRS_Brief) #McDonald’s omega (ω)
  rho(CRS_Brief, 1) #Tarkkonen’s rho (ρ)
  
  #CRS full version without CRS_5 and CRS_20
  CRS_Full <- CRS_dec %>% dplyr::select(CRS_6, CRS_9r, CRS_11r, CRS_15r, CRS_2, CRS_17, CRS_24, CRS_28r, CRS_30, CRS_31r, CRS_32r, CRS_33r, CRS_34r, CRS_35r, CRS_3, CRS_10, CRS_19, CRS_25, CRS_26, CRS_27, CRS_8r, CRS_12r, CRS_13r, CRS_16r, CRS_21r, CRS_22r, CRS_1, CRS_4, CRS_7r, CRS_14, CRS_18, CRS_23, CRS_29r)
  ###3: Realiability measures
  psych::alpha(CRS_Full) #Cronbach’s alpha coefficient
  psych::omega(CRS_Full) #McDonald’s omega (ω)
  rho(CRS_Full, 1) #Tarkonnen' rho