MS 2 - inhibitory control

1 Data processing

1.1 Packages

pacman::p_load(tidyverse, readxl, effsize, janitor, irr, dplyr, tidyr, ggplot2, openxlsx, RColorBrewer, car, semPlot,reshape2, mediation, afex, ez, lavaan, tidyLPA, compareGroups, car, emmeans, arsenal)

1.2 Data - ASD group

df_asd <- read_excel("C:/Users/louis/OneDrive/Área de Trabalho/Artigo 2 - Mestrado/df_asd.xlsx")

1.3 Check NA’s

variaveis <- c("go_no_go", "go_no_go_z", "stroop_a", "stroop_a_z", "stroop_a_errors", "stroop_a_errors_z", "stroop_b", "stroop_b_z",                "stroop_b_errors", "stroop_b_errors_z", "stroop_c", "stroop_c_z", "stroop_c_errors",
               "stroop_c_errors_z", "interference", "fdt_reading", "fdt_reading_z",
               "fdt_reading_errors", "fdt_reading_errors_z", "fdt_counting", "fdt_counting_z",
               "fdt_counting_errors", "fdt_counting_errors_z", "fdt_choosing", "fdt_choosing_z",
               "fdt_choosing_errors", "fdt_choosing_errors_z", "fdt_shifting", "fdt_shifting_z",
               "fdt_shifting_errors", "fdt_shifting_errors_z", "fdt_inhibition", "fdt_inhibition_z",
               "fdt_flexibility", "fdt_flexibility_z")


sapply(df_asd[variaveis], function(x) sum(is.na(x)))

             go_no_go            go_no_go_z              stroop_a            stroop_a_z       stroop_a_errors 
                    8                     8                     5                     8                     8 
    stroop_a_errors_z              stroop_b            stroop_b_z       stroop_b_errors     stroop_b_errors_z 
                    8                     8                    11                    11                    11 
             stroop_c            stroop_c_z       stroop_c_errors     stroop_c_errors_z          interference 
                    8                    11                    11                    11                     9 
          fdt_reading         fdt_reading_z    fdt_reading_errors  fdt_reading_errors_z          fdt_counting 
                   24                    24                    24                    24                    24 
       fdt_counting_z   fdt_counting_errors fdt_counting_errors_z          fdt_choosing        fdt_choosing_z 
                   24                    24                    24                    24                    24 
  fdt_choosing_errors fdt_choosing_errors_z          fdt_shifting        fdt_shifting_z   fdt_shifting_errors 
                   24                    24                    24                    24                    24 
fdt_shifting_errors_z        fdt_inhibition      fdt_inhibition_z       fdt_flexibility     fdt_flexibility_z 
                   24                    24                    24                    24                    24 

1.4 Imputation of missing values using the mean

df_asd[variaveis] <- lapply(df_asd[variaveis], function(x) {
  ifelse(is.na(x), mean(x, na.rm = TRUE), x)
})

1.5 Transform - Sociodemographic variables

df_asd <- df_asd %>%
  dplyr::mutate(
    sex = as.factor(sex),
    education_level = as.factor(education_level),
    comorbidities = as.factor(comorbidities),
    diagnosis = as.factor(diagnosis)
  )

1.6 Rename - Sex

df_asd <- df_asd %>%
  dplyr::mutate(
    sex = dplyr::recode_factor(sex, 
                               `1` = "Feminine", 
                               `2` = "Masculine")
  )

1.7 Transform in numeric

vars_to_convert <- c("go_no_go", "go_no_go_z", "stroop_a", "stroop_a_z", "stroop_a_errors", "stroop_a_errors_z",
                     "stroop_b", "stroop_b_z", "stroop_b_errors", "stroop_b_errors_z",
                     "stroop_c", "stroop_c_z", "stroop_c_errors", "stroop_c_errors_z",
                     "interference", "fdt_reading", "fdt_reading_z", "fdt_reading_errors",
                     "fdt_reading_errors_z", "fdt_counting", "fdt_counting_z", "fdt_counting_errors",
                     "fdt_counting_errors_z", "fdt_choosing", "fdt_choosing_z", "fdt_choosing_errors",
                     "fdt_choosing_errors_z", "fdt_shifting", "fdt_shifting_z", "fdt_shifting_errors",
                     "fdt_shifting_errors_z", "fdt_inhibition", "fdt_inhibition_z", "fdt_flexibility",
                     "fdt_flexibility_z")

df_asd <- df_asd %>%
  mutate(across(all_of(vars_to_convert), ~ round(as.numeric(.), 2)))


summary(df_asd[vars_to_convert])

    go_no_go      go_no_go_z         stroop_a        stroop_a_z      stroop_a_errors  stroop_a_errors_z    stroop_b    
 Min.   :0.67   Min.   :-10.330   Min.   :  8.00   Min.   :-9.2100   Min.   :0.0000   Min.   :-13.170   Min.   :14.00  
 1st Qu.:0.94   1st Qu.: -1.340   1st Qu.: 16.25   1st Qu.:-1.6275   1st Qu.:0.0000   1st Qu.: -0.910   1st Qu.:23.25  
 Median :0.95   Median : -1.000   Median : 27.50   Median :-1.2900   Median :0.0000   Median :  0.170   Median :32.48  
 Mean   :0.94   Mean   : -1.336   Mean   : 29.49   Mean   :-1.3212   Mean   :0.6904   Mean   : -0.913   Mean   :32.48  
 3rd Qu.:1.00   3rd Qu.:  0.670   3rd Qu.: 33.50   3rd Qu.: 0.0375   3rd Qu.:0.6900   3rd Qu.:  0.250   3rd Qu.:35.00  
 Max.   :1.00   Max.   :  0.670   Max.   :115.00   Max.   : 2.3000   Max.   :8.0000   Max.   :  0.330   Max.   :77.00  
   stroop_b_z      stroop_b_errors stroop_b_errors_z    stroop_c        stroop_c_z     stroop_c_errors  stroop_c_errors_z
 Min.   :-4.7900   Min.   : 0.00   Min.   :-19.600   Min.   : 21.00   Min.   :-6.490   Min.   : 0.000   Min.   :-7.1400  
 1st Qu.:-1.3900   1st Qu.: 0.00   1st Qu.: -1.870   1st Qu.: 30.25   1st Qu.:-1.718   1st Qu.: 0.000   1st Qu.:-0.6700  
 Median :-0.8500   Median : 0.00   Median :  0.000   Median : 47.50   Median :-1.020   Median : 0.500   Median :-0.5150  
 Mean   :-0.8542   Mean   : 1.41   Mean   : -1.873   Mean   : 48.64   Mean   :-1.020   Mean   : 2.153   Mean   :-0.6722  
 3rd Qu.: 0.1775   3rd Qu.: 1.41   3rd Qu.:  0.250   3rd Qu.: 55.00   3rd Qu.: 0.300   3rd Qu.: 2.150   3rd Qu.: 0.4400  
 Max.   : 1.6300   Max.   :10.00   Max.   :  0.500   Max.   :140.00   Max.   : 1.270   Max.   :17.000   Max.   : 0.8100  
  interference    fdt_reading    fdt_reading_z     fdt_reading_errors fdt_reading_errors_z  fdt_counting    fdt_counting_z  
 Min.   :1.000   Min.   :17.00   Min.   :-7.4200   Min.   :0.000      Min.   :-6.330       Min.   : 19.00   Min.   :-5.580  
 1st Qu.:1.570   1st Qu.:34.00   1st Qu.:-1.0700   1st Qu.:0.000      1st Qu.:-0.550       1st Qu.: 48.75   1st Qu.:-1.130  
 Median :2.025   Median :38.31   Median :-1.0700   Median :0.150      Median :-0.550       Median : 55.96   Median :-1.130  
 Mean   :2.041   Mean   :38.31   Mean   :-1.0696   Mean   :0.152      Mean   :-0.551       Mean   : 55.96   Mean   :-1.128  
 3rd Qu.:2.192   3rd Qu.:38.31   3rd Qu.:-0.3575   3rd Qu.:0.150      3rd Qu.: 0.000       3rd Qu.: 55.96   3rd Qu.:-0.885  
 Max.   :6.000   Max.   :68.00   Max.   : 1.6600   Max.   :2.000      Max.   : 0.330       Max.   :120.00   Max.   : 1.710  
 fdt_counting_errors fdt_counting_errors_z  fdt_choosing    fdt_choosing_z   fdt_choosing_errors fdt_choosing_errors_z
 Min.   :0.0000      Min.   :-3.6000       Min.   : 15.00   Min.   :-4.010   Min.   :0.000       Min.   :-2.700       
 1st Qu.:0.0000      1st Qu.:-0.0900       1st Qu.: 77.50   1st Qu.:-1.280   1st Qu.:1.000       1st Qu.: 0.090       
 Median :0.4200      Median :-0.0900       Median : 88.92   Median :-1.280   Median :1.880       Median : 0.090       
 Mean   :0.4216      Mean   :-0.0874       Mean   : 88.92   Mean   :-1.282   Mean   :1.882       Mean   : 0.093       
 3rd Qu.:0.4200      3rd Qu.: 0.3900       3rd Qu.: 88.92   3rd Qu.:-1.055   3rd Qu.:1.880       3rd Qu.: 0.620       
 Max.   :3.0000      Max.   : 0.5000       Max.   :176.00   Max.   : 2.720   Max.   :9.000       Max.   : 1.000       
  fdt_shifting    fdt_shifting_z   fdt_shifting_errors fdt_shifting_errors_z fdt_inhibition  fdt_inhibition_z
 Min.   : 45.00   Min.   :-3.300   Min.   : 0.000      Min.   :-2.6800       Min.   : -2.0   Min.   :-3.790  
 1st Qu.: 99.75   1st Qu.:-1.155   1st Qu.: 2.000      1st Qu.:-0.2000       1st Qu.: 48.5   1st Qu.:-1.130  
 Median :105.92   Median :-1.050   Median : 3.380      Median :-0.2000       Median : 52.5   Median :-1.130  
 Mean   :105.92   Mean   :-1.051   Mean   : 3.382      Mean   :-0.2020       Mean   : 52.5   Mean   :-1.131  
 3rd Qu.:105.92   3rd Qu.:-1.035   3rd Qu.: 3.380      3rd Qu.: 0.1375       3rd Qu.: 52.5   3rd Qu.:-0.910  
 Max.   :177.00   Max.   : 1.280   Max.   :13.000      Max.   : 0.9000       Max.   :118.0   Max.   : 2.200  
 fdt_flexibility  fdt_flexibility_z
 Min.   : 23.00   Min.   :-3.450   
 1st Qu.: 64.17   1st Qu.:-1.653   
 Median : 67.69   Median :-1.070   
 Mean   : 67.69   Mean   :-1.068   
 3rd Qu.: 67.92   3rd Qu.:-1.070   
 Max.   :130.00   Max.   : 0.920   

2 Data Analysis

2.1 Check - participants characteristics

df_asd %>%
  dplyr::select(age, sex,education_level, comorbidities, diagnosis, iq, verbal_index, performance_index, go_no_go_z,stroop_a_z, stroop_b_z, stroop_c_z,fdt_inhibition_z, fdt_flexibility_z) %>%
  arsenal::tableby(~ ., digits = 2, data = .) %>% summary(text= T) %>% 
  as.data.frame()

2.2 How many participants has a IQ <70?

 df_asd %>%
  dplyr::filter(iq < 70) %>%
  dplyr::summarise(count = n()) # São 8, sendo que 8 tiveram o diagnóstico de TEA e DI

2.3 Select data

df_selected <- df_asd %>%
  dplyr::select("go_no_go_z", "stroop_a_z", "stroop_b_z", "stroop_c_z", "fdt_inhibition_z")

2.4 Clustering (traditional)

2.4.1 How many classes?

library(factoextra)
fviz_nbclust(df_selected, kmeans, method = "wss")+
  geom_vline(xintercept = 4, linetype = 2) # incialmente sugerem quatro 

set.seed(123)
km.res <- kmeans(df_selected, 3, nstart = 25)

palette <- c("#255293", "#db0a16", "#f8c72d")

# Visualize
library(factoextra)

fviz_cluster(km.res, data = df_selected,
             ellipse.type = "convex",
             #palette = palette,
             ggtheme = theme_minimal())

fviz_cluster(km.res, data=df_selected,
             palette = c("#2E9FDF", "#00AFBB", "#E7B800", "#FC4E07"),
             ellipse.type="euclid",
             star.plot=TRUE,
             repel=TRUE,
             ggtheme=theme_minimal()
             )

2.4.2 Final model

#set.seed(123)
#lpa_exploratory = estimate_profiles(df_selected,1:5, package = "MplusAutomation",
                     # variances = c("equal", "varying"),
                     # covariances = c("zero", "varying"))
#compare_solutions(lpa_exploratory,statistics = c("LogLik","AWE", "AIC", "BIC", "SABIC", "BLRT_val", "BLRT_p","Entropy"))

Now I’ll run this solution

mod_lap = estimate_profiles(df_selected, n_profiles = 3)

get_fit(mod_lap)

get_estimates(mod_lap)

df_cluster = get_data(mod_lap)

tabyl(df_cluster$Class)

 df_cluster$Class  n percent
                1 29    0.58
                2  7    0.14
                3 14    0.28

2.4.3 Plot fit

prof_estimates <- get_estimates(mod_lap)

prof_estimates %>%
  filter(Category == "Means") %>%
  mutate(Class = str_replace_all(Class, "^", "Profile "),
         Parameter = str_remove_all(Parameter, "_centered")) %>%
  mutate(Parameter = str_remove_all(Parameter, "_zscore"))%>%
  mutate(Parameter = str_remove_all(Parameter, "(total)"))%>%
  mutate(Parameter = str_remove_all(Parameter,   "\\(.*")) %>%
  ggplot(aes(x = Class , y = Estimate , group=Parameter, fill = Parameter)) +
  geom_bar(stat = "summary", position = "dodge") +
  theme_light() + theme(legend.position = "bottom") +  scale_fill_grey(start = 0, end = .9)

bind_cols(
  df_asd %>% dplyr::select("diagnosis"),
  df_cluster
)

 bind_cols(
  df_asd %>% dplyr::select(diagnosis),
  df_cluster
) %>% 
  tabyl(diagnosis, Class, show_na = FALSE) %>% 
  adorn_percentages("row") %>% 
  adorn_pct_formatting(digits = 1) %>% 
  adorn_ns()

 diagnosis          1         2          3
       ASD 62.2% (23)  8.1% (3) 29.7% (11)
  ASD+ADHD 40.0%  (2) 20.0% (1) 40.0%  (2)
   ASD+IDD 50.0%  (4) 37.5% (3) 12.5%  (1)

3 Predictors

df_asd <- bind_cols(
  df_asd,
  df_cluster %>% dplyr::select("Class")
)

rm(res)
res<-compareGroups(Class ~ age + sex + iq + verbal_index + performance_index + fdt_flexibility_z + working_memory_index + processing_speed_index, data=df_asd)
createTable(res)

--------Summary descriptives table by 'Class'---------

_______________________________________________________________________ 
                            1            2            3       p.overall 
                           N=29         N=7          N=14               
¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯ 
age                    10.1 (2.71)  8.86 (2.04)  10.6 (2.28)    0.346   
sex:                                                            0.132   
    Feminine            6 (20.7%)    3 (42.9%)    1 (7.14%)             
    Masculine           23 (79.3%)   4 (57.1%)    13 (92.9%)            
iq                     90.0 (20.1)  79.6 (12.6)  83.3 (10.7)    0.249   
verbal_index           86.3 (21.9)  78.6 (12.3)  85.3 (12.1)    0.615   
performance_index      95.4 (17.6)  83.3 (15.7)  87.2 (12.2)    0.119   
fdt_flexibility_z      -1.09 (0.85) -1.52 (0.88) -0.79 (0.92)   0.201   
working_memory_index   85.6 (16.9)  76.5 (10.3)  86.3 (20.0)    0.467   
processing_speed_index 92.6 (15.6)  76.8 (11.0)  90.5 (12.2)    0.065   
¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯ 

3.1 Getting to Know the People in the Clusters

df_asd %>% 
  dplyr::select(age, sex,education_level, comorbidities, diagnosis, iq, verbal_index, performance_index, go_no_go_z,stroop_a_z, stroop_b_z, stroop_c_z,fdt_inhibition_z, fdt_flexibility_z, working_memory_index, processing_speed_index) %>%
  arsenal::tableby(df_asd $Class ~ .,  digits=2, data = .) %>% summary(text= T)%>% as.data.frame()

3.2 ANOVA - Stroop

# ANOVA para cada variável
anova_stroop_a <- aov(stroop_a_z ~ Class, data = df_asd)
anova_stroop_b <- aov(stroop_b_z ~ Class, data = df_asd)
anova_stroop_c <- aov(stroop_c_z ~ Class, data = df_asd)

# Resumo dos resultados
summary(anova_stroop_a)

            Df Sum Sq Mean Sq F value Pr(>F)  
Class        1  17.36  17.359   3.872 0.0549 .
Residuals   48 215.19   4.483                 
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

summary(anova_stroop_b) # diferença significativa

            Df Sum Sq Mean Sq F value  Pr(>F)   
Class        1  16.58  16.582   9.865 0.00288 **
Residuals   48  80.68   1.681                   
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

summary(anova_stroop_c) # diferença significativa 

            Df Sum Sq Mean Sq F value Pr(>F)  
Class        1  15.53   15.53   5.586 0.0222 *
Residuals   48 133.45    2.78                 
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

# Post-hoc usando Tukey HSD para uma das variáveis
posthoc_a <- emmeans(anova_stroop_a, pairwise ~ Class)
posthoc_a

$emmeans
 Class emmean    SE df lower.CL upper.CL
   1.7  -1.32 0.299 48    -1.92   -0.719

Confidence level used: 0.95 

$contrasts
 contrast  estimate SE df z.ratio p.value
 (nothing)   nonEst NA NA      NA      NA

# ANOVA para cada variável
anova_go_no_go <- aov(go_no_go_z ~ Class, data = df_asd)

# Resumo dos resultados
summary(anova_go_no_go)

            Df Sum Sq Mean Sq F value Pr(>F)
Class        1  10.02  10.016     1.6  0.212
Residuals   48 300.41   6.258               

# Post-hoc usando Tukey HSD para uma das variáveis
posthoc_go <- emmeans(anova_go_no_go, pairwise ~ Class)
posthoc_go

$emmeans
 Class emmean    SE df lower.CL upper.CL
   1.7  -1.34 0.354 48    -2.05   -0.625

Confidence level used: 0.95 

$contrasts
 contrast  estimate SE df z.ratio p.value
 (nothing)   nonEst NA NA      NA      NA