---
title: "R Codes"
author: "Bianca Silva"
date: "29 de abril de 2022"
output:
  word_document: default
  html_document: default
---

```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = TRUE)
```
##Libraries, Database, Items and Models
```{r}
# Libraries/Packages
library(lavaan)
library(semTools)
library(foreign)
library(haven)
library(readxl)

## read database



#Amostras teste e validacao

bd <- read_sav ("C:/Users/samsung/Dropbox/Bianca Odonto/MESTRADO/Dissertação/Artigo B/BD_trabalho.sav")

## Define items analysis

it1 = c("ECOHIS1", "ECOHIS2", "ECOHIS3", "ECOHIS4", "ECOHIS5", "ECOHIS6", "ECOHIS7" , "ECOHIS8", "ECOHIS9" , "ECOHIS10", "ECOHIS11", "ECOHIS12", "ECOHIS13")


## Specify Factorial Model to CFA - Confirmatoria
#Modelos de medida
modelo1 = '
ECOHIS =~ ECOHIS1 + ECOHIS2 + ECOHIS3 + ECOHIS4 + ECOHIS5 + ECOHIS6 + ECOHIS7 + ECOHIS8 + ECOHIS9 + ECOHIS10 + ECOHIS11 + ECOHIS12 + ECOHIS13
'

```

#Sensibilidade psicometrica
```{r sensibilidade, echo=F, message=F, warning=F, results="asis"}
skim_with_defaults()
fun <- list(mean=mean,
            median=median,
            sd=sd,
            sk = function(x) round(psych::skew(x, na.rm=TRUE), 2),
            ku = function(x) round(psych::kurtosi(x, na.rm=TRUE),2))
#skim_with(numeric = list(hist = NULL))
skim_with(numeric = fun, append = T)
#detach("package:knitr", unload = T)
bd[,it2t] %>%
  skim() %>%  
  kable(digits=2, caption = "Descriptives with histogram")
```

##Confirmatory Factor Analysis
```{r}
#fit the model CFA
fit1 = cfa(modelo1, data=bd,std.lv=T, ordered = it1, estimator = "WLSMV")

#factor loadings
summary(fit1, fit.measures=TRUE, standardized=TRUE)
fitmeasures(fit1, fit.measures = "all", baseline.model = NULL)





```

##Modification indexes (MI)
```{r}
#MI
im1 = modindices(fit4, sort.=TRUE, minimum.value=4)


#MI for fator loadings
#selIM=im1[im1$op == "~~",]
#knitr::kable(selIM)

selIM1=im1[im1$op == "~~",]
knitr::kable(selIM1)
```

##Reliability
```{r}
#Reliability
#fiab <- round(reliability(fit1),3)#First order only
#knitr::kable(fiab)

#fiab <- round(reliability(fit1),3)#First order only
#knitr::kable(fiab)
reliability(fit1)
```

# Analise de Invariancia
###Grupo = Amostra (0 - Teste, 1 - Validacao)
```{r new meas invariance, include=T, echo=F, eval=T}
## compare several invariance models

grupo ="Amostra"
test.seq <- c("thresholds","loadings","residuals")
meq.list <- list()
for (i in 0L:length(test.seq)) {
  if (i == 0L) {
    meq.label <- "configural"
    group.equal <- ""
  } else {
    meq.label <- test.seq[i]
    group.equal <- test.seq[1:i]
  }
  meq.list[[meq.label]] <- measEq.syntax(configural.model = modelo1,
                                         ridge=1e-05,
                                         data = bd,
                                         ordered = it3r,
                                         parameterization = "theta",
                                         #ID.fac = "ul",
                                         ID.cat = "Wu.Estabrook.2016",
                                         group = grupo,
                                         group.equal = group.equal,
                                         return.fit = TRUE)
}

compareFit(meq.list)
```


##Structural Model Analysis
```{r}
# Fit the Structural Model
#fit3 = sem(modelo3, data=bd,std.lv=T, ordered = it, estimator="WLSMV")
#fit4 = sem(modelo4, data=bd,std.lv=T, ordered = it, estimator="WLSMV")
fit10 = sem(modelo10, data=bd,std.lv=T, ordered = it10, estimator="WLSMV")
#Factor loadings
#summary(fit3, fit.measures=TRUE, standardized=TRUE)
#inspect(fit3, 'r2')
#summary(fit4, fit.measures=TRUE, standardized=TRUE)
#inspect(fit4, 'r2')
summary(fit10, fit.measures=TRUE, standardized=TRUE)
inspect(fit10, 'r2')
```