Cognitive Strategies In Older Autistic Adults

Code explanation

This markdown shows all output for the project on Cognitive Strategies in Older/Autistic Adults. Code is written by CT. Authors are CT, APG, EvdB, RC, CH, HMG.

Set-up

Required packages are loaded, workspace is defined, data is loaded and transformed into a final format. For more information on the loaded data check separate Qual/VPT/JND Rcodes.

### some basic settings
# Check packages installed:
required_packages <- c("rstudioapi","plyr", "dplyr", "stringr", 
                       "readr", "gridExtra", "moments", "readxl",
                       "tidyr", "ggplot2","ggpubr","extrafont",
                       "kableExtra", "broom", "knitr","broom.mixed",
                       "haven", "lubridate", "magrittr", "rstatix",
                       "ggbeeswarm", "ggstatsplot", "stats", "effsize",
                       "lmboot", "BayesFactor")
need_install <- required_packages[!(required_packages) %in% installed.packages()]

# Install packages that were not:
if(length(need_install)> 0) {
  install.packages(need_install)
}

# Load packages
invisible(lapply(required_packages, require, character.only=TRUE))

# Set working directory and font
dir_qual <- dirname(rstudioapi::getSourceEditorContext()$path)
setwd(dir_qual)
windowsFonts(sans=windowsFont("Helvetica")) 

### Load data
VPT_df <- read.csv2("Data/VPT_aggregated_wide.csv")
JND_df <- read.csv2("Data/JND_aggregated_wide.csv")
Qual_df <- read.csv2("Data/Qualtrics_merged_coderesponses.csv")

colnames(VPT_df)[10:length(VPT_df)] <- paste0("VPT_", colnames(VPT_df[10:length(VPT_df)]))
colnames(JND_df)[2:length(JND_df)] <- paste0("JND_", colnames(JND_df[2:length(JND_df)]))

### Merge data
Total_df <- merge(Qual_df, VPT_df, by="ppn", all.x=TRUE)
Total_df <- merge(Total_df, JND_df, by="ppn", all.x=TRUE) 
Total_df <- Total_df %>% mutate_all(na_if,"") ## replace empty cells by NA

### Select relevant data

Total_df <- Total_df %>% select(ppn, Exclusion, Gender_bio, Gender_id, Age_test, Age_cat, Education, Autism,
                                Autism_year_1_1, ADHD, Psych, Psych_past, Neuro, Alcohol, Drugs, Meds, Sleep_qual_1,
                                Sleep_quant, AQ_10_sum, MMQCon_sum, MMQAb_sum, MMQStrat_sum, Strategy_VPT, Coded_specific,  
                                Coded_aggregated,Coded_summary, Strategy_JND,Concentration_VPT, Effort_VPT, Performance_VPT, Concentration_JND, 
                                Effort_JND, Performance_JND, VPT_span_traditional_high, VPT_span_traditional_low,
                                VPT_span_traditional_mean_high, VPT_span_traditional_mean_low,
                                VPT_mrt_wo_high, VPT_mrt_wo_low,VPT_mrt_cor_wo_high,VPT_mrt_cor_wo_low, VPT_acc_wo_high, 
                                VPT_acc_wo_low, VPT_n_outlier_high, VPT_n_outlier_low, VPT_n_FaR_high, VPT_n_FaR_low,
                                JND_mrt_wo_full, JND_acc_wo_full, JND_mrt_cor_wo_full,JND_score_wo_full,
                                JND_mrt_wo_final, JND_acc_wo_final,JND_mrt_cor_wo_final, JND_score_wo_final,
                                JND_n_outlier_full, JND_n_FaR_full, JND_n_min_full, JND_n_max_full)

# Add concatenated group variable

Total_df$Groups <- paste0(Total_df$Autism, "_", Total_df$Age_cat)
Total_df$Groups <- factor(Total_df$Groups, levels=c("autism_younger", "autism_older", "no-autism_younger", "no-autism_older"))
levels(Total_df$Groups) <- list(aut_yng = "autism_younger" , aut_old = "autism_older",
                                naut_yng = "no-autism_younger", naut_old = "no-autism_older")

Total_df$Age_cat <- factor(Total_df$Age_cat, levels=c("younger", "older"))
levels(Total_df$Age_cat) <- list(yng = "younger" , old = "older")

Total_df$Autism <- factor(Total_df$Autism, levels=c("autism", "no-autism"))
levels(Total_df$Autism) <- list(aut = "autism", naut= "no-autism")

Total_df$Gender_bio <- factor(Total_df$Gender_bio, levels = c("man", "woman", "other"))
Total_df$Edu_cat <- ifelse(Total_df$Education < 5, "low",
                           ifelse(Total_df$Education < 7, "middle", "high"))
Total_df$Edu_cat <- factor(Total_df$Edu_cat, levels = c("low", "middle", "high"))

Exclude participants

#general exclude

Total_df <- Total_df %>% filter(Exclusion == 0) #N exclusion = 4 (all based on autism family)

#VPT exclude
VPT_long_df <- Total_df %>% pivot_longer(c(VPT_span_traditional_high:VPT_n_FaR_low), 
                                         names_to = c(".value","type"),
                                         names_pattern = "(.*)_(.*)")

VPT_exclusion_df <- VPT_long_df %>% select(ppn, Exclusion, Age_cat, Autism, type, VPT_mrt_wo)

VPT_avg <- VPT_exclusion_df %>% 
  group_by(type, Age_cat, Autism) %>%
  summarise(mean_mrt = mean(VPT_mrt_wo, na.rm=TRUE),
            sd_mrt = sd(VPT_mrt_wo, na.rm=TRUE)) %>%
  ungroup()

VPT_exclusion_df <- merge(VPT_exclusion_df, VPT_avg, by=c("type", "Age_cat", "Autism"))
VPT_exclusion_df$Exclusion_VPT <- ifelse(VPT_exclusion_df$VPT_mrt_wo > (VPT_exclusion_df$mean_mrt+(3*VPT_exclusion_df$sd_mrt)), 1,
                                         ifelse(VPT_exclusion_df$VPT_mrt_wo < (VPT_exclusion_df$mean_mrt-(3*VPT_exclusion_df$sd_mrt)),1,0))
VPT_exclusion_df <- pivot_wider(VPT_exclusion_df, names_from = type, names_sep= "_", values_from = c(VPT_mrt_wo:Exclusion_VPT))
VPT_exclusion_df$Exclusion_VPT <- ifelse(VPT_exclusion_df$Exclusion_VPT_high == 1 | VPT_exclusion_df$Exclusion_VPT_low == 1, 1, 0) 
VPT_exclusion_df <- select(VPT_exclusion_df, c(ppn, Exclusion_VPT))
Total_df <- merge(Total_df, VPT_exclusion_df, by = "ppn") #add to total_df
Total_df$Exclusion_VPT <- ifelse(is.na(Total_df$Exclusion_VPT), -9, Total_df$Exclusion_VPT) #set NA to -9, because then no VPT data was present
VPT_long_df <- merge(VPT_long_df, VPT_exclusion_df, by = "ppn") 
VPT_long_df$Exclusion_VPT <- ifelse(is.na(VPT_long_df$Exclusion_VPT), -9, VPT_long_df$Exclusion_VPT) # here too.

#JND exclude
JND_long_df <- Total_df %>% pivot_longer(c(JND_mrt_wo_full:JND_score_wo_final),
                                         names_to = c(".value", "dataset"),
                                         names_pattern = "(.*)_(.*)")

JND_exclusion_df <- JND_long_df %>% select(ppn, Exclusion, Age_cat, Autism, dataset, JND_mrt_wo)

JND_avg <- JND_exclusion_df %>% 
  group_by(dataset, Age_cat, Autism) %>%
  summarise(mean_mrt = mean(JND_mrt_wo, na.rm=TRUE),
            sd_mrt = sd(JND_mrt_wo, na.rm=TRUE)) %>%
  ungroup()

JND_exclusion_df <- merge(JND_exclusion_df, JND_avg, by=c("dataset", "Age_cat", "Autism"))
JND_exclusion_df$Exclusion_JND <- ifelse(JND_exclusion_df$JND_mrt_wo > (JND_exclusion_df$mean_mrt+(3*JND_exclusion_df$sd_mrt)), 1,
                                         ifelse(JND_exclusion_df$JND_mrt_wo < (JND_exclusion_df$mean_mrt-(3*JND_exclusion_df$sd_mrt)),1,0))

JND_exclusion_df <- pivot_wider(JND_exclusion_df, names_from = dataset, names_sep= "_", values_from = c(JND_mrt_wo:Exclusion_JND))
## even the outlier responder (N=1) is really fast (<2 sec) with no indication of taking brakes, so are not removed.
JND_exclusion_df$Exclusion_JND <- ifelse(JND_exclusion_df$Exclusion_JND_full, 0, 0) ## no one is excluded.
JND_exclusion_df <- select(JND_exclusion_df, c(ppn, Exclusion_JND))
Total_df <- merge(Total_df, JND_exclusion_df, by = "ppn") #add to total_df
Total_df$Exclusion_JND <- ifelse(is.na(Total_df$Exclusion_JND), -9, Total_df$Exclusion_JND)
JND_long_df <- merge(JND_long_df, JND_exclusion_df, by = "ppn")
JND_long_df$Exclusion_JND <- ifelse(is.na(JND_long_df$Exclusion_JND), -9, JND_long_df$Exclusion_JND)

write_csv2(Total_df, "Datasets_merged.csv") # save total dataset with excluded ppns indicated.

Descriptives (Aut/nAut)

We create basic information tables for N and means/sds for all relevant variables.

get_N = function(data, var, title){
  
  Ns <- list()
  
  for (i in 1:length(var)){
    
  var1 = sym(var[i])
  groups = c("Autism", "Age_cat", var[i])
  
  N_df <- data %>% #getting some info on categorical data
  group_by(.dots = groups) %>% 
  summarise(N = n()) %>%
  pivot_wider(names_from = c("Age_cat", "Autism"), values_from = N)
  colnames(N_df)[1] <- "category"

  Chidf <- data %>%
           group_by(Age_cat) %>%
           summarise(Chi = chisq.test(Autism, !!var1)$statistic,
                     p = chisq.test(Autism, !!var1)$p.value) %>%
           mutate_if(is.numeric, round, 3)  %>%
           pivot_wider(names_from = "Age_cat", values_from = c(Chi, p))
  
  Var <- var[i]
  
  N_df <- cbind(N_df, Chidf)
  N_df <- cbind(Var, N_df)
  Ns[[i]] <- N_df
  
  }
  
  finalNs <- do.call(rbind,Ns)
  
return(kbl(finalNs, caption=paste0(title), digits = 3,booktabs = TRUE) %>%  
       kable_classic(full_width = F, font_size = 11) %>%
       kable_styling(latex_options = "HOLD_position"))
}

nvars <- c("Gender_bio","Edu_cat", "ADHD", "Psych", "Psych_past", "Neuro", "Drugs", "Alcohol")
get_N(Total_df, nvars, title="N demographs")

N demographs

Var	category	yng_aut	old_aut	yng_naut	old_naut	Chi_yng	Chi_old	p_yng	p_old
Gender_bio	man	9	31	11	21	0.127	1.790	0.721	0.409
Gender_bio	woman	23	15	20	16	0.127	1.790	0.721	0.409
Gender_bio	other		1			0.127	1.790	0.721	0.409
Edu_cat	middle	22	30	15	19	3.309	1.635	0.191	0.441
Edu_cat	high	10	14	15	16	3.309	1.635	0.191	0.441
Edu_cat	low		3	1	2	3.309	1.635	0.191	0.441
ADHD	ADHD	4	12	4		0.000	9.035	1.000	0.003
ADHD	no-ADHD	28	35	27	37	0.000	9.035	1.000	0.003
Psych	no	26	39	27	37	0.084	5.126	0.772	0.024
Psych	yes	6	8	4		0.084	5.126	0.772	0.024
Psych_past	no	17	34	28	34	8.931	3.943	0.003	0.047
Psych_past	yes	15	13	3	3	8.931	3.943	0.003	0.047
Neuro	no	28	42	31	34	2.303	0.000	0.129	0.986
Neuro	yes	4	5		3	2.303	0.000	0.129	0.986
Drugs	no	31	47	31	36	0.000	0.015	1.000	0.904
Drugs	yes	1			1	0.000	0.015	1.000	0.904
Alcohol	no	31	34	24	25	3.765	0.055	0.052	0.814
Alcohol	yes	1	13	7	12	3.765	0.055	0.052	0.814

nvarsVPT <- c("Gender_bio","Edu_cat", "ADHD", "Psych", "Psych_past", "Neuro", "Alcohol")
get_N(subset(Total_df, Exclusion_VPT == 0), nvarsVPT, title= "N demographics VPT sample")

N demographics VPT sample

Var	category	yng_aut	old_aut	yng_naut	old_naut	Chi_yng	Chi_old	p_yng	p_old
Gender_bio	man	6	26	10	20	0.517	1.602	0.472	0.449
Gender_bio	woman	21	12	19	14	0.517	1.602	0.472	0.449
Gender_bio	other		1			0.517	1.602	0.472	0.449
Edu_cat	middle	18	25	13	16	3.239	2.772	0.198	0.250
Edu_cat	high	9	11	15	16	3.239	2.772	0.198	0.250
Edu_cat	low		3	1	2	3.239	2.772	0.198	0.250
ADHD	ADHD	3	12	4		0.000	10.380	1.000	0.001
ADHD	no-ADHD	24	27	25	34	0.000	10.380	1.000	0.001
Psych	no	22	34	25	34	0.014	2.886	0.907	0.089
Psych	yes	5	5	4		0.014	2.886	0.907	0.089
Psych_past	no	14	29	26	31	8.026	2.455	0.005	0.117
Psych_past	yes	13	10	3	3	8.026	2.455	0.005	0.117
Neuro	no	23	35	29	31	2.663	0.000	0.103	1.000
Neuro	yes	4	4		3	2.663	0.000	0.103	1.000
Alcohol	no	26	27	22	23	3.245	0.000	0.072	1.000
Alcohol	yes	1	12	7	11	3.245	0.000	0.072	1.000

get_N(subset(Total_df, Exclusion_JND == 0), c("Gender_bio"), title= "N (sex) JND sample")

N (sex) JND sample

Var	category	yng_aut	old_aut	yng_naut	old_naut	Chi_yng	Chi_old	p_yng	p_old
Gender_bio	man	5	25	10	20	0.935	2.167	0.334	0.338
Gender_bio	woman	22	10	20	14	0.935	2.167	0.334	0.338
Gender_bio	other		1			0.935	2.167	0.334	0.338

get_chi_agg <- function(data, var, title){
  
  Chi <- list()
  
  for (i in 1:length(var)){
    
  var1 = sym(var[i])
  
  Chidf <- data %>%
           summarise(Chi_autism = chisq.test(Autism,!!var1)$statistic,
                     p_autism = chisq.test(Autism,!!var1)$p.value,
                     Chi_age = chisq.test(Age_cat,!!var1)$statistic,
                     p_age = chisq.test(Age_cat,!!var1)$p.value) %>%
           mutate_if(is.numeric, round, 3) %>%
           pivot_longer(everything(),names_to = c("statistic", "predictor"), names_pattern = "(.*)_(.*)")
  Chidf$outcome <- var[i]
  
  Chi[[i]] <- Chidf
  }
  
  finalchi <- do.call(rbind, Chi)
  
  return(kbl(finalchi, caption=paste0(title), digits = 3,booktabs = TRUE) %>%  
       kable_classic(full_width = F, font_size = 11) %>%
       kable_styling(latex_options = "HOLD_position"))
}

get_chi_agg(Total_df, nvars, title="main effects demographics")

main effects demographics

statistic	predictor	value	outcome
Chi	autism	1.126	Gender_bio
p	autism	0.569	Gender_bio
Chi	age	14.463	Gender_bio
p	age	0.001	Gender_bio
Chi	autism	3.857	Edu_cat
p	autism	0.145	Edu_cat
Chi	age	1.833	Edu_cat
p	age	0.400	Edu_cat
Chi	autism	5.256	ADHD
p	autism	0.022	ADHD
Chi	age	0.001	ADHD
p	age	0.972	ADHD
Chi	autism	3.729	Psych
p	autism	0.053	Psych
Chi	age	0.824	Psych
p	age	0.364	Psych
Chi	autism	13.106	Psych_past
p	autism	0.000	Psych_past
Chi	age	1.340	Psych_past
p	age	0.247	Psych_past
Chi	autism	1.535	Neuro
p	autism	0.215	Neuro
Chi	age	0.153	Neuro
p	age	0.696	Neuro
Chi	autism	0.000	Drugs
p	autism	1.000	Drugs
Chi	age	0.000	Drugs
p	age	1.000	Drugs
Chi	autism	1.645	Alcohol
p	autism	0.200	Alcohol
Chi	age	5.081	Alcohol
p	age	0.024	Alcohol

get_chi_agg(subset(Total_df, Exclusion_VPT == 0), nvarsVPT, title= "main effects demographics VPT sample")

main effects demographics VPT sample

statistic	predictor	value	outcome
Chi	autism	0.995	Gender_bio
p	autism	0.608	Gender_bio
Chi	age	16.533	Gender_bio
p	age	0.000	Gender_bio
Chi	autism	5.028	Edu_cat
p	autism	0.081	Edu_cat
Chi	age	2.027	Edu_cat
p	age	0.363	Edu_cat
Chi	autism	5.642	ADHD
p	autism	0.018	ADHD
Chi	age	0.141	ADHD
p	age	0.708	ADHD
Chi	autism	1.752	Psych
p	autism	0.186	Psych
Chi	age	1.914	Psych
p	age	0.167	Psych
Chi	autism	10.453	Psych_past
p	autism	0.001	Psych_past
Chi	age	1.534	Psych_past
p	age	0.215	Psych_past
Chi	autism	1.394	Neuro
p	autism	0.238	Neuro
Chi	age	0.031	Neuro
p	age	0.861	Neuro
Chi	autism	0.947	Alcohol
p	autism	0.331	Alcohol
Chi	age	4.248	Alcohol
p	age	0.039	Alcohol

Descriptives <- c("Age_test", "AQ_10_sum", "Sleep_qual_1", "Sleep_quant")
Describe_long <- Total_df %>% 
  pivot_longer(cols = all_of(Descriptives), names_to ="var", values_to = "value") 

get_ss_df = function (data,title) {
  ss_df <- data %>%
    group_by(var, Age_cat, Autism)%>%
    summarise(mean = round(mean(value, na.rm=TRUE),2),
              sd = round(sd(value, na.rm=TRUE),2),
              min = round(min(value, na.rm=TRUE),2), 
              max = round(max(value, na.rm=TRUE),2),
              final = paste0(mean, ",", sd, "(", min, "-", max, ")")) %>%
    pivot_wider(names_from = c(Autism, Age_cat), values_from = c(mean, sd, min,max, final))
  
  ss_df <- ss_df[,c(1,(length(ss_df)-3):length(ss_df))]
  
  tval <- data %>%
    group_by(var, Age_cat) %>%
    summarise(tval = t.test(value[Autism=="aut"], value[Autism=="naut"])$statistic,
              pval = t.test(value[Autism=="aut"], value[Autism=="naut"])$p.value,
              n1 = sum(Autism=="aut"),
              n2 = sum(Autism=="naut"),
              dval = cohen.d(value[Autism=="aut"], value[Autism=="naut"], na.rm=TRUE)$estimate,
              BF = exp(ttest.tstat(t = tval, n1 = n1, n2=n2, rscale = 0.707)$bf)) %>%
    select(!c(n1,n2)) %>%
    mutate_if(is.numeric, round, 3) %>%
    pivot_wider(names_from = Age_cat,values_from = c(tval, pval, dval, BF))
  
  ss_df <- merge(ss_df, tval, by = "var")
  ss_df <- ss_df %>% relocate(c(tval_yng, pval_yng, dval_yng, BF_yng), .after=`final_naut_yng`)
  
  names(ss_df) <- sub('final_', '', names(ss_df))
  
  return(kbl(ss_df, caption = paste0(title), digits = 3, booktabs = TRUE, row.names=FALSE) 
         %>% kable_classic(full_width = F, font_size = 13) %>%
           kable_styling(latex_options = "HOLD_position"))
}

get_ss_df(Describe_long, "Descriptives per group") 

Descriptives per group

var	aut_yng	naut_yng	tval_yng	pval_yng	dval_yng	BF_yng	aut_old	naut_old	tval_old	pval_old	dval_old	BF_old
Age_test	36.63,3.69(24.25-40.82)	32.45,4.57(25.06-40.6)	3.987	0.000	1.008	1.325810e+02	71.83,4.67(65.19-82.23)	73.69,5.07(65.13-87.75)	-1.736	0.087	-0.385	8.430000e-01
AQ_10_sum	6.84,1.95(1-9)	2.52,1.79(0-6)	9.182	0.000	2.311	1.231554e+10	7,1.98(2-10)	2.32,1.62(0-6)	11.850	0.000	2.560	1.694087e+16
Sleep_qual_1	5.22,2.6(0-10)	6.32,1.89(3-10)	-1.933	0.058	-0.485	1.221000e+00	6.61,2.14(0-10)	7.03,1.62(3-10)	-1.011	0.315	-0.217	3.570000e-01
Sleep_quant	7.03,2.19(1-10)	6.84,1.04(5-9)	0.448	0.656	0.112	2.800000e-01	7.07,1.68(2-10)	7.11,1.02(5-10)	-0.137	0.891	-0.029	2.310000e-01

get_ss_df(subset(Describe_long, Exclusion_VPT==0), "Descriptives per group VPT sample") 

Descriptives per group VPT sample

var	aut_yng	naut_yng	tval_yng	pval_yng	dval_yng	BF_yng	aut_old	naut_old	tval_old	pval_old	dval_old	BF_old
Age_test	36.99,3.06(28.87-40.82)	32.66,4.59(25.06-40.6)	4.175	0.000	1.101	2.048990e+02	71.72,4.77(65.19-82.23)	73.96,5.12(65.13-87.75)	-1.924	0.059	-0.454	1.1660e+00
AQ_10_sum	6.74,1.99(1-9)	2.48,1.82(0-6)	8.322	0.000	2.233	2.087548e+08	7.26,1.95(2-10)	2.24,1.63(0-6)	11.881	0.000	2.777	2.1491e+15
Sleep_qual_1	5.07,2.56(0-8)	6.34,1.88(3-10)	-2.108	0.040	-0.570	1.659000e+00	6.39,2.3(0-10)	7.24,1.5(3-10)	-1.857	0.068	-0.429	1.0500e+00
Sleep_quant	7.15,2.13(1-10)	6.86,0.99(5-9)	0.638	0.527	0.175	3.200000e-01	6.89,1.74(2-10)	7.18,1(5-10)	-0.854	0.396	-0.196	3.3100e-01

### anova function
get_anova_df <- function(Data, Variables, Title){
  
  AOV_sum_list <- list()
  
  for (i in 1:length(Variables)){
    
    Describe_temp <- Data[Data['var']==Variables[i],]
    
    AOV_descr <- aov(value ~ Age_cat*Autism, data=Describe_temp)
    AOV_sum <- data.frame(matrix(ncol = 0, nrow = 3))
    rownames(AOV_sum) <- c("Age_cat", "Autism", "Age_cat*Autism")
    AOV_sum$Fval <- summary(AOV_descr)[[1]][["F value"]][1:3]
    AOV_sum$Pval <- summary(AOV_descr)[[1]][["Pr(>F)"]][1:3]
    AOV_sum$Eta <- (summary(AOV_descr)[[1]][["Sum Sq"]]/sum(summary(AOV_descr)[[1]][["Sum Sq"]]))[1:3]
    AOV_sum <- round(AOV_sum,3)
    
    Predictors <- rownames(AOV_sum)
    AOV_sum <- cbind(Predictors, AOV_sum)

    Outcome <- Variables[i] #AOV_sum <- AOV_sum %>% pivot_wider(names_from ="Predictors", values_from = c(Fval, Pval, Eta))
    AOV_sum <- cbind(Outcome, AOV_sum)
    
    AOV_sum_list[[i]] <- AOV_sum
    
  }
  
  Anovas <- do.call(rbind, AOV_sum_list)
  return(kbl(Anovas, caption = paste0(Title), digits = 3, booktabs = TRUE, row.names=FALSE)  %>% kable_classic(full_width = F, font_size = 13) %>%
           kable_styling(latex_options = "HOLD_position"))
}

get_anova_df(Describe_long, Descriptives, "Anova Descriptives")

Anova Descriptives

Outcome	Predictors	Fval	Pval	Eta
Age_test	Age_cat	2507.117	0.000	0.940
Age_test	Autism	0.973	0.326	0.000
Age_test	Age_cat*Autism	15.721	0.000	0.006
AQ_10_sum	Age_cat	0.426	0.515	0.001
AQ_10_sum	Autism	217.698	0.000	0.604
AQ_10_sum	Age_cat*Autism	0.317	0.575	0.001
Sleep_qual_1	Age_cat	8.784	0.004	0.056
Sleep_qual_1	Autism	4.270	0.041	0.027
Sleep_qual_1	Age_cat*Autism	0.961	0.328	0.006
Sleep_quant	Age_cat	0.324	0.570	0.002
Sleep_quant	Autism	0.055	0.816	0.000
Sleep_quant	Age_cat*Autism	0.199	0.656	0.001

get_anova_df(subset(Describe_long, Exclusion_VPT ==0), Descriptives, "Anova Descriptives VPT sample")

Anova Descriptives VPT sample

Outcome	Predictors	Fval	Pval	Eta
Age_test	Age_cat	2236.429	0.000	0.940
Age_test	Autism	0.594	0.442	0.000
Age_test	Age_cat*Autism	16.601	0.000	0.007
AQ_10_sum	Age_cat	1.144	0.287	0.003
AQ_10_sum	Autism	204.517	0.000	0.618
AQ_10_sum	Age_cat*Autism	1.356	0.246	0.004
Sleep_qual_1	Age_cat	8.156	0.005	0.058
Sleep_qual_1	Autism	7.797	0.006	0.056
Sleep_qual_1	Age_cat*Autism	0.335	0.564	0.002
Sleep_quant	Age_cat	0.010	0.919	0.000
Sleep_quant	Autism	0.015	0.903	0.000
Sleep_quant	Age_cat*Autism	1.086	0.299	0.009

MMQ

Descriptives and analyses for MMQ are provided below. Including figures.

### some figures and prelim analyses MMQ
MMQVars <- c("MMQCon_sum", "MMQAb_sum", "MMQStrat_sum")
MMQ_long <- Total_df %>% 
  pivot_longer(cols = all_of(MMQVars), names_to ="var", values_to = "value") 

get_ss_df(MMQ_long, "MMQ per group")

MMQ per group

var	aut_yng	naut_yng	tval_yng	pval_yng	dval_yng	BF_yng	aut_old	naut_old	tval_old	pval_old	dval_old	BF_old
MMQAb_sum	55.22,12.67(30-78)	64.1,8.14(37-74)	-3.320	0.002	-0.831	21.685	56.3,13(12-73)	63.49,7.57(43-76)	-3.142	0.002	-0.657	14.720
MMQCon_sum	48.19,13.16(23-69)	55.65,9.31(34-71)	-2.603	0.012	-0.653	4.162	44.04,15.02(9-72)	52.73,9.54(32-72)	-3.202	0.002	-0.675	17.190
MMQStrat_sum	27.47,12.37(6-66)	24.58,9.66(7-45)	1.034	0.305	0.260	0.404	24.04,10.39(6-47)	19.41,7.14(7-34)	2.404	0.019	0.510	2.721

get_anova_df(MMQ_long, MMQVars, "Anova MMQ")

Anova MMQ

Outcome	Predictors	Fval	Pval	Eta
MMQCon_sum	Age_cat	3.694	0.057	0.023
MMQCon_sum	Autism	15.996	0.000	0.099
MMQCon_sum	Age_cat*Autism	0.089	0.766	0.001
MMQAb_sum	Age_cat	0.002	0.964	0.000
MMQAb_sum	Autism	19.402	0.000	0.120
MMQAb_sum	Age_cat*Autism	0.219	0.641	0.001
MMQStrat_sum	Age_cat	5.923	0.016	0.039
MMQStrat_sum	Autism	5.439	0.021	0.035
MMQStrat_sum	Age_cat*Autism	0.272	0.603	0.002

# and for only those with VPT data:
get_ss_df(subset(MMQ_long, Exclusion_VPT==0), "MMQ VPT sample per group")

MMQ VPT sample per group

var	aut_yng	naut_yng	tval_yng	pval_yng	dval_yng	BF_yng	aut_old	naut_old	tval_old	pval_old	dval_old	BF_old
MMQAb_sum	57.19,11.33(35-78)	65.1,6.6(51-74)	-3.166	0.003	-0.862	14.332	55.89,13.38(12-73)	63.68,7.53(43-76)	-3.081	0.003	-0.706	12.401
MMQCon_sum	50.48,12.35(23-69)	56.34,9.2(34-71)	-2.004	0.051	-0.541	1.397	43,14.82(9-72)	52.62,9.54(32-72)	-3.308	0.002	-0.763	21.993
MMQStrat_sum	25,9.98(6-47)	24.48,9.93(7-45)	0.194	0.847	0.052	0.274	25.05,10.72(6-47)	19.24,6.92(7-34)	2.762	0.007	0.637	5.873

get_anova_df(subset(MMQ_long, Exclusion_VPT==0), MMQVars, "Anova MMQ VPT sample")

Anova MMQ VPT sample

Outcome	Predictors	Fval	Pval	Eta
MMQCon_sum	Age_cat	7.992	0.005	0.054
MMQCon_sum	Autism	14.419	0.000	0.098
MMQCon_sum	Age_cat*Autism	0.787	0.377	0.005
MMQAb_sum	Age_cat	0.882	0.350	0.006
MMQAb_sum	Autism	18.651	0.000	0.130
MMQAb_sum	Age_cat*Autism	0.001	0.970	0.000
MMQStrat_sum	Age_cat	2.056	0.154	0.015
MMQStrat_sum	Autism	4.325	0.040	0.033
MMQStrat_sum	Age_cat*Autism	2.446	0.120	0.018

plotting <- function(data, y,lims,ylab,title){
  plotting <- ggplot(data, aes(x=Groups, y=y))+
    stat_boxplot(geom = "errorbar", width = 0.05) + 
    geom_boxplot(size=.75, width=.3, outlier.alpha = 0)+
    geom_beeswarm(cex = 1.2, size=1,alpha=.6, aes(color=Groups))+
    stat_summary(fun=mean, aes(group=1), size=2, geom="point", color="darkgrey", alpha=.7) +
    scale_color_manual(name="Groups",values = c("#059D05", "#07D607", 
                                                "#9933FF", "#BF80FF"))+
    theme_classic(base_size = 7)+
    coord_cartesian(ylim = lims)+
    scale_x_discrete(labels =  c("aut_old" = "AutOld",
                                 "aut_yng" = "AutYoung",
                                 "naut_old" = "nAutOld",
                                 "naut_yng" = "nAutYoung"))+
    ylab(paste0(ylab))+
    xlab("")+
    theme(text=element_text(family="sans", size=12))+
    ggtitle(" ")+
    theme(legend.position = "bottom",legend.justification = "bottom", legend.title = element_blank())
  
  ggsave(filename=paste0(title), plot = plotting, device = NULL, path = NULL,
         scale = 1, units = c("in"), width=5, height =6, dpi = 600, limitsize = TRUE)
  plotting
}

plotting(Total_df, y=Total_df$MMQCon_sum,lims = c(0,72),ylab="MMQ Contentment", title="MMQ_contentment.png")

plotting(Total_df, y=Total_df$MMQAb_sum,lims = c(0,80),ylab="MMQ Ability", title="MMQ_ability.png")

plotting(Total_df, y=Total_df$MMQStrat_sum,lims = c(0,76),ylab="MMQ Strategies", title="MMQ_strategy.png")

# only vpt
plotting(subset(Total_df, Exclusion_VPT ==0), y=subset(Total_df, Exclusion_VPT ==0)$MMQStrat_sum,lims = c(0,76),ylab="MMQ Strategies", title="MMQ_strategy_VPTSample.png")

## VPT

Descriptives and analyses for VPT are provided below. Including figures. New functions are needed due to the extra “within” level for these data. Motivation stats/checks can be done using generic functions.

VPT_long_df <- subset(VPT_long_df, Exclusion_VPT==0)
VPT_vars <- c("VPT_span_traditional", "VPT_span_traditional_mean", "VPT_acc_wo", "VPT_mrt_cor_wo")
VPT_analyses_df <-  VPT_long_df %>% 
  pivot_longer(cols = all_of(VPT_vars), names_to ="vars", values_to = "value") 

get_vpt_descriptives <- function(data, title) {
VPT_df <- data %>%
    group_by(vars,type, Age_cat, Autism)%>%
    summarise(mean = round(mean(value, na.rm=TRUE),2),
              sd = round(sd(value, na.rm=TRUE),2),
              min = round(min(value, na.rm=TRUE),2), 
              max = round(max(value, na.rm=TRUE),2),
              final = paste0(mean, ",", sd, "(", min, "-", max, ")")) %>%
    pivot_wider(names_from = type, values_from = c(mean, sd, min,max, final))
  
  VPT_df <- VPT_df[,c(1:3,(length(VPT_df)-1):length(VPT_df))]
  
  Pairs_tot <- VPT_analyses_df %>%
  group_by(vars, Autism, Age_cat) %>%
  summarise(tval = t.test(value[type=="high"], value[type=="low"], paired=TRUE)$statistic,
            pval = t.test(value[type=="high"], value[type=="low"], paired=TRUE)$p.value,
            dval = cohen.d(value[type=="high"], value[type=="low"], na.rm=TRUE, paired=TRUE)$estimate,
            n1 = n()/2,
            BF = exp(ttest.tstat(t = tval, n1 = n1, rscale = 0.707)$bf)) %>%
  select(!n1) %>%
  mutate_if(is.numeric, round, 3)
  
  VPT_df <- merge(VPT_df, Pairs_tot, by = c("vars","Age_cat", "Autism"))

  names(VPT_df) <- sub('final_', '', names(VPT_df))
  
return(kbl(VPT_df, caption = paste0(title), digits = 3, booktabs = TRUE, row.names=FALSE)%>% kable_classic(full_width = F, font_size = 13) %>%
kable_styling(latex_options = "HOLD_position"))
}

get_anova_VPT <- function (data){
  
  AOV_VPT_three = get_anova_table(anova_test(data = data, dv = value, wid = ppn, between = c(Autism, Age_cat), within = type),
                                   correction="auto")
  AOV_VPT_three = AOV_VPT_three[AOV_VPT_three$Effect == "Autism:Age_cat:type",]     
  AOV_VPT_btw = get_anova_table(anova_test(data = data, dv = value, between = c(Autism, Age_cat)))
  AOV_VPT_btw = AOV_VPT_btw[AOV_VPT_btw$Effect == "Autism:Age_cat",]  
  AOV_VPT_wth1 = get_anova_table(anova_test(data = data, dv = value, wid = ppn, between = c(Age_cat), within = type),
                                  correction="auto")
  AOV_VPT_wth1 = AOV_VPT_wth1[AOV_VPT_wth1$Effect == "Age_cat:type",]
  AOV_VPT_wth2 = get_anova_table(anova_test(data = data, dv = value, wid = ppn, between = c(Autism), within = type),
                                  correction="auto")
  AOV_VPT_wth2 = AOV_VPT_wth2[AOV_VPT_wth2$Effect == "Autism:type",]
  
  AOV_tot = rbind(AOV_VPT_three, AOV_VPT_btw, AOV_VPT_wth1, AOV_VPT_wth2)
  
  AOV_VPT_Aut = get_anova_table(anova_test(data = data, dv = value, between = Autism))
  AOV_VPT_Aut = AOV_VPT_Aut[AOV_VPT_Aut$Effect == "Autism",] 
  
  AOV_VPT_Age = get_anova_table(anova_test(data = data, dv = value, between = Age_cat))
  AOV_VPT_Age = AOV_VPT_Age[AOV_VPT_Age$Effect == "Age_cat",] 
  
  AOV_VPT_type = get_anova_table(anova_test(data = data, dv = value, within = type, wid=ppn))
  AOV_VPT_type = AOV_VPT_type[AOV_VPT_type$Effect == "type",] 
  
  AOV_tot = rbind(AOV_tot, AOV_VPT_Aut, AOV_VPT_Age, AOV_VPT_type)

}

get_vpt_descriptives(VPT_analyses_df, "VPT differences high/low per group")

VPT differences high/low per group

vars	Age_cat	Autism	high	low	tval	pval	dval	BF
VPT_acc_wo	old	aut	0.74,0.12(0.52-0.95)	0.69,0.11(0.41-0.86)	2.669	0.011	0.471	3.754
VPT_acc_wo	old	naut	0.73,0.1(0.5-0.91)	0.68,0.09(0.5-0.86)	2.593	0.014	0.493	3.233
VPT_acc_wo	yng	aut	0.78,0.1(0.57-0.95)	0.72,0.13(0.38-0.91)	2.527	0.018	0.498	2.863
VPT_acc_wo	yng	naut	0.77,0.09(0.59-0.95)	0.76,0.12(0.5-1)	0.463	0.647	0.116	0.218
VPT_mrt_cor_wo	old	aut	3326.61,869.02(1776.74-5708.74)	3640.6,934.69(2086.05-6257.5)	-3.906	0.000	-0.345	74.444
VPT_mrt_cor_wo	old	naut	3050.84,780.79(1994.45-5763.42)	3481.81,1195.85(2004.57-8219.81)	-3.083	0.004	-0.384	9.247
VPT_mrt_cor_wo	yng	aut	2938.73,821.3(1875.02-4687.39)	3285.52,1046.41(1805.78-5525.51)	-2.714	0.012	-0.352	4.099
VPT_mrt_cor_wo	yng	naut	2887.9,692.5(1760.21-4712.96)	3243.48,879.94(2181.65-5484.79)	-2.870	0.008	-0.437	5.670
VPT_span_traditional	old	aut	12.44,3.09(6-15)	11.16,3.03(6-15)	2.266	0.029	0.448	1.671
VPT_span_traditional	old	naut	12.55,2.72(7-15)	10.32,2.74(6-15)	3.384	0.002	0.795	18.553
VPT_span_traditional	yng	aut	12.96,2.62(7-15)	11.85,2.97(7-15)	1.938	0.064	0.395	1.032
VPT_span_traditional	yng	naut	13.21,2.38(7-15)	12.45,3.12(7-15)	0.933	0.359	0.274	0.294
VPT_span_traditional_mean	old	aut	11.41,3.28(5.33-14.67)	9.99,3.1(5.33-14.67)	2.365	0.024	0.512	2.021
VPT_span_traditional_mean	old	naut	11.3,2.96(5.33-14.67)	9,2.87(5.33-14.67)	3.235	0.003	0.760	13.090
VPT_span_traditional_mean	yng	aut	11.83,3(5.33-14.67)	10.56,3.19(5.33-14.67)	2.013	0.055	0.410	1.163
VPT_span_traditional_mean	yng	naut	12.01,2.63(5.67-14.67)	11.39,3.42(5.67-14.67)	0.707	0.486	0.204	0.248

ANOVA_vpt <- VPT_analyses_df %>% group_by(vars) %>% do(data.frame(val=get_anova_VPT(.)))
kbl(ANOVA_vpt, caption = "VPT ANOVAs", digits = 3, booktabs = TRUE, row.names=FALSE )%>% 
kable_classic(full_width = F, font_size = 13) %>% kable_styling(latex_options = "HOLD_position")

VPT ANOVAs

vars	val.Effect	val.DFn	val.DFd	val.F	val.p	val.p..05	val.ges
VPT_acc_wo	Autism:Age_cat:type	1	125	0.833	0.363		0.002
VPT_acc_wo	Autism:Age_cat	1	254	0.859	0.355		0.003
VPT_acc_wo	Age_cat:type	1	127	0.500	0.481		0.001
VPT_acc_wo	Autism:type	1	127	1.279	0.260		0.003
VPT_acc_wo	Autism	1	256	0.056	0.813		0.000
VPT_acc_wo	Age_cat	1	256	11.827	0.001		0.044
VPT_acc_wo	type	1	128	16.043	0.000		0.038
VPT_mrt_cor_wo	Autism:Age_cat:type	1	125	0.209	0.649		0.000
VPT_mrt_cor_wo	Autism:Age_cat	1	254	0.535	0.465		0.002
VPT_mrt_cor_wo	Age_cat:type	1	127	0.021	0.884		0.000
VPT_mrt_cor_wo	Autism:type	1	127	0.350	0.555		0.000
VPT_mrt_cor_wo	Autism	1	256	1.827	0.178		0.007
VPT_mrt_cor_wo	Age_cat	1	256	6.371	0.012		0.024
VPT_mrt_cor_wo	type	1	128	38.652	0.000		0.037
VPT_span_traditional	Autism:Age_cat:type	1	123	0.759	0.385		0.003
VPT_span_traditional	Autism:Age_cat	1	252	1.192	0.276		0.005
VPT_span_traditional	Age_cat:type	1	125	1.511	0.221		0.005
VPT_span_traditional	Autism:type	1	125	0.147	0.702		0.000
VPT_span_traditional	Autism	1	254	0.003	0.960		0.000
VPT_span_traditional	Age_cat	1	254	7.354	0.007		0.028
VPT_span_traditional	type	1	126	17.520	0.000		0.055
VPT_span_traditional_mean	Autism:Age_cat:type	1	120	0.776	0.380		0.003
VPT_span_traditional_mean	Autism:Age_cat	1	249	1.799	0.181		0.007
VPT_span_traditional_mean	Age_cat:type	1	122	1.804	0.182		0.006
VPT_span_traditional_mean	Autism:type	1	122	0.000	0.989		0.000
VPT_span_traditional_mean	Autism	1	251	0.008	0.929		0.000
VPT_span_traditional_mean	Age_cat	1	251	6.390	0.012		0.025
VPT_span_traditional_mean	type	1	123	16.383	0.000		0.054

VPT_analyses_df$var <- paste0(VPT_analyses_df$vars, "_", VPT_analyses_df$type)
get_ss_df(VPT_analyses_df, "VPT group differences AUT/NAUT")

VPT group differences AUT/NAUT

var	aut_yng	naut_yng	tval_yng	pval_yng	dval_yng	BF_yng	aut_old	naut_old	tval_old	pval_old	dval_old	BF_old
VPT_acc_wo_high	0.78,0.1(0.57-0.95)	0.77,0.09(0.59-0.95)	0.305	0.762	0.082	0.281	0.74,0.12(0.52-0.95)	0.73,0.1(0.5-0.91)	0.525	0.601	0.122	0.273
VPT_acc_wo_low	0.72,0.13(0.38-0.91)	0.76,0.12(0.5-1)	-1.158	0.252	-0.311	0.471	0.69,0.11(0.41-0.86)	0.68,0.09(0.5-0.86)	0.306	0.761	0.071	0.252
VPT_mrt_cor_wo_high	2938.73,821.3(1875.02-4687.39)	2887.9,692.5(1760.21-4712.96)	0.249	0.804	0.067	0.277	3326.61,869.02(1776.74-5708.74)	3050.84,780.79(1994.45-5763.42)	1.428	0.158	0.333	0.579
VPT_mrt_cor_wo_low	3285.52,1046.41(1805.78-5525.51)	3243.48,879.94(2181.65-5484.79)	0.162	0.872	0.044	0.273	3640.6,934.69(2086.05-6257.5)	3481.81,1195.85(2004.57-8219.81)	0.625	0.534	0.149	0.287
VPT_span_traditional_high	12.96,2.62(7-15)	13.21,2.38(7-15)	-0.363	0.718	-0.098	0.285	12.44,3.09(6-15)	12.55,2.72(7-15)	-0.160	0.873	-0.037	0.245
VPT_span_traditional_low	11.85,2.97(7-15)	12.45,3.12(7-15)	-0.732	0.467	-0.196	0.338	11.16,3.03(6-15)	10.32,2.74(6-15)	1.228	0.224	0.288	0.462
VPT_span_traditional_mean_high	11.83,3(5.33-14.67)	12.01,2.63(5.67-14.67)	-0.244	0.808	-0.065	0.277	11.41,3.28(5.33-14.67)	11.3,2.96(5.33-14.67)	0.147	0.883	0.035	0.244
VPT_span_traditional_mean_low	10.56,3.19(5.33-14.67)	11.39,3.42(5.67-14.67)	-0.945	0.349	-0.252	0.392	9.99,3.1(5.33-14.67)	9,2.87(5.33-14.67)	1.387	0.170	0.331	0.552

# in figures:
GroupNames <- c(
  `aut_old` = "AutOld",
  `aut_yng` = "AutYoung",
  `naut_old` = "nAutOld",
  `naut_yng` = "nAutYoung")

plot_agg <- function(y,lims,ylab,title){
  plot_agg_wo <- ggplot(VPT_long_df, aes(x=type, y=y))+
    stat_boxplot(geom = "errorbar", width = 0.05) + 
    geom_boxplot(size=.75, width=.3, outlier.alpha = 0)+
    geom_beeswarm(cex = 1.2, size=1,alpha=.6, aes(color=type))+
    geom_line(aes(group=ppn),alpha=.7, size=.1, color="grey")+
    stat_summary(fun=mean, aes(group=1), size=1, geom="line", color="darkgrey", alpha=.7) +
    stat_summary(fun=mean, aes(group=1), size=2, geom="point", color="darkgrey", alpha=.7) +
    scale_color_manual(name="type",values = c("#059D05", "#9933FF"), labels=c("high","low"))+
    theme_classic(base_size = 7)+
    facet_wrap(~Groups, labeller=as_labeller(GroupNames))+
    coord_cartesian(ylim = lims)+
    ylab(paste0(ylab))+
    theme(text=element_text(family="sans", size=12))+
    ggtitle(" ")+
    theme(legend.position = "bottom",legend.justification = "bottom", legend.title = element_blank())
  
  ggsave(filename=paste0(title), plot = plot_agg_wo, device = NULL, path = NULL,
         scale = 1, units = c("in"), width=5, height =6, dpi = 600, limitsize = TRUE)
  
  plot_agg_wo
}

plot_agg(y=VPT_long_df$VPT_acc_wo,lims = c(0,1),ylab="accuracy", title="VPT_accuracy_pergroup.png")

plot_agg(y=VPT_long_df$VPT_span_traditional,lims = c(4,15),ylab="VPT span traditional", title="VPT_span_traditional_pergroup.png")

plot_agg(y=VPT_long_df$VPT_span_traditional_mean,lims = c(4,15),ylab="VPT mean span traditional", title="VPT_span_traditional_mean_pergroup.png")

plot_agg(y=VPT_long_df$VPT_mrt_cor_wo,lims = c(0,7500),ylab="mean correct RT", title="VPT_mrtcor_pergroup.png")

## Self-reported responses to the task:
VPTself_vars <- c("Effort_VPT", "Performance_VPT",  "Concentration_VPT")
VPTself_long <- Total_df %>% 
  pivot_longer(cols = all_of(VPTself_vars), names_to ="var", values_to = "value") 

get_ss_df(subset(VPTself_long, Exclusion_VPT==0), "VPT self reported performance")

VPT self reported performance

var	aut_yng	naut_yng	tval_yng	pval_yng	dval_yng	BF_yng	aut_old	naut_old	tval_old	pval_old	dval_old	BF_old
Concentration_VPT	7.65,1.39(4-10)	6.96,1.85(4-10)	1.455	0.153	0.411	0.648	7.34,1.66(4-10)	7.59,1.25(5-10)	-0.657	0.514	-0.168	0.291
Effort_VPT	8.24,0.94(7-10)	7.89,1.55(3-10)	0.962	0.341	0.264	0.397	8.47,1.44(5-10)	8.11,1.15(5-10)	1.060	0.294	0.272	0.392
Performance_VPT	6.38,1.63(3-9)	6.7,1.3(3-8)	-0.744	0.462	-0.223	0.340	6.41,1.36(3-9)	6.62,1.68(2-10)	-0.513	0.610	-0.138	0.271

get_anova_df(subset(VPTself_long, Exclusion_VPT==0), VPTself_vars, "ANOVA VPT self reported performance")

ANOVA VPT self reported performance

Outcome	Predictors	Fval	Pval	Eta
Effort_VPT	Age_cat	1.056	0.307	0.010
Effort_VPT	Autism	1.905	0.170	0.018
Effort_VPT	Age_cat*Autism	0.000	0.987	0.000
Performance_VPT	Age_cat	0.047	0.829	0.000
Performance_VPT	Autism	0.806	0.371	0.008
Performance_VPT	Age_cat*Autism	0.038	0.846	0.000
Concentration_VPT	Age_cat	0.435	0.511	0.004
Concentration_VPT	Autism	0.281	0.598	0.003
Concentration_VPT	Age_cat*Autism	2.289	0.133	0.022

## self reported strategy use after VPT:
vars = c("Strategy_VPT", "Coded_summary")
titles = c("Strategy_use", "Strategy_type_sum")
legends = c("Strategy use (y/n/?)", "Strategy type")
colors = list(c("#059D05","#7300E6","#9933FF"),c("#024e02","#036d03","#059D05","#36b036","#9933FF","#7300E6")) #better than palette 2.

Total_df$Coded_summary <- factor(Total_df$Coded_summary, levels=c("no", "?", "other", "mixed", "visual", "verbal"))

for (i in 1:length(vars)){
  
  data_temp <- Total_df[complete.cases(Total_df[,vars[i]]),]
  
  strategy_plot <- ggbarstats(
  data             = data_temp,
  x                = !!vars[i], # !! to unquote
  y                = Groups,
  label            = "both",
  title            = titles[i],
  xlab             = "",
  legend.title     = legends[i],
  palette          = "Set2"
)+
scale_fill_manual(values = colors[[i]])

ggsave(filename=paste0(titles[i], ".png"), plot = strategy_plot, device = NULL, path = NULL,
         scale = 1, units = c("in"), width=10, height =8, dpi = 600, limitsize = TRUE)
  
  strategy_plot

}

#simpler table in each age-group:
for (i in 1:length(vars)){
  
  data_temp <- subset(Total_df, Exclusion_VPT == 0 & !is.na(Total_df[,vars[i]]))
  
  strategy_plot <- ggbarstats(
  data             = data_temp,
  x                = !!vars[i], # !! to unquote
  y                = Groups,
  label            = "both",
  title            = titles[i],
  xlab             = "",
  legend.title     = legends[i],
  palette          = "Set2"
)+
scale_fill_manual(values = colors[[i]])

ggsave(filename=paste0(titles[i], ".png"), plot = strategy_plot, device = NULL, path = NULL,
         scale = 1, units = c("in"), width=10, height =8, dpi = 600, limitsize = TRUE)
  
  strategy_plot

}

#simpler table in each age-group:
get_N(subset(Total_df, Exclusion_VPT == 0 & !is.na(Total_df$Strategy_VPT)),"Strategy_VPT", title="N per strategy use")

N per strategy use

Var	category	yng_aut	old_aut	yng_naut	old_naut	Chi_yng	Chi_old	p_yng	p_old
Strategy_VPT	?	2	5	3	4	0.291	2.554	0.865	0.279
Strategy_VPT	no	6	8	5	14	0.291	2.554	0.865	0.279
Strategy_VPT	yes	19	18	19	13	0.291	2.554	0.865	0.279

get_N(subset(Total_df, Exclusion_VPT == 0 & !is.na(Total_df$Coded_aggregated)), "Coded_aggregated", title="N per strategy type")

N per strategy type

Var	category	yng_aut	old_aut	yng_naut	old_naut	Chi_yng	Chi_old	p_yng	p_old
Coded_aggregated	mixed	6	3	6	1	1.059	4.474	0.787	0.215
Coded_aggregated	other	1	1		4	1.059	4.474	0.787	0.215
Coded_aggregated	verbal	8	8	9	3	1.059	4.474	0.787	0.215
Coded_aggregated	visual	4	6	4	5	1.059	4.474	0.787	0.215

get_chi_agg(subset(Total_df, Exclusion_VPT == 0 & !is.na(Total_df$Strategy_VPT)), "Strategy_VPT", title="main effects strategy use")

main effects strategy use

statistic	predictor	value	outcome
Chi	autism	1.120	Strategy_VPT
p	autism	0.571	Strategy_VPT
Chi	age	4.992	Strategy_VPT
p	age	0.082	Strategy_VPT

get_chi_agg(subset(Total_df, Exclusion_VPT == 0 & !is.na(Total_df$Coded_aggregated)), "Coded_aggregated", title="main effects strategy type")

main effects strategy type

statistic	predictor	value	outcome
Chi	autism	1.185	Coded_aggregated
p	autism	0.757	Coded_aggregated
Chi	age	7.796	Coded_aggregated
p	age	0.050	Coded_aggregated

JND

JND_long_df <- subset(JND_long_df, Exclusion_JND == 0)
JND_vars <- c("JND_acc_wo", "JND_score_wo", "JND_mrt_cor_wo")
JND_analyses_df <- JND_long_df %>% pivot_longer(cols = all_of(JND_vars), names_to ="var", values_to = "value")

get_ss_df(subset(JND_analyses_df, dataset=="full"), "JND per group all data")

JND per group all data

var	aut_yng	naut_yng	tval_yng	pval_yng	dval_yng	BF_yng	aut_old	naut_old	tval_old	pval_old	dval_old	BF_old
JND_acc_wo	0.86,0.06(0.73-0.95)	0.86,0.06(0.73-1)	-0.549	0.585	-0.145	0.304	0.84,0.06(0.75-1)	0.83,0.07(0.58-1)	0.669	0.506	0.161	0.298
JND_mrt_cor_wo	819.42,292.94(444.78-1520)	783.29,278.92(380.64-1286.46)	0.476	0.636	0.127	0.295	923,310.64(493.31-1625)	887.89,403.71(333.52-1918.33)	0.406	0.686	0.098	0.264
JND_score_wo	1.16,0.03(1.12-1.26)	1.16,0.02(1.12-1.19)	0.232	0.818	0.062	0.274	1.16,0.03(1.13-1.28)	1.16,0.03(1.11-1.24)	0.994	0.324	0.237	0.375

get_anova_df(subset(JND_analyses_df, dataset=="full"), JND_vars, "JND per group all data")

JND per group all data

Outcome	Predictors	Fval	Pval	Eta
JND_acc_wo	Age_cat	4.266	0.041	0.033
JND_acc_wo	Autism	0.029	0.865	0.000
JND_acc_wo	Age_cat*Autism	0.730	0.394	0.006
JND_score_wo	Age_cat	0.055	0.814	0.000
JND_score_wo	Autism	0.855	0.357	0.007
JND_score_wo	Age_cat*Autism	0.293	0.589	0.002
JND_mrt_cor_wo	Age_cat	3.250	0.074	0.026
JND_mrt_cor_wo	Autism	0.372	0.543	0.003
JND_mrt_cor_wo	Age_cat*Autism	0.000	0.993	0.000

## but need bootstrapped results (for main analysis) because the data is quite skewed:
AOVboot <- ANOVA.boot(JND_score_wo_final ~ Autism*Age_cat, B = 1000, type = "residual", wild.dist = "normal",seed = NULL, data = subset(Total_df, Exclusion_JND == 0), keep.boot.resp = FALSE)
preds <- c(AOVboot$terms[1:3])
pvals <- c(AOVboot$`p-values`)
AOVboot <- data.frame(preds, pvals)
kbl(AOVboot, caption = "JND score bootstrapped p-values", digits = 3, booktabs = TRUE, row.names=FALSE )%>% 
kable_classic(full_width = F, font_size = 13) %>% kable_styling(latex_options = "HOLD_position")

JND score bootstrapped p-values

preds	pvals
Autism	0.751
Age_cat	0.674
Autism:Age_cat	0.515

plotting(subset(JND_long_df, dataset=="full"), y=subset(JND_long_df, dataset=="full")$JND_acc_wo,lims = c(0,1),ylab="JND accuracy full data", title="JND_acc_fulldata.png")

plotting(subset(JND_long_df, dataset=="full"), y=subset(JND_long_df, dataset=="full")$JND_score_wo,lims = c(1.01, 1.3),ylab="JND score full data", title="JND_score_fulldata.png")

plotting(subset(JND_long_df, dataset=="full"), y=subset(JND_long_df, dataset=="full")$JND_mrt_cor_wo,lims = c(0, 3000),ylab="JND correct RT full data", title="JND_mrt_fulldata.png")

get_ss_df(subset(JND_analyses_df, dataset=="final"), "JND per group final reversals set")

JND per group final reversals set

var	aut_yng	naut_yng	tval_yng	pval_yng	dval_yng	BF_yng	aut_old	naut_old	tval_old	pval_old	dval_old	BF_old
JND_acc_wo	0.72,0.13(0.33-0.92)	0.73,0.13(0.38-1)	-0.341	0.735	-0.090	0.282	0.7,0.12(0.43-1)	0.72,0.13(0.39-1)	-0.522	0.603	-0.125	0.277
JND_mrt_cor_wo	899.94,826.98(326-4338.8)	790.11,322.87(304.67-1360.83)	0.647	0.522	0.179	0.320	982.9,501.26(413.09-2421.8)	846.06,493.1(229.8-2722.97)	1.151	0.254	0.275	0.433
JND_score_wo	1.06,0.04(1.03-1.19)	1.06,0.03(1.02-1.15)	0.337	0.738	0.090	0.281	1.05,0.03(1.02-1.18)	1.06,0.04(1.02-1.26)	-0.663	0.510	-0.159	0.297

get_anova_df(subset(JND_analyses_df, dataset=="final"), JND_vars, "JND per group final reversals set")

JND per group final reversals set

Outcome	Predictors	Fval	Pval	Eta
JND_acc_wo	Age_cat	0.447	0.505	0.004
JND_acc_wo	Autism	0.377	0.540	0.003
JND_acc_wo	Age_cat*Autism	0.007	0.935	0.000
JND_score_wo	Age_cat	0.185	0.668	0.001
JND_score_wo	Autism	0.080	0.778	0.001
JND_score_wo	Age_cat*Autism	0.492	0.484	0.004
JND_mrt_cor_wo	Age_cat	0.567	0.453	0.005
JND_mrt_cor_wo	Autism	1.612	0.207	0.013
JND_mrt_cor_wo	Age_cat*Autism	0.019	0.891	0.000

plotting(subset(JND_long_df, dataset=="final"), y=subset(JND_long_df, dataset=="final")$JND_acc_wo,lims = c(0,1),ylab="JND accuracy final data", title="JND_acc_finaldata.png")

plotting(subset(JND_long_df, dataset=="final"), y=subset(JND_long_df, dataset=="final")$JND_score_wo,lims = c(1.01,1.3),ylab="JND score final data", title="JND_score_finaldata.png")

plotting(subset(JND_long_df, dataset=="final"), y=subset(JND_long_df, dataset=="final")$JND_mrt_cor_wo,lims = c(0, 3000),ylab="JND correct RT final data", title="JND_mrt_finaldata.png")

## Self-reported responses to the task:
JNDself_vars <- c("Effort_JND", "Performance_JND",  "Concentration_JND")
JNDself_long <- Total_df %>% 
  pivot_longer(cols = all_of(JNDself_vars), names_to ="var", values_to = "value") 

get_ss_df(subset(JNDself_long, Exclusion_JND==0), "JND self reported performance")

JND self reported performance

var	aut_yng	naut_yng	tval_yng	pval_yng	dval_yng	BF_yng	aut_old	naut_old	tval_old	pval_old	dval_old	BF_old
Concentration_JND	7.24,2.05(2-10)	6.85,1.92(4-10)	0.667	0.509	0.196	0.323	7.31,1.71(4-10)	7.75,1.29(5-10)	-1.063	0.293	-0.286	0.399
Effort_JND	8.14,1.01(6-10)	7.7,1.56(3-10)	1.175	0.246	0.325	0.477	8.21,1.5(5-10)	8,1.32(5-10)	0.535	0.595	0.146	0.278
Performance_JND	5.95,1.53(2-9)	6,1.67(1-9)	-0.102	0.919	-0.030	0.269	6.36,1.54(3-10)	6.54,1.67(2-10)	-0.411	0.683	-0.115	0.265

get_anova_df(subset(JNDself_long, Exclusion_JND==0), JNDself_vars, "ANOVA JND self reported performance")

ANOVA JND self reported performance

Outcome	Predictors	Fval	Pval	Eta
Effort_JND	Age_cat	0.617	0.434	0.006
Effort_JND	Autism	1.300	0.257	0.013
Effort_JND	Age_cat*Autism	0.174	0.677	0.002
Performance_JND	Age_cat	2.054	0.155	0.021
Performance_JND	Autism	0.136	0.713	0.001
Performance_JND	Age_cat*Autism	0.044	0.834	0.000
Concentration_JND	Age_cat	1.945	0.166	0.019
Concentration_JND	Autism	0.019	0.891	0.000
Concentration_JND	Age_cat*Autism	1.372	0.244	0.014

get_N(subset(Total_df, Exclusion_JND == 0 & !is.na(Total_df$Strategy_JND)),"Strategy_JND", title="N strategy use JND")

N strategy use JND

Var	category	yng_aut	old_aut	yng_naut	old_naut	Chi_yng	Chi_old	p_yng	p_old
Strategy_JND	?	4	2	2	4	1.62	1.481	0.445	0.477
Strategy_JND	no	14	24	22	18	1.62	1.481	0.445	0.477
Strategy_JND	yes	3	4	4	5	1.62	1.481	0.445	0.477

Exploratory

### relating performance to strategy use (MMQ/VPT)

cor.test(subset(Total_df,Autism =="naut" & Age_cat == "yng")$VPT_acc_wo_high, subset(Total_df, Autism=="naut" & Age_cat == "yng")$JND_score_wo_final)

## 
##  Pearson's product-moment correlation
## 
## data:  subset(Total_df, Autism == "naut" & Age_cat == "yng")$VPT_acc_wo_high and subset(Total_df, Autism == "naut" & Age_cat == "yng")$JND_score_wo_final
## t = -0.54816, df = 27, p-value = 0.5881
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  -0.4539617  0.2720552
## sample estimates:
##        cor 
## -0.1049111

Total_df$VPT_acc_diff <- Total_df$VPT_acc_wo_high - Total_df$VPT_acc_wo_low #higher acc_diff = more benefit from semantics

MMQstrat <- tidy(cor.test(Total_df$MMQStrat_sum, Total_df$VPT_acc_diff)) #broader indication (continuous) of strategy use (MMQ)
kbl(MMQstrat[1:6], caption = "MMQ strategies ~ acc diff", digits = 3, booktabs = TRUE, row.names=TRUE )%>% 
kable_classic(full_width = F, font_size = 13) %>% kable_styling(latex_options = "HOLD_position")

MMQ strategies ~ acc diff

	estimate	statistic	p.value	parameter	conf.low	conf.high
1	-0.018	-0.201	0.841	131	-0.187	0.153

VPTstrat <- summary(aov(VPT_acc_diff~Strategy_VPT, data = subset(Total_df, Exclusion_VPT ==0)))[[1]] ## narrow indication of strategy use (VPT)
kbl(VPTstrat, caption = "VPT strategies (y/n/?) ~ acc diff", digits = 3, booktabs = TRUE, row.names=TRUE )%>% 
kable_classic(full_width = F, font_size = 13) %>% kable_styling(latex_options = "HOLD_position")

VPT strategies (y/n/?) ~ acc diff

	Df	Sum Sq	Mean Sq	F value	Pr(>F)
Strategy_VPT	2	0.001	0.000	0.029	0.971
Residuals	113	1.737	0.015

### effort?
VPT_long_df$Effort_centered <- VPT_long_df$Effort_VPT - mean(VPT_long_df$Effort_VPT, na.rm=TRUE) ## center covariate (!)
ancova <- get_anova_table(anova_test(VPT_long_df, dv = VPT_acc_wo, wid = ppn, between = c(Age_cat, Autism), covariate = Effort_centered, within = type),correction="auto")
kbl(ancova,caption = "effort ancova", digits = 3, booktabs = TRUE)%>% 
kable_classic(full_width = F, font_size = 13) %>% kable_styling(latex_options = "HOLD_position") 

effort ancova

Effect	DFn	DFd	F	p	p<.05	ges
Effort_centered	1	102	0.197	0.658		0.001
Age_cat	1	102	10.379	0.002		0.063
Autism	1	102	0.001	0.981		0.000
type	1	102	8.008	0.006		0.026
Age_cat:Autism	1	102	0.227	0.634		0.001
Effort_centered:type	1	102	0.147	0.702		0.000
Age_cat:type	1	102	2.151	0.146		0.007
Autism:type	1	102	0.155	0.695		0.001
Age_cat:Autism:type	1	102	0.031	0.860		0.000

# homogeneity of regression slopes seems assumed...:
seps <- VPT_long_df %>% group_by(Groups) %>% anova_test(dv = VPT_acc_wo, wid = ppn, covariate = Effort_centered, within = type)
kbl(seps,caption = "hmgn regression slopes ancova", digits = 3, booktabs = TRUE)%>% 
kable_classic(full_width = F, font_size = 13) %>% kable_styling(latex_options = "HOLD_position") 

hmgn regression slopes ancova

Groups	Effect	DFn	DFd	F	p	p<.05	ges
aut_yng	Effort_centered	1	19	0.014	0.908		0.001
aut_yng	type	1	19	1.140	0.299		0.015
aut_yng	Effort_centered:type	1	19	0.066	0.800		0.001
aut_old	Effort_centered	1	30	0.279	0.602		0.006
aut_old	type	1	30	5.924	0.021		0.055
aut_old	Effort_centered:type	1	30	0.116	0.736		0.001
naut_yng	Effort_centered	1	25	0.049	0.827		0.001
naut_yng	type	1	25	0.084	0.774		0.002
naut_yng	Effort_centered:type	1	25	0.851	0.365		0.015
naut_old	Effort_centered	1	25	0.189	0.667		0.005
naut_old	type	1	25	5.975	0.022		0.072
naut_old	Effort_centered:type	1	25	0.136	0.716		0.002

#plot
effortplot <- ggplot(VPT_long_df, aes(x=VPT_acc_wo, y = Effort_VPT, grp=type))+ 
  geom_smooth(method = "lm", se=TRUE, alpha = .3, size = 1.5, aes(x=VPT_acc_wo, color=type, linetype=type))+
  geom_point(alpha = 1, size = 1.5, aes(x=VPT_acc_wo, color=type))+
  scale_linetype_manual(values=c(1,2,1,2)) +
  scale_color_manual(name="type", values = c("#059D05","#9933FF"),labels=c("high", "low"))+
  facet_wrap(vars(Groups), scales = "free", labeller = as_labeller(GroupNames), drop=FALSE)+
  ylab("Effort") +
  coord_cartesian(ylim = c(2,10.4),xlim = c(0.4,1.05), expand=FALSE)+
  xlab("Accuracy")+
  theme_classic(base_size = 10)+
  theme(text=element_text(family="sans", size=10))+
  ggtitle(" ")+
  theme(legend.position = "bottom",legend.justification = "bottom", legend.title = element_blank())
effortplot

ggsave("effortplot.png", plot = effortplot, device = NULL, path = NULL,
       scale = 1, units = c("in"), width=15, height = 13, dpi = 600, limitsize = TRUE)
# in no-autism group:
naut <- get_anova_table(anova_test(subset(VPT_long_df, Autism == "naut"), dv = VPT_acc_wo, wid = ppn, between = c(Age_cat), within = type),correction="auto")
naut

## ANOVA Table (type III tests)
## 
##         Effect DFn DFd     F     p p<.05   ges
## 1      Age_cat   1  61 9.146 0.004     * 0.084
## 2         type   1  61 3.552 0.064       0.022
## 3 Age_cat:type   1  61 1.211 0.275       0.008

naut2 <- get_anova_table(anova_test(subset(VPT_long_df, Autism == "naut"), dv = VPT_acc_wo, wid = ppn, between = c(Age_cat), covariate = Effort_centered, within = type),correction="auto")
naut2

## ANOVA Table (type II tests)
## 
##                 Effect DFn DFd     F     p p<.05      ges
## 1      Effort_centered   1  51 0.004 0.951       4.35e-05
## 2              Age_cat   1  51 8.355 0.006     * 8.90e-02
## 3                 type   1  51 2.661 0.109       2.10e-02
## 4 Effort_centered:type   1  51 0.471 0.496       4.00e-03
## 5         Age_cat:type   1  51 1.194 0.280       9.00e-03