Dateien hochladen nach 'Behavioral data + R-script'

Behavioral data + R-script/R-Skript_FrameBetween_publication.R

@@ -0,0 +1,397 @@
+
+#### FrameBetween ####
+
+
+# # Sind alle notwendigen Pakete installiert?
+# install.packages("car")
+# install.packages("ez")
+# install.packages("ggplot2")
+# install.packages("pastecs")
+# install.packages("tidyr")
+# install.packages("Hmisc")
+# install.packages("dplyr")
+# install.packages("tidyverse")
+# install.packages("gdata")
+# install.packages("lme4")
+# install.packages("psych")
+# install.packages("psychReport")
+# install.packages("rcompanion")
+# install.packages("data.table")
+# install.packages("hrbrthemes")
+# install.packages("viridis")
+library(ggplot2)
+library(psychReport)
+library(dplyr)
+library(apa)
+library(afex)
+library(dplyr)
+library(forcats)
+library(data.table)
+library(psych)
+library(doBy)
+library(rcompanion)
+library(effectsize)
+
+
+
+#### Read data ####
+
+
+setwd("Your Directory")    
+
+list_of_files <- list.files(path = "Your Directory", recursive = TRUE,
+                            pattern = "\\.dat$", 
+                            full.names = TRUE)
+
+FrameBetween <- vroom(list_of_files)
+
+FrameBetween <- subset(FrameBetween, FrameBetween$subjID !=2) # MEG artifacts
+
+# stim: 1 = anodal /excitatory    2 = Sham/Placebo
+FrameBetween$stim <- ifelse(FrameBetween$stim == 1, "Excitatory", "Sham")
+
+
+#### Check of Programming ####
+table  (FrameBetween$InitRew)
+table  (FrameBetween$RisktoLoose)
+table  (FrameBetween$GainLossTr)
+table  (FrameBetween$subjID)
+table  (FrameBetween$choice)
+describe(FrameBetween$choice)
+table  (FrameBetween$stim)
+kreuztabelle <- xtabs (~ FrameBetween$RisktoLoose + FrameBetween$GainLossTr)
+xtabs (~ FrameBetween$stim + FrameBetween$subjID)
+print(kreuztabelle)
+
+table(FrameBetween$subjID)
+
+# Half of the trials risky choices?
+observed_proportion <- c(10808, 12232)
+theoretical_proportion <- c(1/2, 1/2)
+chisq.test(observed_proportion, p=theoretical_proportion)
+
+FrameBetween = FrameBetween %>%
+  dplyr::group_by(subjID) %>%
+  dplyr::mutate(TrialNumber = (row_number()-1) %% 320 +1) %>%
+  ungroup()
+
+
+
+
+#### Add expected value ####
+FrameBetween <- FrameBetween %>% 
+  mutate(expectation_valueRisk = case_when(RisktoLoose == 20 & InitRew == 100 ~ 80,
+                                           RisktoLoose == 40 & InitRew == 100 ~ 60,
+                                           RisktoLoose == 60 & InitRew == 100 ~ 40,
+                                           RisktoLoose == 80 & InitRew == 100 ~ 20,
+                                           RisktoLoose == 20 & InitRew == 75 ~ 60,
+                                           RisktoLoose == 40 & InitRew == 75 ~ 45,
+                                           RisktoLoose == 60 & InitRew == 75 ~ 30,
+                                           RisktoLoose == 80 & InitRew == 75 ~ 15,
+                                           RisktoLoose == 20 & InitRew == 50 ~ 40,
+                                           RisktoLoose == 40 & InitRew == 50 ~ 30,
+                                           RisktoLoose == 60 & InitRew == 50 ~ 20,
+                                           RisktoLoose == 80 & InitRew == 50 ~ 10,
+                                           RisktoLoose == 20 & InitRew == 25 ~ 20,
+                                           RisktoLoose == 40 & InitRew == 25 ~ 15,
+                                           RisktoLoose == 60 & InitRew == 25 ~ 10,
+                                           RisktoLoose == 80 & InitRew == 25 ~ 5,
+  ))
+
+FrameBetween <- FrameBetween %>% 
+  mutate(expectation_valueSafe = case_when(InitRew == 100 ~ 40,
+                                           InitRew == 75 ~ 30,
+                                           InitRew == 50 ~ 20,
+                                           InitRew == 25 ~ 10,
+  ))
+
+FrameBetween <- FrameBetween %>% 
+  mutate(expectation_valueChoice = case_when(choice == 1 ~ expectation_valueRisk,
+                                             choice == 0 ~ expectation_valueSafe,
+  ))
+
+
+FrameBetween <- FrameBetween %>% 
+  mutate(RationalDecision = case_when(choice == 1 & expectation_valueRisk >= expectation_valueSafe  ~ 1,
+                                      choice == 0 & expectation_valueSafe >= expectation_valueRisk ~ 1,
+                                      choice == 1 & expectation_valueRisk < expectation_valueSafe  ~ 0,
+                                      choice == 0 & expectation_valueSafe < expectation_valueRisk ~ 0,
+  ))
+
+
+FrameBetween <- FrameBetween %>% 
+  mutate(Outcome = case_when(choice == 1 & RiskFeedback == 22 | choice == 1 &  RiskFeedback == 52 | choice == 1 & RiskFeedback == 72 | choice == 1 & RiskFeedback == 102  ~ 1,
+                             choice == 1 & RiskFeedback == 21 | choice == 1 & RiskFeedback == 51 | choice == 1 & RiskFeedback == 71 | choice == 1 & RiskFeedback == 101  ~ 0,
+                             choice == 0 & GainLossTr == 0 & RiskFeedback == 0  ~ 1, # 99 oder 1
+                             choice == 0 & GainLossTr == 1 & RiskFeedback == 0 ~ 0 # 99 oder 0
+  ))
+
+FrameBetween$choicepercent <- FrameBetween$choice*100
+
+#FrameBetween <- FrameBetween %>% 
+#  mutate(Outcome4 = case_when(choice == 1 & RiskFeedback == 22 | choice == 1 &  RiskFeedback == 52 | choice == 1 & RiskFeedback == 72 | choice == 1 & RiskFeedback == 102  ~ 4,
+#                             choice == 1 & RiskFeedback == 21 | choice == 1 & RiskFeedback == 51 | choice == 1 & RiskFeedback == 71 | choice == 1 & RiskFeedback == 101  ~ 1,
+#                             choice == 0 & GainLossTr == 0 & RiskFeedback == 0  ~ 3, 
+#                             choice == 0 & GainLossTr == 1 & RiskFeedback == 0 ~ 2 
+#  ))
+
+FrameBetween <- data.table(FrameBetween)
+FrameBetween[ , previousTr := shift(Outcome)]
+#FrameBetween[ , previousTr4 := shift(Outcome4)]
+
+PrevTrial <- subset(FrameBetween, FrameBetween$previousTr !=99)
+PrevTrial <- subset(PrevTrial, PrevTrial$previousTr !="NA")
+PrevTrial <- subset(PrevTrial, PrevTrial$TrialNumber !=1)
+#PrevTrial4 <- subset(FrameBetween, FrameBetween$previousTr4 !=99)
+#PrevTrial4 <- subset(PrevTrial, PrevTrial$previousTr4 !="NA")
+PrevTrial$previousTr <- ifelse(PrevTrial$previousTr == 1, "Gain", "Loss")
+
+
+# Choice = 0 --> no risk, Choice = 1 --> risk
+
+# GainLossTr: 0 = GainTrial    1 = LossTrial
+
+FrameRat <- subset(FrameBetween, FrameBetween$RisktoLoose !=60)
+
+GLMRat <- glm(formula = RationalDecision ~ stim, family = binomial(logit), data = FrameRat)
+summary(GLMRat)
+
+GLMFrame <- glm(formula = choice ~ GainLossTr, family = binomial(logit), data = FrameBetween)
+summary(GLMFrame)
+effectsizes <- standardize_parameters(GLMFrame, exp = TRUE)
+
+GLMIntFrame <- glm(formula = choice ~ stim * GainLossTr, family = binomial(logit), data = FrameBetween)
+summary(GLMIntFrame)
+effectsizes <- standardize_parameters(GLMIntFrame, exp = TRUE)
+
+GLMALL <- glm(formula = choice ~ stim * GainLossTr * RisktoLoose, family = binomial(logit), data = FrameBetween)
+summary(GLMALL)
+effectsizes <- standardize_parameters(GLMALL, exp = TRUE)
+
+
+GLMInitRew <- glm(formula = choice ~ stim * InitRew, family = binomial(logit), data = FrameBetween)
+summary(GLMInitRew)
+
+
+Frame_Sham <- subset(FrameBetween, FrameBetween$stim !="Excitatory")
+Frame_Ano <- subset(FrameBetween, FrameBetween$stim !="Sham")
+
+GLMInitRew <- glm(formula = choice ~ InitRew, family = binomial(logit), data = Frame_Sham)
+summary(GLMInitRew)
+effectsizes <- standardize_parameters(GLMInitRew, exp = TRUE)
+
+GLMInitRew <- glm(formula = choice ~ InitRew, family = binomial(logit), data = Frame_Ano)
+summary(GLMInitRew)
+
+GLMprevTr <- glm(formula = choice ~ stim * previousTr, family = binomial(logit), data = PrevTrial)
+summary(GLMprevTr)
+effectsizes <- standardize_parameters(GLMprevTr, exp = TRUE)
+
+Frame_Sham <- subset(PrevTrial, PrevTrial$stim !="Excitatory")
+Frame_Ano <- subset(PrevTrial, PrevTrial$stim !="Sham")
+
+Frame_prevLoss <- subset(PrevTrial, PrevTrial$previousTr !="Gain")
+Frame_prevGain <- subset(PrevTrial, PrevTrial$previousTr !="Loss")
+
+
+GLMprevTr <- glm(formula = choice ~ previousTr, family = binomial(logit), data = Frame_Ano)
+summary(GLMprevTr)
+effectsizes <- standardize_parameters(GLMprevTr, exp = TRUE)
+
+GLMprevTr <- glm(formula = choice ~ previousTr, family = binomial(logit), data = Frame_Sham)
+summary(GLMprevTr)
+effectsizes <- standardize_parameters(GLMprevTr, exp = TRUE)
+
+GLMprevTr <- glm(formula = choice ~ stim, family = binomial(logit), data = Frame_prevLoss)
+summary(GLMprevTr)
+
+GLMprevTr <- glm(formula = choice ~ stim, family = binomial(logit), data = Frame_prevGain)
+summary(GLMprevTr)
+
+GLMNull <- glm(formula = choice ~ 1, family = binomial(logit), data = FrameBetween)
+summary(GLMNull)
+
+anova(GLMFrame, GLMNull, test = "Chisq")
+anova(GLMRat, GLMNull, test = "Chisq")
+anova(GLMALL, GLMNull, test = "Chisq")
+
+
+
+#### Post-Hocs Interaction Risk and Expected Value ####
+
+
+apa(t.test(expectation_valueChoice ~ stim, data = FrameBetween, var.equal = TRUE, paired = FALSE,
+           alternative = "greater")) # two.sided
+
+groupwiseMean(expectation_valueChoice ~ stim,
+              data = FrameBetween,
+              conf = 0.95,
+              digits = 3)
+
+FrameAgg <- aggregate(choice ~ subjID + stim + RisktoLoose,
+                      data = FrameBetween, FUN = mean)
+
+
+Frame20 <- subset(FrameBetween, FrameBetween$RisktoLoose !=80)
+Frame20 <- subset(Frame20, Frame20$RisktoLoose !=60)
+Frame20 <- subset(Frame20, Frame20$RisktoLoose !=40)
+
+apa(t.test(expectation_valueChoice ~ stim, data = Frame20, var.equal = TRUE, paired = FALSE,
+                       alternative = "greater")) # two.sided
+
+GLM20 <- glm(formula = choice ~ stim, family = binomial(logit), data = Frame20)
+summary(GLM20)
+effectsizes <- standardize_parameters(GLM20, exp = TRUE)
+
+xtabs (~ Frame20$choice + Frame20$stim)
+
+prop.test(x = c(1089, 1001), n = c(1280, 1280),alternative = "two.sided")
+
+Frame40 <- subset(FrameBetween, FrameBetween$RisktoLoose !=80)
+Frame40 <- subset(Frame40, Frame40$RisktoLoose !=60)
+Frame40 <- subset(Frame40, Frame40$RisktoLoose !=20)
+
+apa(t.test(expectation_valueChoice ~ stim, data = Frame40, var.equal = TRUE, paired = FALSE,
+           alternative = "greater")) # two.sided
+
+xtabs (~ Frame40$choice + Frame40$stim)
+
+GLM40 <- glm(formula = choice ~ stim, family = binomial(logit), data = Frame40)
+summary(GLM40)
+
+prop.test(x = c(852, 867), n = c(1280, 1280),alternative = "two.sided")
+
+Frame60 <- subset(FrameBetween, FrameBetween$RisktoLoose !=80)
+Frame60 <- subset(Frame60, Frame60$RisktoLoose !=20)
+Frame60 <- subset(Frame60, Frame60$RisktoLoose !=40)
+
+apa(t.test(expectation_valueChoice ~ stim, data = Frame60, var.equal = TRUE, paired = FALSE,
+           alternative = "greater"))
+
+xtabs (~ Frame60$choice + Frame60$stim)
+
+GLM60 <- glm(formula = choice ~ stim, family = binomial(logit), data = Frame60)
+summary(GLM60)
+effectsizes <- standardize_parameters(GLM60, exp = TRUE)
+
+prop.test(x = c(468, 637), n = c(1280, 1280),alternative = "two.sided")
+
+Frame80 <- subset(FrameBetween, FrameBetween$RisktoLoose !=60)
+Frame80 <- subset(Frame80, Frame80$RisktoLoose !=20)
+Frame80 <- subset(Frame80, Frame80$RisktoLoose !=40)
+
+apa(t.test(expectation_valueChoice ~ stim, data = Frame80, var.equal = TRUE, paired = FALSE,
+           alternative = "greater"))
+
+xtabs (~ Frame80$choice + Frame80$stim)
+
+GLM80 <- glm(formula = choice ~ stim, family = binomial(logit), data = Frame80)
+summary(GLM80)
+effectsizes <- standardize_parameters(GLM80, exp = TRUE)
+interpret_oddsratio(1.88, rules = "cohen1988")
+
+prop.test(x = c(223, 361), n = c(1280, 1280),alternative = "two.sided")
+
+apa(t.test(expectation_valueChoice ~ stim, data = FrameRat, var.equal = TRUE, paired = FALSE,
+           alternative = "greater")) # two.sided
+
+#### Post-hoc Tests Interaction Frame ####
+LossFrame <- subset(FrameBetween, FrameBetween$GainLossTr !="GainTr")
+GainFrame <- subset(FrameBetween, FrameBetween$GainLossTr !="LossTr")
+
+GainFrame$Diff <- GainFrame$choice - LossFrame$choice
+
+bar <- ggplot(GainFrame, aes(stim, Diff, fill = stim))
+bar + 
+  stat_summary(fun = mean, geom = "bar", position = "dodge") +
+  stat_summary(fun.data = mean_cl_normal, geom = "errorbar", position = position_dodge(width = 0.90),width = 0.2) +
+  labs (x = "Frame", y = "Risky Choice", fill = "Stim") +
+  scale_fill_brewer(palette = "Set1")
+  theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(),
+      panel.background = element_blank(), axis.line = element_line(colour = "black"))
+
+xtabs (~ GainFrame$choice + GainFrame$stim)
+
+prop.test(x = c(1293, 1384), n = c(2560, 2560),alternative = "two.sided")
+
+xtabs (~ LossFrame$choice + LossFrame$stim)
+
+prop.test(x = c(1339, 1482), n = c(2560, 2560),alternative = "two.sided")
+
+GainFrameexc <- subset(GainFrame, GainFrame$stim !="Sham")
+GainFrameSham <- subset(GainFrame, GainFrame$stim !="Excitatory")
+
+table(GainFrameexc$Diff)
+table(GainFrameSham$Diff)
+
+prop.test(x = c(640, 680), n = c(2560, 2560),alternative = "less")
+
+#### Post-Hoc Interaction InitRew ####
+
+Frame25 <- subset(FrameBetween, FrameBetween$InitRew !=50)
+Frame25 <- subset(Frame25, Frame25$InitRew !=75)
+Frame25 <- subset(Frame25, Frame25$InitRew !=100)
+
+xtabs (~ Frame25$choice + Frame25$stim)
+
+prop.test(x = c(652, 647), n = c(1280, 1280),alternative = "two.sided")
+
+Frame50 <- subset(FrameBetween, FrameBetween$InitRew !=25)
+Frame50 <- subset(Frame50, Frame50$InitRew !=75)
+Frame50 <- subset(Frame50, Frame50$InitRew !=100)
+
+xtabs (~ Frame50$choice + Frame50$stim)
+
+prop.test(x = c(658, 692), n = c(1280, 1280),alternative = "two.sided")
+
+Frame75 <- subset(FrameBetween, FrameBetween$InitRew !=50)
+Frame75 <- subset(Frame75, Frame75$InitRew !=25)
+Frame75 <- subset(Frame75, Frame75$InitRew !=100)
+
+xtabs (~ Frame75$choice + Frame75$stim)
+
+prop.test(x = c(652, 728), n = c(1280, 1280),alternative = "two.sided")
+
+GLM75 <- glm(formula = choice ~ stim, family = binomial(logit), data = Frame75)
+summary(GLM75)
+effectsizes <- standardize_parameters(GLM75, exp = TRUE)
+
+Frame100 <- subset(FrameBetween, FrameBetween$InitRew !=25)
+Frame100 <- subset(Frame100, Frame100$InitRew !=50)
+Frame100 <- subset(Frame100, Frame100$InitRew !=75)
+
+xtabs (~ Frame100$choice + Frame100$stim)
+
+prop.test(x = c(670, 799), n = c(1280, 1280),alternative = "two.sided")
+
+GLM100 <- glm(formula = choice ~ stim, family = binomial(logit), data = Frame100)
+summary(GLM100)
+effectsizes <- standardize_parameters(GLM100, exp = TRUE)
+
+
+#### Post-Hoc previous Trial ####
+prvLoss <- subset(PrevTrial, PrevTrial$previousTr !="Gain")
+prvGain <- subset(PrevTrial, PrevTrial$previousTr !="Loss")
+Exc <- subset(PrevTrial, PrevTrial$stim !="Sham")
+Sham <- subset(PrevTrial, PrevTrial$stim !="Excitatory")
+
+xtabs (~ prvLoss$choice + prvLoss$stim)
+
+prop.test(x = c(1126, 1262), n = c(2217, 2181),alternative = "two.sided")
+
+xtabs (~ prvGain$choice + prvGain$stim)
+
+prop.test(x = c(1454, 1495), n = c(2787, 2758),alternative = "two.sided")
+
+xtabs (~ Exc$choice + Exc$previousTr)
+
+prop.test(x = c(1454, 1126), n = c(2787, 2217),alternative = "two.sided")
+
+xtabs (~ Sham$choice + Sham$previousTr)
+
+prop.test(x = c(1496, 1260), n = c(2758, 2178),alternative = "two.sided")
+
+
+
+
+