MS 2 - inhibitory control

1 ASD Group

2 Data processing

2.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, Hmisc, factoextra, boot, lsr)

2.2 Data - ASD group

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

2.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 

2.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)
})

2.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)
  )

2.6 Rename - Sex

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

2.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)))

2.8 Creating medication subgroups

df_asd <- df_asd %>%
  dplyr::mutate(
    medication_class = case_when(
      !is.na(initial_medication) & initial_medication %in% c("Ritalina") ~ "Psychostimulants",
      !is.na(initial_medication) & initial_medication %in% c("Atentah") ~ "Non-stimulants (ADHD)",
      !is.na(initial_medication) & initial_medication %in% c("Risperdal", "Risperidona", "Risperidona e Escitalopram", 
                                                             "Risperidona e Quetiapina", "Neuleptil", "neuleptil", "Neozine") ~ "Antipsychotics",
      !is.na(initial_medication) & initial_medication %in% c("Fluoxetina", "Sertralina", "Daforin", "Escitalopram") ~ "Antidepressants",
      !is.na(initial_medication) & initial_medication %in% c("Topiramato") ~ "Mood stabilizers",
      !is.na(initial_medication) & initial_medication %in% c("without medication", "N-Miss") ~ "No medication",
      FALSE ~ "Others" 
    )
  )

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

df_asd <- df_asd %>%
  dplyr::mutate(
    psychomotor_agitation = as.factor(psychomotor_agitation),
    restlessness = as.factor(restlessness),
    impulsiveness = as.factor(impulsiveness),
    inattention = as.factor(inattention),
    verbal_aggressiveness = as.factor(verbal_aggressiveness),
    physical_aggressiveness = as.factor(physical_aggressiveness),
    difficulty_with_frustration = as.factor(difficulty_with_frustration),
    restricted_patterns = as.factor(restricted_patterns),
    behavioral_inflexibility = as.factor(behavioral_inflexibility),
    repetitive_patterns = as.factor(repetitive_patterns),
    social_difficulty = as.factor(social_difficulty),
    speech_delay = as.factor(speech_delay),
    depression_symptoms = as.factor(depression_symptoms),
    anxiety_symptoms = as.factor(anxiety_symptoms)
  )

3 Data Analysis

3.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()

3.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

3.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")

3.4 Clustering (traditional)

3.4.1 How many classes?

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()
             )

3.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 28 0.5283019
                2 11 0.2075472
                3 14 0.2641509

3.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) +
  scale_y_continuous(limits = c(NA, 2), breaks = scales::pretty_breaks(n = 5))

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 60.5% (23) 13.2% (5) 26.3% (10)
  ASD+ADHD 40.0%  (2) 20.0% (1) 40.0%  (2)
   ASD+IDD 30.0%  (3) 50.0% (5) 20.0%  (2)

4 Predictors

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

rm(res)

Warning: objeto 'res' não encontrado

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=28         N=11         N=14               
¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯ 
age                    10.1 (2.73)  8.91 (1.87)  10.6 (2.28)    0.233   
sex:                                                            0.360   
    Feminine            7 (25.0%)    3 (27.3%)    1 (7.14%)             
    Masculine           21 (75.0%)   8 (72.7%)    13 (92.9%)            
iq                     89.1 (19.5)  83.0 (18.7)  83.3 (10.7)    0.479   
verbal_index           85.9 (21.0)  79.7 (19.2)  85.3 (12.1)    0.644   
performance_index      94.1 (17.4)  88.7 (18.6)  87.2 (12.2)    0.396   
fdt_flexibility_z      -1.02 (0.76) -1.75 (0.81) -0.81 (0.93)   0.016   
working_memory_index   85.6 (17.2)  79.9 (11.3)  86.3 (20.0)    0.640   
processing_speed_index 93.0 (14.8)  80.8 (15.1)  90.5 (12.2)    0.111   
¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯ 

4.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()

4.2 Medication use

df_asd %>% 
  dplyr::select(medication_class) %>%
  arsenal::tableby(df_asd $Class ~ .,  digits=2, data = .) %>% summary(text= T)%>% as.data.frame()

4.3 Plot - Medication use

df_pizza <- df_asd %>%
  filter(!is.na(medication_class) & medication_class != "Others") %>%  # Remove NAs e "Others"
  group_by(Class, medication_class) %>%
  summarise(count = n()) %>%
  mutate(percentage = count / sum(count) * 100)

`summarise()` has grouped output by 'Class'. You can override using the `.groups` argument.

ggplot(df_pizza, aes(x = "", y = percentage, fill = medication_class)) +
  geom_bar(stat = "identity", width = 1, color = "white") +
  coord_polar(theta = "y") +
  facet_wrap(~ Class, ncol = 3) +  # Coloca os gráficos lado a lado (3 colunas)
  labs(
    title = "Distribution of Medication Classes by Class",
    fill = "Medication Class",
    y = "",
    x = ""
  ) +
  geom_text(aes(label = paste0(round(percentage, 1), "%")), 
            position = position_stack(vjust = 0.5), 
            size = 4, color = "black") + # Rótulos de porcentagem no centro das fatias
  theme_minimal() +
  theme(
    plot.title = element_text(size = 16, face = "bold", hjust = 0.5), # Centraliza e aumenta o título
    strip.text = element_text(size = 14, face = "bold"), # Aumenta o texto das facetas
    axis.text = element_blank(),
    axis.ticks = element_blank(),
    panel.grid = element_blank(),
    legend.title = element_text(size = 12, face = "bold"), # Ajusta o tamanho do título da legenda
    legend.text = element_text(size = 10),
    legend.position = "bottom"  # Move a legenda para baixo
  ) +
  guides(fill = guide_legend(title.position = "top", title.hjust = 0.5)) # Centraliza o título da legenda

NA
NA

4.4 Behavioral information

df_asd %>%
  dplyr::select(
    psychomotor_agitation, restlessness, impulsiveness, inattention, verbal_aggressiveness, 
    physical_aggressiveness, difficulty_with_frustration, restricted_patterns, behavioral_inflexibility, 
    repetitive_patterns, social_difficulty, speech_delay, depression_symptoms, anxiety_symptoms
  ) %>%
  arsenal::tableby(df_asd$Class ~ ., digits = 2, data = .) %>%
  summary(text = TRUE) %>%
  as.data.frame()

df_pizza_2 <- df_asd %>%
  dplyr::select(
    psychomotor_agitation, restlessness, impulsiveness, inattention, verbal_aggressiveness, 
    physical_aggressiveness, difficulty_with_frustration, restricted_patterns, behavioral_inflexibility, 
    repetitive_patterns, social_difficulty, speech_delay, depression_symptoms, anxiety_symptoms, Class
  ) %>%
  gather(key = "behavior", value = "score", -Class) %>%
  mutate(score = as.numeric(score)) %>%  # Garantir que 'score' seja numérico
  group_by(Class, behavior) %>%
  summarise(total_score = sum(score, na.rm = TRUE)) %>%
  group_by(Class) %>%
  mutate(percentage = total_score / sum(total_score) * 100)

`summarise()` has grouped output by 'Class'. You can override using the `.groups` argument.

custom_colors <- c(
  "#66c2a5", "#fc8d62", "#8da0cb", "#e78ac3", "#a6d854",
  "#ffd92f", "#e5c494", "#b3b3b3", "#a6cee3", "#1f78b4",
  "#b2df8a", "#33a02c", "#fb9a99", "#e31a1c", "#fdbf6f"
)

options(repr.plot.width = 12, repr.plot.height = 8)

ggplot(df_pizza_2, aes(x = Class, y = percentage, fill = behavior)) +
  geom_bar(stat = "identity", position = "dodge", color = "gray80") +
  labs(
    title = "Distribution of Symptoms by Groups",
    fill = "Symptom",
    y = "Percentage",
    x = "Class"
  ) +
  geom_text(aes(label = paste0(round(percentage, 1), "%")), 
            position = position_dodge(width = 0.8), 
            size = 3.5, color = "black") +
  theme_minimal(base_size = 14) +
  theme(
    plot.title = element_text(size = 18, face = "bold", hjust = 0.5),
    axis.text.x = element_text(size = 12, face = "bold"),
    axis.text.y = element_text(size = 12),
    legend.title = element_text(size = 14, face = "bold"),
    legend.text = element_text(size = 12),
    legend.position = "bottom",
    legend.key.size = unit(0.8, "cm")
  ) +
  scale_fill_manual(values = custom_colors) + # Aplicando cores personalizadas
  guides(fill = guide_legend(nrow = 2, byrow = TRUE))

4.5 ANOVA - Stroop

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  21.11  21.114   4.987   0.03 *
Residuals   51 215.95   4.234                 
---
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   8.54   8.536   4.452 0.0398 *
Residuals   51  97.79   1.918                 
---
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  11.17  11.169   4.124 0.0475 *
Residuals   51 138.11   2.708                 
---
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.74  -1.37 0.283 51    -1.94   -0.802

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   19.1  19.091   2.754  0.103
Residuals   51  353.5   6.932               

# 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.74  -1.45 0.362 51    -2.18   -0.729

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_fdt  <- aov(fdt_inhibition_z ~ Class, data = df_asd)

# Resumo dos resultados
summary(anova_fdt)

            Df Sum Sq Mean Sq F value Pr(>F)  
Class        1   4.83   4.835   3.729  0.059 .
Residuals   51  66.12   1.297                 
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

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

$emmeans
 Class emmean    SE df lower.CL upper.CL
  1.74  -1.26 0.156 51    -1.58   -0.948

Confidence level used: 0.95 

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

4.6 Correlation

# Selecionar variáveis
correlation_data <- df_asd %>%
  dplyr::select(go_no_go_z, stroop_a_z, stroop_b_z, stroop_c_z, fdt_inhibition_z, 
                fdt_flexibility_z, performance_index, working_memory_index)

# Matriz de correlação
correlation_matrix <- cor(correlation_data, use = "pairwise.complete.obs", method = "pearson")

# Visualização da matriz
library(corrplot)
corrplot::corrplot(correlation_matrix, method = "color", type = "upper", addCoef.col = "black", tl.col = "black")

# Obter matriz de correlações e valores-p
cor_results <- Hmisc::rcorr(as.matrix(correlation_data), type = "pearson")

# Matriz de correlações
cor_matrix <- cor_results$r

# Matriz de valores-p
p_matrix <- cor_results$P

# Criar matriz combinada com correlação e valores-p
combined_matrix <- paste0(
  format(round(cor_matrix, 2), nsmall = 2),  # Coeficiente de correlação
  " (p=", format(round(p_matrix, 3), nsmall = 3), ")"  # Valor-p
)

# Reformatar em formato de matriz
dim(combined_matrix) <- dim(cor_matrix)
rownames(combined_matrix) <- rownames(cor_matrix)
colnames(combined_matrix) <- colnames(cor_matrix)

# Visualizar matriz combinada
print(combined_matrix, quote = FALSE)

                     go_no_go_z      stroop_a_z      stroop_b_z      stroop_c_z      fdt_inhibition_z fdt_flexibility_z
go_no_go_z            1.00 (p=   NA)  0.09 (p=0.501)  0.11 (p=0.449)  0.18 (p=0.187)  0.29 (p=0.036)   0.21 (p=0.123)  
stroop_a_z            0.09 (p=0.501)  1.00 (p=   NA)  0.49 (p=0.000)  0.46 (p=0.001)  0.22 (p=0.121)   0.09 (p=0.522)  
stroop_b_z            0.11 (p=0.449)  0.49 (p=0.000)  1.00 (p=   NA)  0.72 (p=0.000)  0.14 (p=0.302)   0.31 (p=0.026)  
stroop_c_z            0.18 (p=0.187)  0.46 (p=0.001)  0.72 (p=0.000)  1.00 (p=   NA)  0.14 (p=0.301)   0.29 (p=0.034)  
fdt_inhibition_z      0.29 (p=0.036)  0.22 (p=0.121)  0.14 (p=0.302)  0.14 (p=0.301)  1.00 (p=   NA)   0.40 (p=0.003)  
fdt_flexibility_z     0.21 (p=0.123)  0.09 (p=0.522)  0.31 (p=0.026)  0.29 (p=0.034)  0.40 (p=0.003)   1.00 (p=   NA)  
performance_index    -0.11 (p=0.430)  0.10 (p=0.483)  0.16 (p=0.260) -0.12 (p=0.420) -0.04 (p=0.800)  -0.14 (p=0.325)  
working_memory_index  0.21 (p=0.151)  0.09 (p=0.551)  0.24 (p=0.099)  0.18 (p=0.217) -0.03 (p=0.813)   0.02 (p=0.901)  
                     performance_index working_memory_index
go_no_go_z           -0.11 (p=0.430)    0.21 (p=0.151)     
stroop_a_z            0.10 (p=0.483)    0.09 (p=0.551)     
stroop_b_z            0.16 (p=0.260)    0.24 (p=0.099)     
stroop_c_z           -0.12 (p=0.420)    0.18 (p=0.217)     
fdt_inhibition_z     -0.04 (p=0.800)   -0.03 (p=0.813)     
fdt_flexibility_z    -0.14 (p=0.325)    0.02 (p=0.901)     
performance_index     1.00 (p=   NA)    0.51 (p=0.000)     
working_memory_index  0.51 (p=0.000)    1.00 (p=   NA)     

5 Non-Clinical Group

5.1 Data

df_non_clinical <- read_excel("C:/Users/louis/Downloads/Controle - escolinha (Louise).xlsx")

5.2 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_non_clinical[variaveis], function(x) sum(is.na(x)))

             go_no_go            go_no_go_z              stroop_a            stroop_a_z       stroop_a_errors 
                    2                     3                     1                     1                     1 
    stroop_a_errors_z              stroop_b            stroop_b_z       stroop_b_errors     stroop_b_errors_z 
                    1                     1                     1                     1                     1 
             stroop_c            stroop_c_z       stroop_c_errors     stroop_c_errors_z          interference 
                    1                     1                     1                     1                     1 
          fdt_reading         fdt_reading_z    fdt_reading_errors  fdt_reading_errors_z          fdt_counting 
                    0                     0                     0                     0                     0 
       fdt_counting_z   fdt_counting_errors fdt_counting_errors_z          fdt_choosing        fdt_choosing_z 
                    0                     0                     0                     0                     0 
  fdt_choosing_errors fdt_choosing_errors_z          fdt_shifting        fdt_shifting_z   fdt_shifting_errors 
                    0                     0                     1                     1                     1 
fdt_shifting_errors_z        fdt_inhibition      fdt_inhibition_z       fdt_flexibility     fdt_flexibility_z 
                    1                     0                     0                     1                     1 

5.3 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_non_clinical <- df_non_clinical %>%
  mutate(across(all_of(vars_to_convert), ~ round(as.numeric(.), 2)))

5.4 Imputation of missing values using the mean

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

5.5 Rename - Sex

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

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

df_non_clinical <- df_non_clinical %>%
  mutate_all(~ replace(., . == 999, NA))

5.6 Check - participants characteristics

df_non_clinical %>%
  dplyr::select(age, sex,education_level, 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()

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

df_selected_2 <- df_selected_2 %>%
  mutate_all(~ ifelse(is.na(.), mean(., na.rm = TRUE), .))

5.7 Clustering (traditional)

5.7.1 How many classes?

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

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

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


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

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

mod_lap_2 = estimate_profiles(df_selected_2, n_profiles = 2)

get_fit(mod_lap)

get_estimates(mod_lap_2)

df_cluster_2 = get_data(mod_lap_2)

tabyl(df_cluster_2$Class)

 df_cluster_2$Class  n   percent
                  1 34 0.7391304
                  2 12 0.2608696

5.7.2 Plot fit

prof_estimates <- get_estimates(mod_lap_2)

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) +
  scale_y_continuous(limits = c(NA, 2), breaks = scales::pretty_breaks(n = 5))

df_non_clinical<- bind_cols(
  df_non_clinical,
  df_cluster_2 %>% 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_non_clinical)

Warning: Aproximação do qui-quadrado pode estar incorreta

createTable(res)

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

__________________________________________________________ 
                            1            2       p.overall 
                           N=34         N=12               
¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯ 
age                    8.09 (1.31)  9.08 (1.31)    0.036   
sex:                                               0.513   
    Feminine            13 (38.2%)   6 (50.0%)             
    Masculine           21 (61.8%)   6 (50.0%)             
iq                     87.9 (13.6)  89.0 (16.7)    0.833   
verbal_index           82.8 (16.1)  80.6 (21.7)    0.796   
performance_index      93.2 (13.4)  99.8 (16.8)    0.333   
fdt_flexibility_z      -1.26 (1.75) -0.82 (1.06)   0.321   
working_memory_index   91.5 (13.6)  88.2 (10.2)    0.604   
processing_speed_index  100 (13.3)   101 (8.81)    0.954   
¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯ 

5.8 Getting to Know the People in the Clusters

df_non_clinical %>% 
  dplyr::select(age, sex,education_level, 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, cdi_classification,z_inattention, z_hyperactivity_impulsivity, z_general_anxiety) %>%
  arsenal::tableby(df_non_clinical $Class ~ .,  digits=2, data = .) %>% summary(text= T)%>% as.data.frame()

6 Comparing Groups (ASD vs Non Clinical Profiles)

6.1 Fixing the database

6.1.1 Create Group variable

df_non_clinical <- df_non_clinical %>%
  mutate(group = "non clinical")

df_asd <- df_asd %>%
  mutate(group = "asd")

df_combined <- bind_rows(df_non_clinical, df_asd)

df_combined %>% 
  dplyr::select(age, sex,education_level, 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, stroop_a_errors_z, stroop_b_errors_z, stroop_c_errors_z) %>%
  arsenal::tableby(df_combined $group ~ .,  digits=2, data = .) %>% summary(text= T)%>% as.data.frame()

6.2 ANCOVA

ancova_go_nc <- aov(go_no_go_z ~ group + iq + age, data = df_combined)
ancova_stroop_a_nc <- aov(stroop_a_z ~ group + iq + age + education_level, data = df_combined)
ancova_stroop_b_nc <- aov(stroop_b_z ~ group + iq + age + education_level, data = df_combined)
ancova_stroop_c_nc <- aov(stroop_c_z ~ group + iq + age + education_level, data = df_combined)
ancova_fdt_nc <- aov(fdt_inhibition_z ~ group + iq + age + education_level, data = df_combined)

# Resumo dos resultados
summary(ancova_go_nc)

            Df Sum Sq Mean Sq F value Pr(>F)   
group        1    5.2    5.20   1.230 0.2704   
iq           1    0.7    0.75   0.177 0.6750   
age          1   33.3   33.33   7.886 0.0061 **
Residuals   91  384.6    4.23                  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
4 observations deleted due to missingness

summary(ancova_stroop_a_nc)

                Df Sum Sq Mean Sq F value Pr(>F)
group            1   2.52   2.522   0.735  0.394
iq               1   1.39   1.391   0.405  0.526
age              1   5.49   5.488   1.599  0.210
education_level 10  30.20   3.020   0.880  0.556
Residuals       80 274.65   3.433               
5 observations deleted due to missingness

summary(ancova_stroop_b_nc)

                Df Sum Sq Mean Sq F value Pr(>F)
group            1   0.53  0.5342   0.266  0.608
iq               1   0.45  0.4541   0.226  0.636
age              1   0.06  0.0581   0.029  0.865
education_level 10  20.82  2.0815   1.036  0.422
Residuals       80 160.75  2.0094               
5 observations deleted due to missingness

summary(ancova_stroop_c_nc)

                Df Sum Sq Mean Sq F value Pr(>F)
group            1   0.98   0.979   0.389  0.535
iq               1   0.12   0.116   0.046  0.830
age              1   4.96   4.960   1.970  0.164
education_level 10  23.35   2.335   0.927  0.513
Residuals       80 201.46   2.518               
5 observations deleted due to missingness

summary(ancova_fdt_nc)

                Df Sum Sq Mean Sq F value Pr(>F)  
group            1   0.02   0.021   0.007 0.9337  
iq               1   5.40   5.398   1.779 0.1860  
age              1   9.91   9.908   3.265 0.0745 .
education_level 10  28.17   2.817   0.928 0.5118  
Residuals       81 245.78   3.034                 
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
4 observations deleted due to missingness

6.3 ETA

eta_squared_go_nc <- etaSquared(ancova_go_nc)
eta_squared_stroop_a_nc <- etaSquared(ancova_stroop_a_nc)
eta_squared_stroop_b_nc <- etaSquared(ancova_stroop_b_nc)
eta_squared_stroop_c_nc <- etaSquared(ancova_stroop_c_nc)
eta_squared_fdt_nc <- etaSquared(ancova_fdt_nc)

# Exibir os resultados
eta_squared_go_nc

           eta.sq eta.sq.part
group 0.043201204  0.04544874
iq    0.009285187  0.01012968
age   0.078628401  0.07974684

eta_squared_stroop_a_nc

                      eta.sq  eta.sq.part
group           0.0110564752 0.0124926747
iq              0.0063720345 0.0072380512
age             0.0001255695 0.0001436548
education_level 0.0961047500 0.0990683862

eta_squared_stroop_b_nc

                      eta.sq  eta.sq.part
group           6.121838e-03 0.0069063647
iq              7.568578e-04 0.0008590502
age             8.782064e-05 0.0000997541
education_level 1.139865e-01 0.1146434307

eta_squared_stroop_c_nc

                      eta.sq  eta.sq.part
group           1.909123e-02 0.0214098140
iq              9.447366e-05 0.0001082535
age             9.810698e-03 0.0111178928
education_level 1.011577e-01 0.1038824582

eta_squared_fdt_nc

                      eta.sq  eta.sq.part
group           1.254345e-04 0.0001476127
iq              3.837007e-02 0.0432095896
age             9.344139e-05 0.0001099670
education_level 9.738935e-02 0.1028378823

df_long <- df_combined %>%
  pivot_longer(cols = c(go_no_go_z, stroop_a_z, stroop_b_z, stroop_c_z, fdt_inhibition_z),
               names_to = "variable",
               values_to = "value")

6.4 Plot - ASD x non clinical

# Reordene os níveis do fator para garantir que "asd" fique na frente
df_long$group <- factor(df_long$group, levels = c("non clinical", "asd"))

# Gráfico de densidade com marcadores em -1 e +1
ggplot(df_long, aes(x = value, fill = group)) +
  geom_density(alpha = 0.6, position = "identity") + # Curvas sobrepostas
  # Linhas verticais nos pontos -1 e +1
  geom_vline(xintercept = c(-1, 1), linetype = "dashed", color = "black", size = 0.8, alpha = 0.7, show.legend = FALSE) + # Linhas em -1 e +1
  facet_wrap(~ variable, scales = "free") +          # Uma curva para cada variável
  labs(title = "Distribution of Variables by Group",
       subtitle = "Score Z",
       x = "Standard-Distribuition",
       y = "Standard-Distribuition",
       fill = "Groups") +
  theme_minimal(base_size = 14) +
  theme(
    plot.title = element_text(hjust = 0.5, face = "bold"),
    plot.subtitle = element_text(hjust = 0.5, size = 12),
    legend.position = "top"
  ) +
  scale_fill_manual(values = c("#1f78b4", "#e31a1c")) + # Ajuste das cores dos grupos
  scale_x_continuous(breaks = seq(-5, 5, by = 1), limits = c(-5, 5)) # Ajuste dos marcadores do eixo X

NA
NA

7 The end