---
title: "R Notebook"
output: html_notebook
---

This is an [R Markdown](http://rmarkdown.rstudio.com) Notebook. When you execute code within the notebook, the results appear beneath the code. 

Try executing this chunk by clicking the *Run* button within the chunk or by placing your cursor inside it and pressing *Cmd+Shift+Enter*. 


```{r}
##### install and load necessary packages

# installs packages if necessary

if(!"emmeans" %in% rownames(installed.packages())) install.packages("emmeans")
if(!"tidyverse" %in% rownames(installed.packages())) install.packages("tidyverse")
if(!"splithalf" %in% rownames(installed.packages())) devtools::install_github("sdparsons/splithalf")
if(!"pwr" %in% rownames(installed.packages())) install.packages("pwr")
if(!"dplyr" %in% rownames(installed.packages())) install.packages("dplyr")
if(!"tidyr" %in% rownames(installed.packages())) install.packages("tidyr")
if(!"psych" %in% rownames(installed.packages())) install.packages("psych")
if(!"rstudioapi" %in% rownames(installed.packages())) install.packages("rstudioapi")


library("splithalf")   # internal consistency
library("psych")       # test-retest
library("tidyverse")
library("dplyr")
library("cvcqv")       # coefficient of variation
library("ggpubr")
library("pwr")
library("tidyr")
library("emmeans")
library("rstudioapi") 
```

```{r}
# clean slate
rm(list=ls())
cat("\014")
```

```{r}
sessionInfo()
```
```{r}
citation('splithalf')
```
```{r}
citation('psych')
```
```{r}
citation('pwr')
```
```{r}
citation('emmeans')
```


```{r}
##### Power Calculation ######
pwr.r.test(n = , r = .6, sig.level = .001, power = 0.9)
```

```{r}
########## Set working directory and load data 

pathwd<-sub("/AffFlexReliability_DataAnalysis_SingleTrials.Rmd","",system("find . -type f -name AffFlexReliability_DataAnalysis_SingleTrials.Rmd",intern=T)[1]) 
setwd(pathwd)
getwd()

df_affect_raw = read.csv("AffFlexReliability_SingleTrials.csv", dec = ".", sep = ";", header = T)
session1 <- df_affect_raw[df_affect_raw$time == "test",] 
session2 <- df_affect_raw[df_affect_raw$time == "retest",] 
```

```{r}
#### assure that relevant variables are treated as factors

session1$task <- factor(session1$task, levels = c("emotion","gender"))
session1$switch_only <- factor(session1$switch_only, levels = c("switch","baseline"))
session1$time <- factor(session1$time, levels = c("test", "retest"))
session1$task_switch <- factor(session1$task_switch, levels = c("switch_TO_emotion", "baseline_emotion", "switch_TO_gender", "baseline_gender"))
session1$is_error <- !as.logical(session1$is_hit)

session2$task <- factor(session2$task, levels = c("emotion","gender"))
session2$switch_only <- factor(session2$switch_only, levels = c("switch","baseline"))
session2$time <- factor(session2$time, levels = c("test", "retest"))
session2$task_switch <- factor(session2$task_switch, levels = c("switch_TO_emotion", "baseline_emotion", "switch_TO_gender", "baseline_gender"))
session2$is_error <- !as.logical(session2$is_hit)
```

```{r}
##### combine the sessions
df_affect_raw <- rbind(session1, session2)
```

```{r}
######## create subset of data for RT analyses
# filtering for only correct trials within the normative RT ranges and excluding errors
df_affect_raw_RT <- df_affect_raw[df_affect_raw$is_hit == "True" & df_affect_raw$is_in == "True",] 
session1_RT <- session1[session1$is_hit == "True" & session1$is_in == "True",] 
session2_RT <- session2[session2$is_hit == "True" & session2$is_in == "True",] 
```

```{r}
### aggregate data to create summary data averaged over trials (for Table 1)
meanRT <- aggregate(response_time~subj_idx+task+switch_only+time+task_switch,df_affect_raw_RT,mean)
meanER <- aggregate(is_error~subj_idx+task+switch_only+time+task_switch,df_affect_raw,mean)
```

```{r}
######## show means and standard deviations of response times and error rates averaged over participants (for Table 1)
descriptivesRT_mean <- aggregate(response_time~task+switch_only+time,meanRT,mean)
descriptivesRT_mean
```

```{r}
descriptivesRT_sd <- aggregate(response_time~task+switch_only+time,meanRT,sd)
descriptivesRT_sd
```

```{r}
descriptivesER_mean <- aggregate(is_error~task+switch_only+time,meanER,mean)
descriptivesER_mean
```
```{r}
descriptivesER_sd <- aggregate(is_error~task+switch_only+time,meanER,sd)
descriptivesER_sd
```

```{r}
##### aggregate data in wide format (averaged over trials but separately for both sessions) 
##### for calculation of switch costs (and later test-retest analysis)

meanRT_T1 <- aggregate(response_time~subj_idx+task_switch,session1_RT,mean)
meanRT_T1_wide <- reshape(meanRT_T1, idvar = "subj_idx", timevar = "task_switch", direction = "wide")

meanRT_T1_wide <- meanRT_T1_wide %>%
  rename(
    RT_Repeat_Emotion_test = response_time.baseline_emotion,
    RT_Switch_Emotion_test = response_time.switch_TO_emotion,
    RT_Repeat_Gender_test = response_time.baseline_gender,
    RT_Switch_Gender_test = response_time.switch_TO_gender
  )

meanRT_T2 <- aggregate(response_time~subj_idx+task_switch,session2_RT,mean)
meanRT_T2_wide <- reshape(meanRT_T2, idvar = "subj_idx", timevar = "task_switch", direction = "wide")

meanRT_T2_wide <- meanRT_T2_wide %>%
  rename(
    RT_Repeat_Emotion_retest = response_time.baseline_emotion,
    RT_Switch_Emotion_retest = response_time.switch_TO_emotion,
    RT_Repeat_Gender_retest = response_time.baseline_gender,
    RT_Switch_Gender_retest = response_time.switch_TO_gender
  )

meanRT_wide <- merge(meanRT_T1_wide,meanRT_T2_wide,by="subj_idx")
```

```{r}
meanRT_wide$RTcost_Emotion_test <- meanRT_wide$RT_Switch_Emotion_test - meanRT_wide$RT_Repeat_Emotion_test
meanRT_wide$RTcost_Gender_test <- meanRT_wide$RT_Switch_Gender_test - meanRT_wide$RT_Repeat_Gender_test
meanRT_wide$RTcost_Emotion_retest <- meanRT_wide$RT_Switch_Emotion_retest - meanRT_wide$RT_Repeat_Emotion_retest
meanRT_wide$RTcost_Gender_retest <- meanRT_wide$RT_Switch_Gender_retest - meanRT_wide$RT_Repeat_Gender_retest
```

```{r}
meanER_T1 <- aggregate(is_error~subj_idx+task_switch,session1,mean)
meanER_T1_wide <- reshape(meanER_T1, idvar = "subj_idx", timevar = "task_switch", direction = "wide")

meanER_T1_wide <- meanER_T1_wide %>%
  rename(
    ER_Repeat_Emotion_test = is_error.baseline_emotion,
    ER_Switch_Emotion_test = is_error.switch_TO_emotion,
    ER_Repeat_Gender_test = is_error.baseline_gender,
    ER_Switch_Gender_test = is_error.switch_TO_gender
  )

meanER_T2 <- aggregate(is_error~subj_idx+task_switch,session2,mean)
meanER_T2_wide <- reshape(meanER_T2, idvar = "subj_idx", timevar = "task_switch", direction = "wide")

meanER_T2_wide <- meanER_T2_wide %>%
  rename(
    ER_Repeat_Emotion_retest = is_error.baseline_emotion,
    ER_Switch_Emotion_retest = is_error.switch_TO_emotion,
    ER_Repeat_Gender_retest = is_error.baseline_gender,
    ER_Switch_Gender_retest = is_error.switch_TO_gender
  )

meanER_wide <- merge(meanER_T1_wide,meanER_T2_wide,by="subj_idx")
```

```{r}
meanER_wide$ERcost_Emotion_test <- meanER_wide$ER_Switch_Emotion_test - meanER_wide$ER_Repeat_Emotion_test
meanER_wide$ERcost_Gender_test <- meanER_wide$ER_Switch_Gender_test - meanER_wide$ER_Repeat_Gender_test
meanER_wide$ERcost_Emotion_retest <- meanER_wide$ER_Switch_Emotion_retest - meanER_wide$ER_Repeat_Emotion_retest
meanER_wide$ERcost_Gender_retest <- meanER_wide$ER_Switch_Gender_retest - meanER_wide$ER_Repeat_Gender_retest
```

```{r}
####### create merged dataframe in wide format (i.e., for characterization of switch costs)
Switchcosts <- merge(meanRT_wide,meanER_wide,by="subj_idx")
```

```{r}
##### subset switch cost variables to show means and standard deviations of RT and ER switch costs averaged over participants (for Table 1)
Switchcosts_descriptives <-  select(Switchcosts ,RTcost_Gender_test
                                                ,RTcost_Gender_retest
                                                ,RTcost_Emotion_test
                                                ,RTcost_Emotion_retest
                                                ,ERcost_Gender_test
                                                ,ERcost_Gender_retest
                                                ,ERcost_Emotion_test
                                                ,ERcost_Emotion_retest)
```


```{r}
sapply(Switchcosts_descriptives, mean, na.rm=TRUE)
```

```{r}
sapply(Switchcosts_descriptives, sd, na.rm=TRUE)
```

```{r}
########### Calculate Effect sizes for switch cost effects (switch vs. repeat) for table 1
effsize::cohen.d(meanRT_wide$RT_Switch_Gender_test, meanRT_wide$RT_Repeat_Gender_test,hedges.correction=TRUE)
effsize::cohen.d(meanRT_wide$RT_Switch_Emotion_test, meanRT_wide$RT_Repeat_Emotion_test,hedges.correction=TRUE)

effsize::cohen.d(meanRT_wide$RT_Switch_Gender_retest, meanRT_wide$RT_Repeat_Gender_retest,hedges.correction=TRUE)
effsize::cohen.d(meanRT_wide$RT_Switch_Emotion_retest, meanRT_wide$RT_Repeat_Emotion_retest,hedges.correction=TRUE)
```

```{r}
effsize::cohen.d(meanER_wide$ER_Switch_Gender_test, meanER_wide$ER_Repeat_Gender_test,hedges.correction=TRUE)
effsize::cohen.d(meanER_wide$ER_Switch_Emotion_test, meanER_wide$ER_Repeat_Emotion_test,hedges.correction=TRUE)

effsize::cohen.d(meanER_wide$ER_Switch_Gender_retest, meanER_wide$ER_Repeat_Gender_retest,hedges.correction=TRUE)
effsize::cohen.d(meanER_wide$ER_Switch_Emotion_retest, meanER_wide$ER_Repeat_Emotion_retest,hedges.correction=TRUE)
```

###### Calculate repeated measures ANOVAs and corresponding post-hoc contrasts for characterization of switch cost effects ########

```{r}
#### calculate repeated measures ANOVA for response times
RT.aov <- aov( response_time~(time*task*switch_only)+Error(subj_idx/(time*task*switch_only)), meanRT)
summary(RT.aov)
```

```{r}
#### post-hoc testing (pairwise comparisons)
RT.emm <- emmeans(RT.aov, ~ time * task * switch_only)
pairs(RT.emm, simple = "each", combine = TRUE, adjust = "mvt")
```

```{r}
###### characterize main effects

##### EMMs repeat vs. switch (main effect condition)
emmeans(RT.aov, pairwise ~ switch_only)
```

```{r}
##### EMMs test vs. retest (main effect session)
emmeans(RT.aov, pairwise ~ time)
```

```{r}
##### EMMs emotion vs. gender task (main effect task)
emmeans(RT.aov, pairwise ~ task)
```

```{r}
#### visualize interaction effects
emmip(RT.aov, switch_only ~ task | time )
```


```{r}
#### EMMs for Condition x Session interaction
emmeans(RT.aov, pairwise ~ time | switch_only)
```

```{r}
#### Task x Condition interaction
emmeans(RT.aov, pairwise ~ task | switch_only)
```

```{r}
#### calculate repeated measures ANOVA for error rates
ER.aov <- aov( is_error~(time*task*switch_only)+Error(subj_idx/(time*task*switch_only)), meanER)
summary(ER.aov)
```
```{r}
#### post-hoc testing (pairwise comparisons)
ER.emm <- emmeans(ER.aov, ~ time * task * switch_only)
pairs(ER.emm, simple = "each", combine = TRUE, adjust = "mvt")
```
```{r}
#### visualize interaction effects
emmip(ER.aov, switch_only ~ task | time )
```

```{r}
###### characterize main effects

##### EMMs repeat vs. switch (main effect condition)
emmeans(ER.aov, pairwise ~ switch_only)
```
```{r}
##### EMMs test vs. retest (main effect session)
emmeans(ER.aov, pairwise ~ time)
```
```{r}
#### EMMs for Condition x Session interaction
emmeans(ER.aov, pairwise ~ time | switch_only)
```
```{r}
#### EMMs for Condition x Task interaction
emmeans(ER.aov, pairwise ~ task | switch_only)
```



```{r}
##### Correlations between RT and ER Switch costs (to judge the amount of shared variance between these two allegedly 
##### interchangeable behavioral indicators of task-switching ability)

#### testing normality of variables (decision, which correlation coefficient)

shapiro.test(Switchcosts$RTcost_Emotion_test) #sig
shapiro.test(Switchcosts$ERcost_Emotion_test) #sig

shapiro.test(Switchcosts$RTcost_Gender_test) #nonsig
shapiro.test(Switchcosts$ERcost_Gender_test) #sig

shapiro.test(Switchcosts$RTcost_Emotion_retest) #nonsig
shapiro.test(Switchcosts$ERcost_Emotion_retest) #sig

shapiro.test(Switchcosts$RTcost_Gender_retest) #sig
shapiro.test(Switchcosts$ERcost_Gender_retest) #sig

```

```{r}
#### Correlate RT and ER switch costs during the Test session

cor.test(Switchcosts$RTcost_Emotion_test,Switchcosts$ERcost_Emotion_test, method="spearman")
cor.test(Switchcosts$RTcost_Gender_test,Switchcosts$ERcost_Gender_test, method="spearman")

```
```{r}
#### Correlate RT and ER switch costs during the Retest session

cor.test(Switchcosts$RTcost_Emotion_retest,Switchcosts$ERcost_Emotion_retest, method="spearman")
cor.test(Switchcosts$RTcost_Gender_retest,Switchcosts$ERcost_Gender_retest, method="spearman")

```

#####################################################
############ Internal consistencies #################

########### RT #############

```{r permutation internal consistency, results = 'hide', eval=FALSE}
#### calculate internal consistency for RT switchcosts during the test session (for each task)

require(splithalf)
splithalf(data = session1_RT,
          outcome = "RT",
          score = "difference",
          permutations = 5000,
          var.trialnum = "X",
          var.condition = "task",
          conditionlist = c("emotion", "gender"),
          var.compare = "switch_only",
          compare1 = "baseline",
          compare2 = "switch",
          var.participant = "subj_idx",
          var.RT = "response_time" )

```

```{r permutation internal consistency, results = 'hide', eval=FALSE}
#### calculate internal consistency for RT switchcosts during the retest session (for each task)

require(splithalf)
splithalf(data = session2_RT,
          outcome = "RT",
          score = "difference",
          permutations = 5000,
          var.trialnum = "X",
          var.condition = "task",
          conditionlist = c("emotion", "gender"),
          var.compare = "switch_only",
          compare1 = "baseline",
          compare2 = "switch",
          var.participant = "subj_idx",
          var.RT = "response_time" )

```

###################################
############# ER ##################

```{r permutation internal consistency, results = 'hide', eval=FALSE}
#### calculate internal consistency for ER switchcosts during the test session (for each task)

require(splithalf)
splithalf(data = session1,
          outcome = "accuracy",
          score = "difference",
          permutations = 5000,
          var.trialnum = "X",
          var.condition = "task",
          conditionlist = c("emotion", "gender"),
          var.compare = "switch_only",
          compare1 = "baseline",
          compare2 = "switch",
          var.participant = "subj_idx",
          var.ACC = "is_error" )

```

```{r permutation internal consistency, results = 'hide', eval=FALSE}
#### calculate internal consistency for ER switchcosts during the retest session (for each task)

require(splithalf)
splithalf(data = session2,
          outcome = "accuracy",
          score = "difference",
          permutations = 5000,
          var.trialnum = "X",
          var.condition = "task",
          conditionlist = c("emotion", "gender"),
          var.compare = "switch_only",
          compare1 = "baseline",
          compare2 = "switch",
          var.participant = "subj_idx",
          var.ACC = "is_error" )

```

###########################################
########## Test-Retest Analysis ###########

```{r}
### calculate test-retest reliability for RT during repeat trials in the gender task
RT_Repeat_Gender <- cbind(Switchcosts$RT_Repeat_Gender_test,Switchcosts$RT_Repeat_Gender_retest)
ICC.RT_Repeat_Gender <- ICC(RT_Repeat_Gender)
ICC.RT_Repeat_Gender
cor.test(Switchcosts$RT_Repeat_Gender_test,Switchcosts$RT_Repeat_Gender_retest) # corresponding Pearson correlation
```

```{r}
### calculate test-retest reliability for RT during switch trials in the gender task
RT_Switch_Gender <- cbind(Switchcosts$RT_Switch_Gender_test,Switchcosts$RT_Switch_Gender_retest)
ICC.RT_Switch_Gender <- ICC(RT_Switch_Gender)
ICC.RT_Switch_Gender
cor.test(Switchcosts$RT_Switch_Gender_test,Switchcosts$RT_Switch_Gender_retest) # corresponding Pearson correlation
```

```{r}
### calculate test-retest reliability for RT switch costs in the gender task
RTcost_Gender <- cbind(Switchcosts$RTcost_Gender_test,Switchcosts$RTcost_Gender_retest)
ICC.RT_RTcost_Gender <- ICC(RTcost_Gender)
ICC.RT_RTcost_Gender
cor.test(Switchcosts$RTcost_Gender_test,Switchcosts$RTcost_Gender_retest) # corresponding Pearson correlation
```

```{r}
### calculate test-retest reliability for RT during repeat trials in the emotion task
RT_Repeat_Emotion <- cbind(Switchcosts$RT_Repeat_Emotion_test,Switchcosts$RT_Repeat_Emotion_retest)
ICC.RT_Repeat_Emotion <- ICC(RT_Repeat_Emotion)
ICC.RT_Repeat_Emotion
cor.test(Switchcosts$RT_Repeat_Emotion_test,Switchcosts$RT_Repeat_Emotion_retest) # corresponding Pearson correlation
```

```{r}
### calculate test-retest reliability for RT during switch trials in the emotion task
RT_Switch_Emotion <- cbind(Switchcosts$RT_Switch_Emotion_test,Switchcosts$RT_Switch_Emotion_retest)
ICC.RT_Switch_Emotion <- ICC(RT_Switch_Emotion)
ICC.RT_Switch_Emotion
cor.test(Switchcosts$RT_Switch_Emotion_test,Switchcosts$RT_Switch_Emotion_retest) # corresponding Pearson correlation
```

```{r}
### calculate test-retest reliability for RT switch costs in the emotion task
RTcost_Emotion <- cbind(Switchcosts$RTcost_Emotion_test,Switchcosts$RTcost_Emotion_retest)
ICC.RT_RTcost_Emotion <- ICC(RTcost_Emotion)
ICC.RT_RTcost_Emotion
cor.test(Switchcosts$RTcost_Emotion_test,Switchcosts$RTcost_Emotion_retest) # corresponding Pearson correlation
```

```{r}
### calculate test-retest reliability for ER during repeat trials in the gender task
ER_Repeat_Gender <- cbind(Switchcosts$ER_Repeat_Gender_test,Switchcosts$ER_Repeat_Gender_retest)
ICC.ER_Repeat_Gender <- ICC(ER_Repeat_Gender)
ICC.ER_Repeat_Gender
cor.test(Switchcosts$ER_Repeat_Gender_test,Switchcosts$ER_Repeat_Gender_retest) # corresponding Pearson correlation
```

```{r}
### calculate test-retest reliability for ER during switch trials in the gender task
ER_Switch_Gender <- cbind(Switchcosts$ER_Switch_Gender_test,Switchcosts$ER_Switch_Gender_retest)
ICC.ER_Switch_Gender <- ICC(ER_Switch_Gender)
ICC.ER_Switch_Gender
cor.test(Switchcosts$ER_Switch_Gender_test,Switchcosts$ER_Switch_Gender_retest) # corresponding Pearson correlation
```

```{r}
### calculate test-retest reliability for ER during switch costs in the gender task
ERcost_Gender <- cbind(Switchcosts$ERcost_Gender_test,Switchcosts$ERcost_Gender_retest)
ICC.ER_ERcost_Gender <- ICC(ERcost_Gender)
ICC.ER_ERcost_Gender
cor.test(Switchcosts$ERcost_Gender_test,Switchcosts$ERcost_Gender_retest) # corresponding Pearson correlation
```

```{r}
### calculate test-retest reliability for ER during repeat trials in the emotion task
ER_Repeat_Emotion <- cbind(Switchcosts$ER_Repeat_Emotion_test,Switchcosts$ER_Repeat_Emotion_retest)
ICC.ER_Repeat_Emotion <- ICC(ER_Repeat_Emotion)
ICC.ER_Repeat_Emotion
cor.test(Switchcosts$ER_Repeat_Emotion_test,Switchcosts$ER_Repeat_Emotion_retest) # corresponding Pearson correlation
```

```{r}
### calculate test-retest reliability for ER during switch trials in the emotion task
ER_Switch_Emotion <- cbind(Switchcosts$ER_Switch_Emotion_test,Switchcosts$ER_Switch_Emotion_retest)
ICC.ER_Switch_Emotion <- ICC(ER_Switch_Emotion)
ICC.ER_Switch_Emotion
cor.test(Switchcosts$ER_Switch_Emotion_test,Switchcosts$ER_Switch_Emotion_retest) # corresponding Pearson correlation
```

```{r}
### calculate test-retest reliability for ER switch costs in the emotion task
ERcost_Emotion <- cbind(Switchcosts$ERcost_Emotion_test,Switchcosts$ERcost_Emotion_retest)
ICC.ER_ERcost_Emotion <- ICC(ERcost_Emotion)
ICC.ER_ERcost_Emotion
cor.test(Switchcosts$ERcost_Emotion_test,Switchcosts$ERcost_Emotion_retest) # corresponding Pearson correlation
```

###############################################
############# Figure 2 ########################

```{r}
###### save variance components of between-sessions, between-subjects and error variance for 
###### stacked barplot (figure 2)

Var_RT_Repeat_Gender <- ICC.RT_Repeat_Gender$lme[1:3,2]
Var_RT_Switch_Gender <- ICC.RT_Switch_Gender$lme[1:3,2]
Var_RT_RTcost_Gender <- ICC.RT_RTcost_Gender$lme[1:3,2]

Var_RT_Repeat_Emotion <- ICC.RT_Repeat_Emotion$lme[1:3,2]
Var_RT_Switch_Emotion <- ICC.RT_Switch_Emotion$lme[1:3,2]
Var_RT_RTcost_Emotion <- ICC.RT_RTcost_Emotion$lme[1:3,2]

Var_ER_Repeat_Gender <- ICC.ER_Repeat_Gender$lme[1:3,2]
Var_ER_Switch_Gender <- ICC.ER_Switch_Gender$lme[1:3,2]
Var_ER_ERcost_Gender <- ICC.ER_ERcost_Gender$lme[1:3,2]

Var_ER_Repeat_Emotion <- ICC.ER_Repeat_Emotion$lme[1:3,2]
Var_ER_Switch_Emotion <- ICC.ER_Switch_Emotion$lme[1:3,2]
Var_ER_ERcost_Emotion <- ICC.ER_ERcost_Emotion$lme[1:3,2]

Type <- rep(c('Between-subjects','Between-sessions','Error variance'),3)
Condition <- c(rep(c('Repeat'),3),rep(c('Switch'),3),rep(c('Switch costs'),3))
```

```{r}
# RT
Percent <- c(Var_RT_Repeat_Gender, Var_RT_Switch_Gender, Var_RT_RTcost_Gender)
Task <- rep(c('Gender'),9)
Measure <- rep(c('RT'),9)
RT_Gender <- data.frame(Percent,Type,Condition,Task,Measure)

Percent <- c(Var_RT_Repeat_Emotion, Var_RT_Switch_Emotion, Var_RT_RTcost_Emotion)
Task <- rep(c('Emotion'),9)
Measure <- rep(c('RT'),9)
RT_Emotion <- data.frame(Percent,Type,Condition,Task,Measure)

RT <- rbind(RT_Gender,RT_Emotion)

# ER
Percent <- c(Var_ER_Repeat_Gender, Var_ER_Switch_Gender, Var_ER_ERcost_Gender)
Task <- rep(c('Gender'),9)
Measure <- rep(c('ER'),9)
ER_Gender <- data.frame(Percent,Type,Condition,Task,Measure)

Percent <- c(Var_ER_Repeat_Emotion, Var_ER_Switch_Emotion, Var_ER_ERcost_Emotion)
Task <- rep(c('Emotion'),9)
Measure <- rep(c('ER'),9)
ER_Emotion <- data.frame(Percent,Type,Condition,Task,Measure)

ER <- rbind(ER_Gender,ER_Emotion)
```

```{r}
#stacked Response Time
p_RT <- ggplot(RT, aes(fill=Type, y=Percent, x=Condition)) +  
  geom_bar(position="stack", stat="identity") +
  ylab("Percentage of variance") +
  facet_wrap(~Task) +
  scale_fill_grey() +
  theme_minimal()

p_RT <- p_RT + guides(fill=guide_legend(title=NULL)) +
            ggtitle("RT based measures") +
            theme(
              legend.text = element_text(size = 16),
              plot.title = element_text(size=16,face="bold",hjust = 0.5),
              text = element_text(size=16),
              axis.title.y = element_text(size=16),
              axis.title.x = element_blank(),
              axis.text.x = element_text(angle=45, hjust=1),
              strip.text = element_text(size=16)
              )

p_RT

```

```{r}
#stacked Error Rate
p_ER <- ggplot(ER, aes(fill=Type, y=Percent, x=Condition)) +  
geom_bar(position="stack", stat="identity") +
  facet_wrap(~Task) +
  scale_fill_grey() +
  theme_minimal()

p_ER <- p_ER + guides(fill=guide_legend(title=NULL)) +
            ggtitle("ER based measures") +
            theme(
              legend.text = element_text(size = 16),
              plot.title = element_text(size=16,face="bold",hjust = 0.5),
              text = element_text(size=16),
              axis.title.y =  element_blank(),
              axis.title.x = element_blank(),
              axis.text.x = element_text(angle=45, hjust=1),
              strip.text = element_text(size=16)
              )

p_ER

```

```{r}
figure <- ggarrange(p_RT, p_ER,
                    labels = c("A)", "B)"),
                    ncol = 2,
                    common.legend = TRUE, legend = "bottom")

plot(figure)
```

###########################################################################################################
###########################################################################################################
############################ Repeat reliability analysis on first run only #################################

```{r}
######## create subsets of data with only first run
session1_run1 <- session1[session1$X < 240.0,]
session2_run1 <- session2[session2$X < 240.0,] 

session1_RT_run1 <- session1_RT[session1_RT$X < 240.0,]
session2_RT_run1 <- session2_RT[session2_RT$X < 240.0,]
```

```{r permutation internal consistency, results = 'hide', eval=FALSE}
#### calculate internal consistency for RT switchcosts during the test session (for each task, first run only)

require(splithalf)
splithalf(data = session1_RT_run1,
          outcome = "RT",
          score = "difference",
          permutations = 5000,
          var.trialnum = "X",
          var.condition = "task",
          conditionlist = c("emotion", "gender"),
          var.compare = "switch_only",
          compare1 = "baseline",
          compare2 = "switch",
          var.participant = "subj_idx",
          var.RT = "response_time" )

```

```{r permutation internal consistency, results = 'hide', eval=FALSE}
#### calculate internal consistency for RT switchcosts during the retest session (for each task, first run only)

require(splithalf)
splithalf(data = session2_RT_run1,
          outcome = "RT",
          score = "difference",
          permutations = 5000,
          var.trialnum = "X",
          var.condition = "task",
          conditionlist = c("emotion", "gender"),
          var.compare = "switch_only",
          compare1 = "baseline",
          compare2 = "switch",
          var.participant = "subj_idx",
          var.RT = "response_time" )

```

```{r permutation internal consistency, results = 'hide', eval=FALSE}
#### calculate internal consistency for ER switchcosts during the test session (for each task, first run only)

require(splithalf)
splithalf(data = session1_run1,
          outcome = "accuracy",
          score = "difference",
          permutations = 5000,
          var.trialnum = "X",
          var.condition = "task",
          conditionlist = c("emotion", "gender"),
          var.compare = "switch_only",
          compare1 = "baseline",
          compare2 = "switch",
          var.participant = "subj_idx",
          var.ACC = "is_error" )

```

```{r permutation internal consistency, results = 'hide', eval=FALSE}
#### calculate internal consistency for ER switchcosts during the retest session (for each task, first run only)

require(splithalf)
splithalf(data = session2_run1,
          outcome = "accuracy",
          score = "difference",
          permutations = 5000,
          var.trialnum = "X",
          var.condition = "task",
          conditionlist = c("emotion", "gender"),
          var.compare = "switch_only",
          compare1 = "baseline",
          compare2 = "switch",
          var.participant = "subj_idx",
          var.ACC = "is_error" )

```

```{r}
##### aggregate data in wide format for calculation of switch costs (and later test-retest analysis, see above)

meanRT_T1_run1 <- aggregate(response_time~subj_idx+task_switch,session1_RT_run1,mean)
meanRT_T1_run1_wide <- reshape(meanRT_T1_run1, idvar = "subj_idx", timevar = "task_switch", direction = "wide")

meanRT_T1_run1_wide <- meanRT_T1_run1_wide %>%
  rename(
    RT_Repeat_Emotion_test = response_time.baseline_emotion,
    RT_Switch_Emotion_test = response_time.switch_TO_emotion,
    RT_Repeat_Gender_test = response_time.baseline_gender,
    RT_Switch_Gender_test = response_time.switch_TO_gender
  )

meanRT_T2_run1 <- aggregate(response_time~subj_idx+task_switch,session2_RT_run1,mean)
meanRT_T2_run1_wide <- reshape(meanRT_T2_run1, idvar = "subj_idx", timevar = "task_switch", direction = "wide")

meanRT_T2_run1_wide <- meanRT_T2_run1_wide %>%
  rename(
    RT_Repeat_Emotion_retest = response_time.baseline_emotion,
    RT_Switch_Emotion_retest = response_time.switch_TO_emotion,
    RT_Repeat_Gender_retest = response_time.baseline_gender,
    RT_Switch_Gender_retest = response_time.switch_TO_gender
  )

meanRT_run1_wide <- merge(meanRT_T1_run1_wide,meanRT_T2_run1_wide,by="subj_idx")
```

```{r}
meanRT_run1_wide$RTcost_Emotion_test <- meanRT_run1_wide$RT_Switch_Emotion_test - meanRT_run1_wide$RT_Repeat_Emotion_test
meanRT_run1_wide$RTcost_Gender_test <- meanRT_run1_wide$RT_Switch_Gender_test - meanRT_run1_wide$RT_Repeat_Gender_test
meanRT_run1_wide$RTcost_Emotion_retest <- meanRT_run1_wide$RT_Switch_Emotion_retest - meanRT_run1_wide$RT_Repeat_Emotion_retest
meanRT_run1_wide$RTcost_Gender_retest <- meanRT_run1_wide$RT_Switch_Gender_retest - meanRT_run1_wide$RT_Repeat_Gender_retest
```


```{r}
meanER_T1_run1 <- aggregate(is_error~subj_idx+task_switch,session1_run1,mean)
meanER_T1_run1_wide <- reshape(meanER_T1_run1, idvar = "subj_idx", timevar = "task_switch", direction = "wide")

meanER_T1_run1_wide <- meanER_T1_run1_wide %>%
  rename(
    ER_Repeat_Emotion_test = is_error.baseline_emotion,
    ER_Switch_Emotion_test = is_error.switch_TO_emotion,
    ER_Repeat_Gender_test = is_error.baseline_gender,
    ER_Switch_Gender_test = is_error.switch_TO_gender
  )

meanER_T2_run1 <- aggregate(is_error~subj_idx+task_switch,session2_run1,mean)
meanER_T2_run1_wide <- reshape(meanER_T2_run1, idvar = "subj_idx", timevar = "task_switch", direction = "wide")

meanER_T2_run1_wide <- meanER_T2_run1_wide %>%
  rename(
    ER_Repeat_Emotion_retest = is_error.baseline_emotion,
    ER_Switch_Emotion_retest = is_error.switch_TO_emotion,
    ER_Repeat_Gender_retest = is_error.baseline_gender,
    ER_Switch_Gender_retest = is_error.switch_TO_gender
  )

meanER_run1_wide <- merge(meanER_T1_run1_wide,meanER_T2_run1_wide,by="subj_idx")
```

```{r}
meanER_run1_wide$ERcost_Emotion_test <- meanER_run1_wide$ER_Switch_Emotion_test - meanER_run1_wide$ER_Repeat_Emotion_test
meanER_run1_wide$ERcost_Gender_test <- meanER_run1_wide$ER_Switch_Gender_test - meanER_run1_wide$ER_Repeat_Gender_test
meanER_run1_wide$ERcost_Emotion_retest <- meanER_run1_wide$ER_Switch_Emotion_retest - meanER_run1_wide$ER_Repeat_Emotion_retest
meanER_run1_wide$ERcost_Gender_retest <- meanER_run1_wide$ER_Switch_Gender_retest - meanER_run1_wide$ER_Repeat_Gender_retest
```

```{r}
Switchcosts_run1 <- merge(meanRT_run1_wide,meanER_run1_wide,by="subj_idx")
```

################## Test-retest Analysis ##################

```{r}
### calculate test-retest reliability for RT during repeat trials in the gender task (first run only)

RT_run1_Repeat_Gender <- cbind(Switchcosts_run1$RT_Repeat_Gender_test,Switchcosts_run1$RT_Repeat_Gender_retest)
ICC.RT_run1_Repeat_Gender <- ICC(RT_run1_Repeat_Gender)
ICC.RT_run1_Repeat_Gender
cor.test(Switchcosts_run1$RT_Repeat_Gender_test,Switchcosts_run1$RT_Repeat_Gender_retest) # corresponding Pearson correlation
```

```{r}
### calculate test-retest reliability for RT during switch trials in the gender task (first run only)

RT_run1_Switch_Gender <- cbind(Switchcosts_run1$RT_Switch_Gender_test,Switchcosts_run1$RT_Switch_Gender_retest)
ICC.RT_run1_Switch_Gender <- ICC(RT_run1_Switch_Gender)
ICC.RT_run1_Switch_Gender
cor.test(Switchcosts_run1$RT_Switch_Gender_test,Switchcosts_run1$RT_Switch_Gender_retest) # corresponding Pearson correlation
```

```{r}
### calculate test-retest reliability for RT switch costs in the gender task (first run only)

RTcost_run1_Gender <- cbind(Switchcosts_run1$RTcost_Gender_test,Switchcosts_run1$RTcost_Gender_retest)
ICC.RT_run1_RTcost_Gender <- ICC(RTcost_run1_Gender)
ICC.RT_run1_RTcost_Gender
cor.test(Switchcosts_run1$RTcost_Gender_test,Switchcosts_run1$RTcost_Gender_retest) # corresponding Pearson correlation
```

```{r}
### calculate test-retest reliability for RT during repeat trials in the emotion task (first run only)

RT_run1_Repeat_Emotion <- cbind(Switchcosts_run1$RT_Repeat_Emotion_test,Switchcosts_run1$RT_Repeat_Emotion_retest)
ICC.RT_run1_Repeat_Emotion <- ICC(RT_run1_Repeat_Emotion)
ICC.RT_run1_Repeat_Emotion
cor.test(Switchcosts_run1$RT_Repeat_Emotion_test,Switchcosts_run1$RT_Repeat_Emotion_retest) # corresponding Pearson correlation
```

```{r}
### calculate test-retest reliability for RT during switch trials in the emotion task (first run only)

RT_run1_Switch_Emotion <- cbind(Switchcosts_run1$RT_Switch_Emotion_test,Switchcosts_run1$RT_Switch_Emotion_retest)
ICC.RT_run1_Switch_Emotion <- ICC(RT_run1_Switch_Emotion)
ICC.RT_run1_Switch_Emotion
cor.test(Switchcosts_run1$RT_Switch_Emotion_test,Switchcosts_run1$RT_Switch_Emotion_retest) # corresponding Pearson correlation
```

```{r}
### calculate test-retest reliability for RT switch costs in the emotion task (first run only)

RTcost_run1_Emotion <- cbind(Switchcosts_run1$RTcost_Emotion_test,Switchcosts_run1$RTcost_Emotion_retest)
ICC.RT_run1_RTcost_Emotion <- ICC(RTcost_run1_Emotion)
ICC.RT_run1_RTcost_Emotion
cor.test(Switchcosts_run1$RTcost_Emotion_test,Switchcosts_run1$RTcost_Emotion_retest) # corresponding Pearson correlation
```

```{r}
### calculate test-retest reliability for ER during repeat trials in the gender task (first run only)

ER_run1_Repeat_Gender <- cbind(Switchcosts_run1$ER_Repeat_Gender_test,Switchcosts_run1$ER_Repeat_Gender_retest)
ICC.ER_run1_Repeat_Gender <- ICC(ER_run1_Repeat_Gender)
ICC.ER_run1_Repeat_Gender
cor.test(Switchcosts_run1$ER_Repeat_Gender_test,Switchcosts_run1$ER_Repeat_Gender_retest) # corresponding Pearson correlation
```

```{r}
### calculate test-retest reliability for ER during switch trials in the gender task (first run only)

ER_run1_Switch_Gender <- cbind(Switchcosts_run1$ER_Switch_Gender_test,Switchcosts_run1$ER_Switch_Gender_retest)
ICC.ER_run1_Switch_Gender <- ICC(ER_run1_Switch_Gender)
ICC.ER_run1_Switch_Gender
cor.test(Switchcosts_run1$ER_Switch_Gender_test,Switchcosts_run1$ER_Switch_Gender_retest) # corresponding Pearson correlation
```

```{r}
### calculate test-retest reliability for ER switch costs in the gender task (first run only)

ERcost_run1_Gender <- cbind(Switchcosts_run1$ERcost_Gender_test,Switchcosts_run1$ERcost_Gender_retest)
ICC.ER_run1_ERcost_Gender <- ICC(ERcost_run1_Gender)
ICC.ER_run1_ERcost_Gender
cor.test(Switchcosts_run1$ERcost_Gender_test,Switchcosts_run1$ERcost_Gender_retest) # corresponding Pearson correlation
```

```{r}
### calculate test-retest reliability for ER during repeat trials in the emotion task (first run only)

ER_run1_Repeat_Emotion <- cbind(Switchcosts_run1$ER_Repeat_Emotion_test,Switchcosts_run1$ER_Repeat_Emotion_retest)
ICC.ER_run1_Repeat_Emotion <- ICC(ER_run1_Repeat_Emotion)
ICC.ER_run1_Repeat_Emotion
cor.test(Switchcosts_run1$ER_Repeat_Emotion_test,Switchcosts_run1$ER_Repeat_Emotion_retest) # corresponding Pearson correlation
```

```{r}
### calculate test-retest reliability for ER during switch trials in the emotion task (first run only)

ER_run1_Switch_Emotion <- cbind(Switchcosts_run1$ER_Switch_Emotion_test,Switchcosts_run1$ER_Switch_Emotion_retest)
ICC.ER_run1_Switch_Emotion <- ICC(ER_run1_Switch_Emotion)
ICC.ER_run1_Switch_Emotion
cor.test(Switchcosts_run1$ER_Switch_Emotion_test,Switchcosts_run1$ER_Switch_Emotion_retest) # corresponding Pearson correlation
```

```{r}
### calculate test-retest reliability for ER switch costs in the emotion task (first run only)

ERcost_run1_Emotion <- cbind(Switchcosts_run1$ERcost_Emotion_test,Switchcosts_run1$ERcost_Emotion_retest)
ICC.ER_run1_ERcost_Emotion <- ICC(ERcost_run1_Emotion)
ICC.ER_run1_ERcost_Emotion
cor.test(Switchcosts_run1$ERcost_Emotion_test,Switchcosts_run1$ERcost_Emotion_retest) # corresponding Pearson correlation
```

Add a new chunk by clicking the *Insert Chunk* button on the toolbar or by pressing *Cmd+Option+I*.

When you save the notebook, an HTML file containing the code and output will be saved alongside it (click the *Preview* button or press *Cmd+Shift+K* to preview the HTML file). 

The preview shows you a rendered HTML copy of the contents of the editor. Consequently, unlike *Knit*, *Preview* does not run any R code chunks. Instead, the output of the chunk when it was last run in the editor is displayed.

