Upload files to 'derivates'

derivates/Behavioral_analysis.Rmd

@@ -0,0 +1,235 @@
+---
+title: "Electrophysiological signatures of temporal context in the bisection task"
+author: "Cemre Baykan, Xiuna Zhu, Artyom Zinchenko, Hermann Mueller, Zhuanghua Shi"
+date: "February 2023"
+---
+
+```{r setup, include=FALSE}
+knitr::opts_chunk$set(echo = TRUE, eval = TRUE)
+source('dataana.R')
+options(mc.cores = parallel::detectCores())
+# flag for save figures 
+saveFigure = TRUE
+```
+
+# Experiment 1- Spacing Effect
+## Fit psychometric functions
+```{r}
+fits_EEG_exp1 = quickpsy(dat_EEG_exp1, x = targetDur, k =RP, prob = .5,
+                         grouping = .(cond, SubName),
+                         fun = logistic_fun,
+                         thresholds = FALSE,bootstrap='none')
+```
+
+
+```{r, fig.cap= 'PSE: Bisection task performance in Experiment 1'}
+# plot fitted function
+plot(fits_EEG_exp1)+ mytheme
+```
+
+```{r}
+par_EEG_exp1 <-fits_EEG_exp1$par %>%  dplyr::select(cond, SubName, parn, par) %>%
+  pivot_wider(names_from = parn, values_from = par) %>% 
+  mutate(jnd = log(3)/p2)
+
+par_EEG_exp1$pse <- par_EEG_exp1$p1
+par_EEG_exp1$p2 <-  NULL
+par_EEG_exp1$p1 <-  NULL
+head(par_EEG_exp1)
+```
+
+
+```{r}
+dat_EEG_exp1$cond= factor(dat_EEG_exp1$cond, levels=c("PS","NS"))
+fits_EEG_exp1_all = quickpsy(dat_EEG_exp1, x = targetDur, k =RP, 
+                             grouping = .(cond),
+                             fun = logistic_fun,
+                             thresholds = FALSE,bootstrap='none')
+fits_EEG_exp1_all
+```
+
+## plot averaged 
+```{r, fig.cap= 'PSE: Bisection task performance in Experiment 1'}
+# plot fitted function
+plot_fit_EEG_exp1 <- plot(fits_EEG_exp1_all) + mytheme+
+  labs(x = "Target interval (ms)", y = "Proportion of 'long' responses") +
+  scale_color_manual(values = colors_plot) +
+  scale_fill_manual(values = colors_plot) +
+  theme(legend.position = "top", legend.title = element_blank()) 
+
+plot_fit_EEG_exp1
+```
+
+```{r}
+plot_pse_EEG_exp1_box <- ggplot(par_EEG_exp1%>% mutate(cond = factor(cond, levels = c("PS", "NS"))),aes(x=cond, y=pse, color = cond))+
+  geom_boxplot(outlier.shape = NA) + 
+  geom_jitter(position=position_jitter(0.2))+
+  geom_signif(comparisons=list(c("PS", "NS")), annotations="*", y_position = 1350, tip_length = 0, vjust=0.4, color = 'black') +
+  labs(x = "Context", y = "PSE (ms)") +
+  coord_cartesian(ylim = c(450,1500)) +
+  scale_color_manual(values = colors_plot) +
+  scale_fill_manual(values = colors_plot) +
+   mytheme + theme(legend.position = "none")
+
+plot_pse_EEG_exp1_box
+```
+
+```{r}
+plot_jnd_EEG_exp1_box <- ggplot(par_EEG_exp1%>% mutate(cond = factor(cond, levels = c("PS", "NS"))),aes(x=cond, y=jnd, color = cond))+
+  geom_boxplot(outlier.shape = NA) + 
+  geom_jitter(position=position_jitter(0.2))+
+  labs(x = "Context", y = "JND (ms)") +
+  coord_cartesian(ylim = c(0,250))+
+  scale_color_manual(values = colors_plot) +
+  scale_fill_manual(values = colors_plot) + mytheme + theme(legend.position = "none")
+
+plot_jnd_EEG_exp1_box
+```
+
+##figure2
+```{r}
+figure1 <- cowplot::plot_grid(plot_fit_EEG_exp1, plot_pse_EEG_exp1_box, plot_jnd_EEG_exp1_box,  nrow = 1, labels = c("a", "b", "c"), rel_widths = c(4,3,3))
+
+if (saveFigure == TRUE){
+  ggsave("figures/figure2.png", figure1, width = 9, height = 4)
+}
+
+figure1
+```
+
+
+## Average fitted paramters and plot
+
+```{r}
+mpse_EEG_exp1 = par_EEG_exp1%>% 
+  group_by(cond) %>%
+  dplyr::summarise(m_pse = mean(pse), m_jnd = mean(jnd), n = n(),
+                   pse_se = sd(pse)/sqrt(n-1), 
+                   jnd_se = sd(jnd)/sqrt(n-1)) 
+
+mpse_EEG_exp1
+```
+
+
+## ANOVA 
+
+```{r}
+Anova_EEG_exp1_pse <- ezANOVA(data = par_EEG_exp1, dv= pse, wid=SubName, within=.(cond))
+Anova_EEG_exp1_pse
+```
+
+
+```{r}
+Anova_EEG_exp1_jnd <- ezANOVA(data = par_EEG_exp1, dv= jnd, wid=SubName, within=.(cond))
+Anova_EEG_exp1_jnd
+```
+
+
+
+
+# Experiment 2- Frequency Effect
+## Fit psychometric functions
+```{r}
+fits_EEG_exp2_v2 = quickpsy(dat_EEG_exp2, x = targetDur, k =RP, prob = .5,
+                            grouping = .(cond, SubName),
+                            fun = logistic_fun,
+                            thresholds = FALSE, bootstrap = 'none')
+```
+
+
+```{r, fig.cap= 'PSE: Bisection task performance in Experiment 2_v2'}
+# plot fitted function
+plot(fits_EEG_exp2_v2)+ mytheme
+```
+
+```{r}
+par_EEG_exp2_v2 <-fits_EEG_exp2_v2$par %>%  dplyr::select(cond, SubName, parn, par) %>%
+  pivot_wider(names_from = parn, values_from = par) %>% 
+  mutate(jnd = log(3)/p2)
+
+par_EEG_exp2_v2$pse <- par_EEG_exp2_v2$p1
+par_EEG_exp2_v2$p2 <-  NULL
+par_EEG_exp2_v2$p1 <-  NULL
+head(par_EEG_exp2_v2)
+```
+
+
+```{r}
+dat_EEG_exp2$cond= factor(dat_EEG_exp2$cond, levels=c("DF","AF"))
+fits_EEG_exp2_v2_all = quickpsy(dat_EEG_exp2, x = targetDur, k =RP, 
+                                grouping = .(cond),
+                                fun = logistic_fun,
+                                thresholds = FALSE, bootstrap = 'none')
+fits_EEG_exp2_v2_all
+```
+
+## plot averaged
+```{r, fig.cap= 'PSE: Bisection task performance in Experiment 2_v2'}
+# plot fitted function
+plot_fit_EEG_exp2_v2 <- plot(fits_EEG_exp2_v2_all) + mytheme+
+  labs(x = "Target interval (ms)", y = "Proportion of 'long' responses") +
+  scale_color_manual(values = colors_plot) +
+  scale_fill_manual(values = colors_plot) +
+  theme(legend.position = "top", legend.title = element_blank()) 
+
+plot_fit_EEG_exp2_v2
+```
+
+```{r}
+plot_pse_EEG_exp2_v2_box <- ggplot(par_EEG_exp2_v2%>% mutate(cond = factor(cond, levels = c("DF", "AF"))),aes(x=cond, y=pse, color = cond))+
+  geom_boxplot() + 
+  geom_jitter(position=position_jitter(0.2))+
+  geom_signif(comparisons=list(c("DF", "AF")), annotations="***", y_position = 1350, tip_length = 0, vjust=0.4, color = 'black') +
+  labs(x = "Context", y = "PSE (ms)") +
+  coord_cartesian(ylim = c(400,1500)) +
+  scale_color_manual(values = colors_plot) +
+  scale_fill_manual(values = colors_plot) +mytheme+theme(legend.position = "none")
+
+plot_pse_EEG_exp2_v2_box
+```
+
+```{r}
+plot_jnd_EEG_exp2_v2_box <- ggplot(par_EEG_exp2_v2%>% mutate(cond = factor(cond, levels = c("DF", "AF"))),aes(x=cond, y=jnd, color = cond))+
+  geom_boxplot() + 
+  geom_jitter(position=position_jitter(0.2))+
+  geom_signif(comparisons=list(c("DF", "AF")), annotations="**", y_position = 215, tip_length = 0, vjust=0.4, color = 'black')+mytheme+
+  labs(x = "Context", y = "JND (ms)") +
+  coord_cartesian(ylim = c(0,250))+
+  scale_color_manual(values = colors_plot) +
+  scale_fill_manual(values = colors_plot) + theme(legend.position = "none")
+
+plot_jnd_EEG_exp2_v2_box
+```
+##figure 4
+
+```{r}
+figure4 <- cowplot::plot_grid(plot_fit_EEG_exp2_v2, plot_pse_EEG_exp2_v2_box, plot_jnd_EEG_exp2_v2_box,  nrow = 1, labels = c("a", "b", "c"), rel_widths = c(4,3,3))
+if (saveFigure == TRUE){
+  ggsave("figures/figure6.png", figure4, width = 9, height = 4)
+}
+figure4
+```
+## Average fitted parameters and plot
+
+```{r}
+mpse_EEG_exp2_v2 = par_EEG_exp2_v2%>% 
+  group_by(cond) %>%
+  dplyr::summarise(m_pse = mean(pse), m_jnd = mean(jnd), n = n(),
+                   pse_se = sd(pse)/sqrt(n-1), 
+                   jnd_se = sd(jnd)/sqrt(n-1)) 
+mpse_EEG_exp2_v2
+```
+
+
+## ANOVA 
+
+```{r}
+Anova_EEG_exp2_v2_pse <- ezANOVA(data = par_EEG_exp2_v2, dv= pse, wid=SubName, within=.(cond))
+Anova_EEG_exp2_v2_pse
+```
+
+
+```{r}
+Anova_EEG_exp2_v2_jnd <- ezANOVA(data = par_EEG_exp2_v2, dv= jnd, wid=SubName, within=.(cond))
+Anova_EEG_exp2_v2_jnd
+```

derivates/EEG_bisection.Rproj

@@ -0,0 +1,13 @@
+Version: 1.0
+
+RestoreWorkspace: Default
+SaveWorkspace: Default
+AlwaysSaveHistory: Default
+
+EnableCodeIndexing: Yes
+UseSpacesForTab: Yes
+NumSpacesForTab: 2
+Encoding: UTF-8
+
+RnwWeave: Sweave
+LaTeX: pdfLaTeX

derivates/dataana.R

@@ -0,0 +1,67 @@
+# load packages and data
+library(tidyverse)
+library(knitr)
+# fitting psychometric function using quickpsy
+library(quickpsy)
+library(ez)
+library(ggplot2)
+library(ggsignif)
+library(PairedData)
+library(cowplot)
+library(effsize)
+library(dplyr)
+library(gtools)
+library(bayesplot)
+library(gridExtra)
+library(plyr)
+library(ggpubr)
+library(lmerTest)        # lmer function
+library(car)             # contr.Sum function
+library(sjPlot)          # tab_model function
+
+theme_set(bayesplot::theme_default())
+
+## define theme for whole report
+mytheme <-  theme_bw()+theme(panel.border = element_blank(),
+                              panel.grid.major = element_blank(), 
+                              panel.grid.minor = element_blank(),
+                              panel.background = element_blank(), 
+                              axis.line = element_line(colour = "black"))+
+  theme(strip.background = element_blank(),
+        axis.text.x = element_text(color = "black", size = 10), 
+        axis.text.y = element_text(color = "black", size = 10), 
+        axis.title.x = element_text(color = "black", size = 11, vjust = -1),
+        axis.title.y = element_text(color = "black", size = 11),
+        plot.title = element_text(color = "black", size = 11),
+        legend.text=element_text(size=9),axis.ticks.length=unit(.10, "cm"))
+
+colors_plot <- c("#185C9B","#D81E3D")
+
+##########################################################################
+# load behavior data
+dat_EEG_exp2 = read.csv('data/allData_exp2_beh.csv') 
+dat_EEG_exp1 = read.csv('data/allData_exp1_beh.csv')
+eeg_dat_all <- data.frame(rbind(dat_EEG_exp1, dat_EEG_exp2))
+
+# load ERP data
+allAverageDat_exp1 <- read.csv(paste0(getwd(), '/data/allAverageDat_cnv_exp1.csv'))
+allAverageDat_exp2 <- read.csv(paste0(getwd(), '/data/allAverageDat_cnv_exp2.csv'))
+
+data_ica_p2_exp1 <- read.csv(paste0(getwd(), '/data/allAverageDat_pos_exp1.csv'))
+data_ica_pos_exp2 <- read.csv(paste0(getwd(), '/data/allAverageDat_pos_exp2.csv'))
+
+# load peak data
+exp1_negPeak_lat <- read.csv(paste0(getwd(), '/data/peak_negative_cnv_exp1.csv'))
+exp2_negPeak_lat <- read.csv(paste0(getwd(), '/data/peak_negative_cnv_exp2.csv'))
+
+exp1_pos_peaks <- read.csv(paste0(getwd(), '/data/peak_positive_offset_exp1.csv'))
+exp2_pos_peaks <- read.csv(paste0(getwd(), '/data/peak_positive_offset_exp2.csv'))
+
+# CNV end-turning point data
+turn_exp1= read.csv("data/width_CNV_exp1.csv")
+turn_exp2= read.csv("data/width_CNV_exp2.csv")
+
+## calulate geometric mean
+gm_mean = function(x, na.rm=TRUE){
+  exp(sum(log(x[x > 0]), na.rm=na.rm) / length(x))
+}