Hippocampal dependence of value-based decisions: Statistical Analyses

Load workspace

load("/Users/zeynepenkavi/Dropbox/CDSPart1/TransitivityOpen/Transitivity_OpenDataOrganizationWS052817.RData")

Group differences in intransitivies

Summary Statistics

both.Intransitive %>%
  group_by(Task, Group) %>%
  summarise(mean_CleanIntr = mean(CleanIntr),
            median_CleanIntr = median(CleanIntr),
            sd_CleanIntr = sd(CleanIntr),
            mean_CleanPercentIntr = mean(CleanPercentIntr),
            median_CleanPercentIntr = median(CleanPercentIntr),
            sd_CleanPercentIntr = sd(CleanPercentIntr))

both.Intransitive %>% 
  group_by(Task, Group) %>% 
  summarise(MeanCleanPercentIntr = mean(CleanPercentIntr, na.rm=T),
            SeCleanPercentIntr = sem(CleanPercentIntr),
            eb.low = MeanCleanPercentIntr - SeCleanPercentIntr,
            eb.high = MeanCleanPercentIntr + SeCleanPercentIntr) %>%
ggplot(aes(x = Group, y=MeanCleanPercentIntr, group = Task))+
  geom_point(aes(shape = Task))+
  geom_line(aes(linetype = Task))+
  geom_errorbar(aes(ymax = eb.high, ymin=eb.low), width=0.25)+
  ylab("Mean Percentage of Intransitivities")+
  xlab("")+
  theme_bw()+
  theme(axis.title.y = element_text(face="bold", size = 14),
        axis.text.x  = element_text(face="bold", size = 14),
        legend.text = element_text(size = 14),
        legend.title = element_text(size = 14))+
  scale_shape_discrete(breaks = c("choice", "numbers"),
                       labels = c("Value-based", "Numbers"))+
  scale_linetype_discrete(breaks = c("choice", "numbers"),
                          labels = c("Value-based", "Numbers"))

ggsave("/Users/zeynepenkavi/Dropbox/CDSPart1/TransitivityOpen/CleanFigure2.png", width=6, height=4, dpi=300)

Looking at the plot we are interested in whether the difference between the two lines is the same for all groups. In other words: Is there a significant intereaction between task and group, i.e. is the effect of changing the task same for the MTL group as it is for the control groups.

To choose the correct model that would check for statistical differences between the groups we checked the assumptions for each possibility. Though all methods yield essentially the same answer we chose to use a linear mixed model with orthogonal contrasts because it was the most appropriate method.

The first thing we could do would be a one-way ANOVA for both tasks checking whether the main dependent measure (percentage of intransitivities) differs depending on the group and task. The distribution of effect sizes within tasks, however, is not normal for either task, which is why this would not be an appropriate test.

both.Intransitive %>%
  ggplot(aes(CleanPercentIntr, fill=Group))+
  geom_histogram(position="identity", alpha=0.5)+
  theme_bw()+
  xlab("Percentage of intransitivity")+
  facet_wrap(~Task)

shapiro.test(both.Intransitive[both.Intransitive$Task == "numbers",]$CleanPercentIntr)

    Shapiro-Wilk normality test

data:  both.Intransitive[both.Intransitive$Task == "numbers", ]$CleanPercentIntr
W = 0.36, p-value <2e-16

shapiro.test(both.Intransitive[both.Intransitive$Task == "choice",]$CleanPercentIntr)

    Shapiro-Wilk normality test

data:  both.Intransitive[both.Intransitive$Task == "choice", ]$CleanPercentIntr
W = 0.69, p-value = 1e-12

To solve this we could try a non-parametric test (e.g. Kruskall-Wallis) but this wouldn’t allow us to check the interaction of group by task and more importantly the variance in the effect size distributions are unequal. That is, as the Bartlett test below and the plot above shows variance in effect size is not independent of task.

bartlett.test(CleanPercentIntr ~ Task, data = both.Intransitive)

    Bartlett test of homogeneity of variances

data:  CleanPercentIntr by Task
Bartlett's K-squared = 58, df = 1, p-value = 3e-14

To deal with this issue of task dependence of the variance we log transform the effect sizes.

bartlett.test(log(CleanPercentIntr + 1) ~ Task, data = both.Intransitive)

    Bartlett test of homogeneity of variances

data:  log(CleanPercentIntr + 1) by Task
Bartlett's K-squared = 3.4, df = 1, p-value = 0.07

Using transformed dependent variables we first confirm that the interactive model is the best fitting model by comparing it to nested simpler models.

gm1a <- lm(log(CleanPercentIntr+1) ~ Task, both.Intransitive)#; summary(gm1a)
gm1b <- lm(log(CleanPercentIntr+1) ~ Group, both.Intransitive)# ; summary(gm1b)
gm2 <- lm(log(CleanPercentIntr+1) ~ Group+Task, both.Intransitive)# ; summary(gm2)
gm3 <- lm(log(CleanPercentIntr+1) ~ Group*Task, both.Intransitive)# ; summary(gm3)

anova(gm1a, gm2, gm3)

anova(gm1b, gm2, gm3)

So what does this interactive model say? To get a better understanding of this we first look at the contrasts used in this model.

contrasts(both.Intransitive$Task)

        numbers
choice        0
numbers       1

contrasts(both.Intransitive$Group)

    ETL MTL
C     0   0
ETL   1   0
MTL   0   1

Based on these contrasts the group differences in this model compare both patient groups to the healthy control alone. Accordingly the parameters of the linear model would have the following interpretations.

summary(gm3)

Call:
lm(formula = log(CleanPercentIntr + 1) ~ Group * Task, data = both.Intransitive)

Residuals:
    Min      1Q  Median      3Q     Max 
-1.2678 -0.3590 -0.0785  0.2261  2.1345 

Coefficients:
                     Estimate Std. Error t value Pr(>|t|)    
(Intercept)            1.2508     0.0839   14.90  < 2e-16 ***
GroupETL               0.1011     0.1187    0.85    0.395    
GroupMTL               0.5441     0.1178    4.62  7.4e-06 ***
Tasknumbers           -1.1322     0.1187   -9.54  < 2e-16 ***
GroupETL:Tasknumbers   0.1422     0.1679    0.85    0.398    
GroupMTL:Tasknumbers  -0.3019     0.1665   -1.81    0.072 .  
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 0.46 on 176 degrees of freedom
Multiple R-squared:  0.655, Adjusted R-squared:  0.645 
F-statistic: 66.7 on 5 and 176 DF,  p-value: <2e-16

b0 = mean percentage of intransitivity for both tasks for the C group is >0 (not interesting) b1 = ETL group doesn’t make significantly more intrans in CHOICE task! (good) - simple effect b2 = MTL group DOES make sig more intrans in CHOICE task! (good) - simple effect b3 = C group makes sig less intrans in NUMBERS task! (not surprising, easy task) b4 = Is the effect of task in ETL group same as C. Yes. b5 = Is the effect of task in MTL group same as C. Marginally no.

As mentioned above, however, these are not orthogonal contrasts and are comparing both groups to Healthy controls only! But we are interested in comparing the MTL group to both control groups. To answer this question and to avoid overlapping variance we orthogonalize the contrasts. This would results in first comparing the ETL group to the C group and then the MTL group to BOTH CONTROL GROUPS!

#Change the contrasts, make them orthogonal
both.Intransitive.contrast <- both.Intransitive
contrasts(both.Intransitive.contrast$Group) <- cbind(c(-1,1, 0), c(-1, -1, 2))
contrasts(both.Intransitive.contrast$Task) <- c(-1, 1)

gm3.con <- lm(log(CleanPercentIntr+1) ~ Group*Task, both.Intransitive.contrast) 
summary(gm3.con)

Call:
lm(formula = log(CleanPercentIntr + 1) ~ Group * Task, data = both.Intransitive.contrast)

Residuals:
    Min      1Q  Median      3Q     Max 
-1.2678 -0.3590 -0.0785  0.2261  2.1345 

Coefficients:
             Estimate Std. Error t value Pr(>|t|)    
(Intercept)    0.8732     0.0341   25.62  < 2e-16 ***
Group1         0.0861     0.0420    2.05    0.042 *  
Group2         0.1023     0.0240    4.27 0.000032 ***
Task1         -0.5927     0.0341  -17.39  < 2e-16 ***
Group1:Task1   0.0356     0.0420    0.85    0.398    
Group2:Task1  -0.0622     0.0240   -2.59    0.010 *  
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 0.46 on 176 degrees of freedom
Multiple R-squared:  0.655, Adjusted R-squared:  0.645 
F-statistic: 66.7 on 5 and 176 DF,  p-value: <2e-16

What does the model say now:

b0 = Unweighted grand mean (mean of group means) (sum(aggregate(CleanPercentIntr ~ Group, both.Intransitive.contrast, mean)[,2])/3) b1 = Total group effect for ETL compared to healthy controls only (quantified) b2 = Total group effect for MTL compared to both control groups (quantified) b3 = Total task effect b4 = Is the effect of task change same for ETL and C groups? (yes.) b5 = Is the effect of task change same for MTL vs both control groups? (No!)

Since there are two measures per subject we check if a mixed model necessary comparing the log-likelihoods of a model with random effects to one without them. We find that that a model with random effects is indeed significantly better fit to data.

gm3.lmer.con <- lmerTest::lmer(log(CleanPercentIntr+1) ~ Group*Task + (1|f.id), both.Intransitive.contrast, REML = F)

x2 = -2*logLik(gm3.con) +2*logLik(gm3.lmer.con)
x2

'log Lik.' 4.39 (df=7)

pchisq(x2, df=1, lower.tail=F)

'log Lik.' 0.03614 (df=7)

summary(gm3.lmer.con) 

Linear mixed model fit by maximum likelihood t-tests use Satterthwaite
  approximations to degrees of freedom [lmerMod]
Formula: log(CleanPercentIntr + 1) ~ Group * Task + (1 | f.id)
   Data: both.Intransitive.contrast

     AIC      BIC   logLik deviance df.resid 
   239.2    264.8   -111.6    223.2      174 

Scaled residuals: 
   Min     1Q Median     3Q    Max 
-2.442 -0.601 -0.152  0.543  4.294 

Random effects:
 Groups   Name        Variance Std.Dev.
 f.id     (Intercept) 0.0444   0.211   
 Residual             0.1601   0.400   
Number of obs: 182, groups:  f.id, 91

Fixed effects:
             Estimate Std. Error      df t value Pr(>|t|)    
(Intercept)    0.8732     0.0370 91.0000   23.61  < 2e-16 ***
Group1         0.0861     0.0455 91.0000    1.89  0.06175 .  
Group2         0.1023     0.0260 91.0000    3.93  0.00016 ***
Task1         -0.5927     0.0297 91.0000  -19.98  < 2e-16 ***
Group1:Task1   0.0356     0.0365 91.0000    0.97  0.33288    
Group2:Task1  -0.0622     0.0209 91.0000   -2.98  0.00369 ** 
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Correlation of Fixed Effects:
            (Intr) Group1 Group2 Task1  Gr1:T1
Group1       0.000                            
Group2      -0.015  0.000                     
Task1        0.000  0.000  0.000              
Group1:Tsk1  0.000  0.000  0.000  0.000       
Group2:Tsk1  0.000  0.000  0.000 -0.015  0.000

The linear mixed model with orthogonal contrasts is the most appropriate analysis and reveals a significant group task interaction in the percentage of intransitivities. Accordingly the MTL group had a higher percentage of intransitivities in the choice task than the numbers task compared to both control groups, while the difference between the percentage of intransitivities did not differ between tasks between the health controls and the ETL group.

Lesion size and intransitivies

symmetry.data <- read.xlsx("/Users/zeynepenkavi/Dropbox/CDSPart1/Transitivity/TransitivityData/Symmetry_index_zeynep.xlsx", 1)

symmetry.data = both.Intransitive %>%
  filter(Task == 'choice') %>%
  mutate(ID = as.numeric(as.character(f.id))) %>%
  right_join(symmetry.data, by= 'ID') %>%
  drop_na()

Spearman correlation between instransitivity percentage and lesion size

cor(symmetry.data$Symmetry.index, symmetry.data$CleanPercentIntr, method = "spearman")

[1] 0.676

Significance test for the Spearman correlation

spearman.test(symmetry.data$Symmetry.index, symmetry.data$CleanPercentIntr)

  Rsquare         F       df1       df2    pvalue         n 
 0.456969 11.781218  1.000000 14.000000  0.004044 16.000000 

Plot of correlation between

ggplot(symmetry.data, aes(Symmetry.index, CleanPercentIntr))+
  geom_point()+
  geom_smooth(method = 'loess', span = 2, color = "black", alpha = 0.1, linetype = "dashed")+
  theme_bw()+
  ylab("Percentage of Intransitve Choice")+
  xlab("Compromised Hippocampal Ratio")+
  theme(axis.title.y = element_text(face="bold", size = 14),
        axis.title.x  = element_text(face="bold", size = 14),
        legend.text = element_text(size = 14),
        legend.title = element_text(size = 14))

ggsave("/Users/zeynepenkavi/Dropbox/CDSPart1/TransitivityOpen/CleanFigure3.png", width=4, height=4, dpi=300)

Simulations

path = '/Users/zeynepenkavi/Dropbox/CDSPart1/TransitivityOpen/sim_output/'

fileList = list.files(path,pattern = "*_hardmax.csv")

read_csv_fn <- function(filename){
    ret <- read.csv(filename)
    ret$Source <- filename 
    ret
}

lapply(paste0(path, fileList), read_csv_fn)  %>% 
  bind_rows() %>%
  separate(Source, c("a", "b", "c", "d", "e", "f","g", "sim"),sep = "/") %>%
  select(pct_Intrans, sim) %>%
  separate(sim, c("a", "b", "c", "d", "noise", "e", "f"), sep = "_") %>%
  select(pct_Intrans, noise) %>%
  filter(noise %in% c(0.05, 0.25) == FALSE) %>%
  group_by(noise) %>%
  summarise(mean_Intrans = mean(pct_Intrans),
            sem_Intrans = sem(pct_Intrans)) %>%
  ggplot(aes(noise, mean_Intrans))+
  geom_point()+
  theme_bw()+
  geom_errorbar(aes(ymin = mean_Intrans-sem_Intrans, ymax = mean_Intrans+sem_Intrans), width=0.25)+
  ylab("Percentage of Intransivity")+
  xlab("Noise level")#+

  # ggtitle("Results of simulation using hardmax")

fileList = list.files(path,pattern = "*_rbinom.csv")

lapply(paste0(path, fileList), read_csv_fn)  %>% 
  bind_rows() %>%
  separate(Source, c("a", "b", "c", "d", "e", "f","g", "sim"),sep = "/") %>%
  select(pct_Intrans, sim) %>%
  separate(sim, c("a", "b", "c", "d", "noise", "e", "f"), sep = "_") %>%
  select(pct_Intrans, noise) %>%
  filter(noise %in% c(0.05, 0.25) == FALSE) %>%
  group_by(noise) %>%
  summarise(mean_Intrans = mean(pct_Intrans),
            sem_Intrans = sem(pct_Intrans)) %>%
  ggplot(aes(noise, mean_Intrans))+
  geom_point()+
  theme_bw()+
  geom_errorbar(aes(ymin = mean_Intrans-sem_Intrans, ymax = mean_Intrans+sem_Intrans), width=0.25)+
  ylab("Percentage of Intransivity")+
  xlab("Noise level")+
  scale_y_continuous(breaks=seq(0,25,5), limits=c(0, 26))+
  ggtitle("Results of simulation flipping coin with noisy choice p")

fileList = list.files(path,pattern = "*_epsilon.csv")

lapply(paste0(path, fileList), read_csv_fn)  %>% 
  bind_rows() %>%
  separate(Source, c("a", "b", "c", "d", "e", "f","g", "sim"),sep = "/") %>%
  select(pct_Intrans, sim) %>%
  separate(sim, c("a", "b", "c", "d", "noise", "e", "f"), sep = "_") %>%
  select(pct_Intrans, noise) %>%
  filter(noise %in% c(0.05, 0.25) == FALSE) %>%
  group_by(noise) %>%
  summarise(mean_Intrans = mean(pct_Intrans),
            sem_Intrans = sem(pct_Intrans)) %>%
  ggplot(aes(noise, mean_Intrans))+
  geom_point()+
  theme_bw()+
  geom_errorbar(aes(ymin = mean_Intrans-sem_Intrans, ymax = mean_Intrans+sem_Intrans), width=0.25)+
  ylab("Percentage of Intransivity")+
  xlab("Noise level")+
  scale_y_continuous(breaks=seq(0,25,5), limits=c(0, 26))+
  ggtitle("Results of simulation with epsilon greedy")

Alternative explanations of intransitivity

Episodic memory

One alternative explanation of the higher level of intransitivities in the MTL group is about episodic memory deficits. Accordingly, the MTL group might be showing higher rates of intransitivity because they are unable to recall choices they have made earlier in the task. This hypothesis would predict a pattern of proportion of intransitive triplets each trial was involved in that increases towards the end of the task (the later the trial the more intransitive triplets it was involved in). The figure below shows the empirical change of intransitive triplet involvement with trial number.

choice2.trial.data$Task <- "choice"
numbers.trial.data$Task <- "numbers"

both.trial.data <- rbind(choice2.trial.data[,c("f.id", "Group", "RT", "Task", "Trialnumber", "IntransTripleCounted")], numbers.trial.data[,c("f.id", "Group", "RT", "Task", "Trialnumber", "IntransTripleCounted")])

both.trial.data = both.Intransitive %>%
  select(f.id, CleanIntr, Task) %>%
  right_join(both.trial.data, by=c("f.id", "Task")) %>%
  mutate(PropInIntransTriplet = IntransTripleCounted/(CleanIntr*3)*100,
         Task = factor(Task, levels = c("choice", "numbers"), labels = c("Value-based", "Perceptual")))%>%
  filter(RT>0, RT<6000)

Note that most of the trials are not involved in any intransitivities so the averaged proportion of intransitivities for each trial number is pretty small (i.e why the y-axis is so small).

both.trial.data %>%
  ggplot(aes(x = Trialnumber, y = PropInIntransTriplet, group = Group, color = Group))+
  geom_smooth(method = "loess", span=5) +
  theme_bw()+
  facet_grid(~Task)+
  ylab("Percent of intransitive Triplets involved") +
  xlab("Order choice was made")

To test if there are any statistical differences both between the groups and throughout each task and between the two tasks we first compare three models where the proportion of intransitive triplets a trial is involved is regressed over the trial number and experimental group. The first model is additive, the second interactive and the third includes a quadratic effect of trial number as well. Model comparison shows that the addition of the interaction terms and the quadratic effect do not lead to significant improvements in fit. Thus we report the effects in the simplest additive model.

mem1 = lmer(PropInIntransTriplet ~ (1|f.id)+scale(Trialnumber)+Group*Task, both.trial.data)

#adding interaction of linear effect of RT between tasks and groups
mem2 = lmer(PropInIntransTriplet ~ (1|f.id)+scale(Trialnumber)*Group*Task, both.trial.data)

#adding the quadratic effect
mem3 = lmer(PropInIntransTriplet ~ (1|f.id)+scale(poly(Trialnumber,2))*Task*Group, both.trial.data)

anova(mem1, mem2, mem3)

The additive model suggests that contrary to the prediction of an episodic memory during the task account there is both no change in the proportion of intransitive triplets a trial in involved in throughout the task and this is true for all groups in both tasks.

summary(mem1)

Linear mixed model fit by REML ['lmerMod']
Formula: PropInIntransTriplet ~ (1 | f.id) + scale(Trialnumber) + Group *  
    Task
   Data: both.trial.data

REML criterion at convergence: 125147

Scaled residuals: 
   Min     1Q Median     3Q    Max 
 -2.22  -0.23  -0.22  -0.11  39.76 

Random effects:
 Groups   Name        Variance Std.Dev.
 f.id     (Intercept) 0.266    0.516   
 Residual             5.676    2.383   
Number of obs: 27302, groups:  f.id, 91

Fixed effects:
                        Estimate Std. Error t value
(Intercept)              0.52975    0.09936    5.33
scale(Trialnumber)      -0.00241    0.01442   -0.17
GroupETL                 0.17915    0.14053    1.27
GroupMTL                 0.01709    0.13946    0.12
TaskPerceptual          -0.00167    0.06236   -0.03
GroupETL:TaskPerceptual -0.02376    0.09162   -0.26
GroupMTL:TaskPerceptual  0.00819    0.07879    0.10

Correlation of Fixed Effects:
            (Intr) scl(T) GrpETL GrpMTL TskPrc GETL:T
scl(Trlnmb)  0.000                                   
GroupETL    -0.707  0.000                            
GroupMTL    -0.712  0.000  0.504                     
TaskPercptl -0.161  0.001  0.114  0.115              
GrpETL:TskP  0.110  0.000 -0.156 -0.078 -0.681       
GrpMTL:TskP  0.128  0.000 -0.090 -0.181 -0.792  0.539

In her comments Elke suggested that based on the episodic memory hypothesis the episodic memory of the choices in the first two trials would only help/influence the decision for the third trial in a triplet. To address this I analyzed the episodic memory hypothesis at the triplet level. Previously we looked at the number of intransitive triplets each trial was involved in depending on when it appeared during the task. This might seem confusing because depending on when they appeared during the task each trial is not equally likely to have lead to an intransitive triplet. That is, as Elke implied not every trial was the third trial in a triplet and looking at all tasks without discriminating the order in which they appeared for a triplet masks this.

The analyses here should take care of this and indeed lead to the same conclusion but might be more intuitive. Here we look at whether a triplet is intransitive depending on when the last trial in that triplet appeared.

both.triplet.data = both.triplet.data %>%
   mutate(CleanIntr = ifelse(Intrans == 1 & Error == 0,1, 0)) %>%
  mutate(Task = factor(Task, levels = c("choice", "numbers"), labels = c("Value-based", "Perceptual")))

both.triplet.data %>%
  ggplot(aes(ThirdTrial, CleanIntr, group=Group, color=Group))+
  geom_smooth(method="glm", method.args = list(family="binomial"))+
  theme_bw()+
  facet_grid(~Task)+
  xlab("Last seen trial in triplet")+
  ylab("Intransitive triplet")

Running model for choice task only because looking at both takes a while. The DV in these regressions is whether the triplet is intransitive and the independent variables are random effects for subjects and fixed effects for patient group and the trial number of the last seen trial in that triplet.

For the choice task there is only a fixed effect of the MTL group suggesting that this group has more intransitive triplets than the control group. The effect of the trial number of the last seen trial in a triplet on whether or not it is intransitive is not significant.

mem_val = glmer(CleanIntr ~ (1|f.id)+Group*scale(ThirdTrial),family=binomial, both.triplet.data[both.triplet.data$Task == "Value-based",])
summary(mem_val)

Generalized linear mixed model fit by maximum likelihood (Laplace
  Approximation) [glmerMod]
 Family: binomial  ( logit )
Formula: CleanIntr ~ (1 | f.id) + Group * scale(ThirdTrial)
   Data: both.triplet.data[both.triplet.data$Task == "Value-based", ]

     AIC      BIC   logLik deviance df.resid 
   32336    32403   -16161    32322   103733 

Scaled residuals: 
   Min     1Q Median     3Q    Max 
-0.553 -0.218 -0.178 -0.136 15.281 

Random effects:
 Groups Name        Variance Std.Dev.
 f.id   (Intercept) 0.466    0.682   
Number of obs: 103740, groups:  f.id, 91

Fixed effects:
                           Estimate Std. Error z value Pr(>|z|)    
(Intercept)                 -3.7174     0.1302  -28.56  < 2e-16 ***
GroupETL                     0.0972     0.1840    0.53  0.59725    
GroupMTL                     0.6825     0.1811    3.77  0.00016 ***
scale(ThirdTrial)            0.0266     0.0335    0.79  0.42713    
GroupETL:scale(ThirdTrial)  -0.0442     0.0452   -0.98  0.32817    
GroupMTL:scale(ThirdTrial)  -0.0319     0.0409   -0.78  0.43417    
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Correlation of Fixed Effects:
            (Intr) GrpETL GrpMTL sc(TT) GETL:(
GroupETL    -0.706                            
GroupMTL    -0.719  0.509                     
scl(ThrdTr) -0.006  0.004  0.004              
GrpETL:(TT)  0.004 -0.001 -0.003 -0.742       
GrpMTL:(TT)  0.005 -0.003 -0.003 -0.821  0.609

The episodic memory hypothesis might make a more specific hypothesis that depends not just on when each trial appears but how far apart each trial in a triplet appears from each other. Specifically, it could be the case that the further apart the trials of a triplet are the more likely they are to form an intransitive triplet because the episodic memory of the choices fade away and can’t help the decision in the thirdpair. For this hypothesis to explain the group difference in intransitivities for the choice task this trend would have to be particularly pronounced for the MTL group.

both.triplet.data %>%
    mutate(TrialNumVar = ((TrialA-((TrialA+TrialB+TrialC)/3))^2+(TrialB-((TrialA+TrialB+TrialC)/3))^2+(TrialC-((TrialA+TrialB+TrialC)/3))^2)/2) %>%
  ggplot(aes(TrialNumVar, CleanIntr, group=Group, color=Group))+
  geom_smooth(method="glm", method.args = list(family="binomial"))+
  theme_bw()+
  facet_grid(~Task)+
  xlab("Variance of the trial numbers involved in the triplet")+
  ylab("Intransitive triplet")

Again, however, this hypothesis does not seem to explain the group difference in intransitivies in the choice task: It is the case that the further apart from each other trials in a triplet appear the more likely is that triplet intransitive but this trend is the same for all groups.

both.triplet.data=both.triplet.data %>%
    mutate(TrialNumVar = ((TrialA-((TrialA+TrialB+TrialC)/3))^2+(TrialB-((TrialA+TrialB+TrialC)/3))^2+(TrialC-((TrialA+TrialB+TrialC)/3))^2)/2)

mem_val_var <- glmer(CleanIntr ~ (1|f.id)+Group*scale(TrialNumVar),family=binomial, both.triplet.data[both.triplet.data$Task == "Value-based",])
summary(mem_val_var)

Generalized linear mixed model fit by maximum likelihood (Laplace
  Approximation) [glmerMod]
 Family: binomial  ( logit )
Formula: CleanIntr ~ (1 | f.id) + Group * scale(TrialNumVar)
   Data: both.triplet.data[both.triplet.data$Task == "Value-based", ]

     AIC      BIC   logLik deviance df.resid 
   32312    32379   -16149    32298   103733 

Scaled residuals: 
   Min     1Q Median     3Q    Max 
-0.608 -0.217 -0.177 -0.134 15.406 

Random effects:
 Groups Name        Variance Std.Dev.
 f.id   (Intercept) 0.466    0.683   
Number of obs: 103740, groups:  f.id, 91

Fixed effects:
                            Estimate Std. Error z value Pr(>|z|)    
(Intercept)                  -3.7228     0.1302  -28.59  < 2e-16 ***
GroupETL                      0.0997     0.1841    0.54  0.58825    
GroupMTL                      0.6861     0.1810    3.79  0.00015 ***
scale(TrialNumVar)            0.1088     0.0320    3.40  0.00068 ***
GroupETL:scale(TrialNumVar)  -0.0271     0.0437   -0.62  0.53493    
GroupMTL:scale(TrialNumVar)  -0.0489     0.0393   -1.25  0.21308    
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Correlation of Fixed Effects:
            (Intr) GrpETL GrpMTL s(TNV) GETL:(
GroupETL    -0.706                            
GroupMTL    -0.719  0.509                     
scl(TrlNmV) -0.026  0.018  0.019              
GrETL:(TNV)  0.019 -0.022 -0.014 -0.732       
GrMTL:(TNV)  0.021 -0.015 -0.019 -0.815  0.597

Response Times

Intransitivities and response times (speed-accuracy tradeoff)

The y-axis is the proportion of intransitive triplets a trial is involved in. The higher lines for the numbers task might be confusing at first look given the few number of intransitivities in this task. It is because there are so few intransitivities that if a trial is involved in even a single intransitivity it translates to a larger proportions. The lower number if intransitivities for the numbers task is captured by the large error bars around the task (i.e. for each response time there are fewer non-zero data points involved in intransitivities compared to the choice task).

both.trial.data %>%
  ggplot(aes(RT, PropInIntransTriplet, group=Group, color=Group))+
  # You don't want to plot the raw measure because then a single subject who has done a lot of intransitivities and is in general faster can make it seem like faster trials are associated with more intransitivities. You want the proportional involvement given the total number of intransitivities a subject makes over response time
  # ggplot(aes(RT, IntransTripleCounted, group=Group, color=Group))+
  theme_bw()+
  geom_smooth(method = "loess", span = 5)+
  facet_grid(~Task)+
  ylim(c(0,NA))+
  xlab("Response Time")+
  ylab("Percent of intransitive Triplets involved")

Is there a statistical difference in the effect of response times on intransitivities between the two tasks?

We’re not necessarily interested in whether the exact slopes are different in the two tasks. We can see they might be visually. The main question of interest is two-fold:

• Is there a speed-accuracy tradeoff? Are shorter RTs associated with more intransitivities (in a linear model: does the linear effect of RT have a negative slope)

• Does this differ between groups and across the two different tasks.

First we do a model comparison and check whether the three-way interactions of RTs with groups and tasks as well as a quadratic effect of RTs lead to significant imporevements. Both do.

rt1 = lmer(PropInIntransTriplet ~ (1|f.id)+scale(RT)+Group*Task, both.trial.data)

#adding interaction of linear effect of RT between tasks and groups
rt2 = lmer(PropInIntransTriplet ~ (1|f.id)+scale(RT)*Group*Task, both.trial.data)

#adding the quadratic effect
rt3 = lmer(PropInIntransTriplet ~ (1|f.id)+scale(poly(RT,2))*Task*Group, both.trial.data)

anova(rt1, rt2, rt3)

So what does this model say? There is a significant linear effect of time where longer RTs are associated with more intransitivies (the opposite of a speed-accuracy tradeoff). There is a main effect of task where the numbers task has larger proportions of intransitivities but this as mentioned above is due to the lower number intransitivities per subject. The linear and quadratic effects of time are different for the ETL group compared to the control group. The increase in the proportion of intransitivities for the MTL group is less than the increase in this proportion for the control group (this isn’t interesting because the control group makes barely any intransitivities in the perceptual task). The change in the linear effect of time in the numbers task compared to the choice task is less this change for the control group (also not interesting because there is a quadratic effect of time for the ETL group and this change is not significant).

As qualitatively suggested in the above figure there is no negative slopes indicating a speed accuracy tradeoff for either task. Also crucially the interactions for the significant effects of time for each group are not significant: The interaction between the linear effect time and task for the control group, the threeway interaction between the linear effect of time and task for the MTL group and the threeway interaction between the quadratic effect of time and task for the ETL group are not significant. That is, the effect of time on intransitivities is the same for both tasks.

summary(rt3)

Linear mixed model fit by REML ['lmerMod']
Formula: PropInIntransTriplet ~ (1 | f.id) + scale(poly(RT, 2)) * Task *  
    Group
   Data: both.trial.data

REML criterion at convergence: 124868

Scaled residuals: 
   Min     1Q Median     3Q    Max 
 -2.52  -0.25  -0.17  -0.06  39.92 

Random effects:
 Groups   Name        Variance Std.Dev.
 f.id     (Intercept) 0.282    0.531   
 Residual             5.611    2.369   
Number of obs: 27302, groups:  f.id, 91

Fixed effects:
                                            Estimate Std. Error t value
(Intercept)                                   0.4000     0.1038    3.85
scale(poly(RT, 2))1                           0.2971     0.0365    8.13
scale(poly(RT, 2))2                          -0.0459     0.0345   -1.33
TaskPerceptual                                0.6372     0.1163    5.48
GroupETL                                      0.1094     0.1467    0.75
GroupMTL                                     -0.0145     0.1473   -0.10
scale(poly(RT, 2))1:TaskPerceptual            0.3847     0.2090    1.84
scale(poly(RT, 2))2:TaskPerceptual           -0.1390     0.1427   -0.97
scale(poly(RT, 2))1:GroupETL                  0.1083     0.0507    2.13
scale(poly(RT, 2))2:GroupETL                 -0.1066     0.0462   -2.31
scale(poly(RT, 2))1:GroupMTL                 -0.0974     0.0503   -1.93
scale(poly(RT, 2))2:GroupMTL                 -0.0239     0.0441   -0.54
TaskPerceptual:GroupETL                      -0.2695     0.1418   -1.90
TaskPerceptual:GroupMTL                      -0.2740     0.1332   -2.06
scale(poly(RT, 2))1:TaskPerceptual:GroupETL  -0.5733     0.2263   -2.53
scale(poly(RT, 2))2:TaskPerceptual:GroupETL   0.0283     0.1603    0.18
scale(poly(RT, 2))1:TaskPerceptual:GroupMTL  -0.1785     0.2228   -0.80
scale(poly(RT, 2))2:TaskPerceptual:GroupMTL   0.1121     0.1548    0.72

We can also analyze the speed accuracy tradeoff on triplet level by checking whether the response times to any of the three trials in each triplet is associated with whether that triplet is intransitive or not. The conclusions don’t change but this analysis is more fine-grained and the dependent measure might be more interpretable.

both.triplet.data = both.triplet.data %>%
  left_join(both.trial.data[,c("Task", "f.id", "Trialnumber", "RT")], by = c("Task", "f.id", "TrialA" = "Trialnumber")) %>%
  rename(RT_A = RT) %>%
    left_join(both.trial.data[,c("Task", "f.id", "Trialnumber", "RT")], by = c("Task", "f.id", "TrialB" = "Trialnumber")) %>%
  rename(RT_B = RT) %>%
    left_join(both.trial.data[,c("Task", "f.id", "Trialnumber", "RT")], by = c("Task", "f.id", "TrialC" = "Trialnumber")) %>%
  rename(RT_C = RT) %>%
  drop_na()

both.triplet.data %>%
  select(RT_A, RT_B, RT_C, Task, CleanIntr, Group, f.id) %>%
  gather(key, value, -Task, -CleanIntr, -f.id, -Group) %>%
  ggplot(aes(value, CleanIntr, group=Group, color=Group))+
  geom_smooth(method = "glm", method.args = list(family = "binomial"))+
  # geom_smooth(method = "lm")+
  theme_bw()+
  facet_grid(key~Task)+
  xlab("Response time")+
  ylab("Intransitive triplet")

In the linear regression I chose to look only at the effect of the response time in the third trial because:
- that’s the trial that would determine whether a trial is intransitive or not - the graph suggests that the relationship between response times and intransitive triplets is the same for the three trials - the response times for each trial in a triplet are higly correlated with each othe - Note: the full interactive model with all three response times does lead to significantly improved model but becomes too uninterpretable

#Note these should all be glmers; ran these because it was faster.
sa1 = lmer(CleanIntr ~ (1|f.id)+Task*Group+scale(RT_C), both.triplet.data)
sa2 = lmer(CleanIntr ~ (1|f.id)+Task*Group*scale(RT_C), both.triplet.data)
sa3 = lmer(CleanIntr ~ (1|f.id)+Task*Group+scale(RT_A)+scale(RT_B)+scale(RT_C), both.triplet.data)
sa4 = lmer(CleanIntr ~ (1|f.id)+Task*Group*scale(RT_C)+scale(RT_B)+scale(RT_A), both.triplet.data)
sa5 = lmer(CleanIntr ~ (1|f.id)+Task*Group+scale(RT_C)*scale(RT_B)*scale(RT_A), both.triplet.data)
sa6 = lmer(CleanIntr ~ (1|f.id)+Task*Group*scale(RT_C)*scale(RT_B)*scale(RT_A), both.triplet.data)

#need the interaction of rt_c
anova(sa1, sa2)
#additive with 3 rt's is better than only the third rt
anova(sa1, sa3)
#interactive with 3'rts is better than additive with 3 rt's
anova(sa3, sa4)
#for interactive do you need the main of the other two  - yes
anova(sa2, sa4)
#when using all three rt's is making them interactive better? - yes
anova(sa3, sa5)
#is the full interactive model with 3 rt's better than adding the other two rt's additively?-yes
anova(sa4, sa6)
#is the full interactive model with 3 rt's better than adding the rt interactions additively?-yes
anova(sa5, sa6)

Looking at the tasks in separate models because the interactive mode including the task takes a while to run: The model confirms that the slope is positive for all groups, suggesting that there is no speed accuracy tradeoff. In fact, we observe the opposite pattern where triplets are more likely to be intransitive if the response time for the last seen trial in that triplet is longer.

sa_val = lme4::glmer(CleanIntr ~ (1|f.id)+Group*scale(RT_C), family = binomial, both.triplet.data[both.triplet.data$Task == 'Value-based',])
summary(sa_val)

Generalized linear mixed model fit by maximum likelihood (Laplace
  Approximation) [glmerMod]
 Family: binomial  ( logit )
Formula: CleanIntr ~ (1 | f.id) + Group * scale(RT_C)
   Data: both.triplet.data[both.triplet.data$Task == "Value-based", ]

     AIC      BIC   logLik deviance df.resid 
   31766    31833   -15876    31752   100470 

Scaled residuals: 
   Min     1Q Median     3Q    Max 
-0.809 -0.222 -0.175 -0.132 15.669 

Random effects:
 Groups Name        Variance Std.Dev.
 f.id   (Intercept) 0.494    0.703   
Number of obs: 100477, groups:  f.id, 91

Fixed effects:
                     Estimate Std. Error z value Pr(>|z|)    
(Intercept)           -3.7406     0.1340  -27.92  < 2e-16 ***
GroupETL               0.1337     0.1893    0.71  0.47987    
GroupMTL               0.7099     0.1865    3.81  0.00014 ***
scale(RT_C)            0.4479     0.0333   13.45  < 2e-16 ***
GroupETL:scale(RT_C)  -0.1470     0.0459   -3.20  0.00136 ** 
GroupMTL:scale(RT_C)  -0.2213     0.0417   -5.31  1.1e-07 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Correlation of Fixed Effects:
            (Intr) GrpETL GrpMTL s(RT_C GETL:(
GroupETL    -0.706                            
GroupMTL    -0.718  0.508                     
scale(RT_C) -0.064  0.043  0.045              
GETL:(RT_C)  0.045 -0.050 -0.033 -0.724       
GMTL:(RT_C)  0.050 -0.035 -0.067 -0.798  0.579

Response times for groups

ggplot(data = both.trial.data, aes(x = Trialnumber, y = RT, group = Group, linetype = Group))+
  geom_smooth(method = "lm", color = "black")+
  theme_bw()+
  facet_grid(~Task)+
  xlab("Trial number")+
  ylab("Response time")

both.trial.data %>%
  group_by(Group, Task) %>%
  summarise(mean_RT = mean(RT),
            sd_RT = sd(RT),
            median_RT = median(RT))

# Full model for above graph. Not run because too long.
rt.lmer1 <- lmerTest::lmer(scale(RT) ~ Group*Task*scale(Trialnumber) + (1|f.id), data = RT.df[RT.df$RT>0,])

summary(rt.lmer1)

Group differences response times (reported)

summary(lmerTest::lmer(scale(RT) ~ Group + (1|f.id), data = choice2.trial.data[choice2.trial.data$RT>0,]))

Linear mixed model fit by REML t-tests use Satterthwaite approximations
  to degrees of freedom [lmerMod]
Formula: scale(RT) ~ Group + (1 | f.id)
   Data: choice2.trial.data[choice2.trial.data$RT > 0, ]

REML criterion at convergence: 42472

Scaled residuals: 
   Min     1Q Median     3Q    Max 
-3.447 -0.614 -0.221  0.374 12.102 

Random effects:
 Groups   Name        Variance Std.Dev.
 f.id     (Intercept) 0.310    0.557   
 Residual             0.687    0.829   
Number of obs: 17084, groups:  f.id, 91

Fixed effects:
            Estimate Std. Error       df t value Pr(>|t|)  
(Intercept) -0.09598    0.10228 87.90000   -0.94    0.351  
GroupETL     0.00416    0.14465 87.90000    0.03    0.977  
GroupMTL     0.30282    0.14349 88.00000    2.11    0.038 *
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Correlation of Fixed Effects:
         (Intr) GrpETL
GroupETL -0.707       
GroupMTL -0.713  0.504

Idiosyncratic stimulus effects

id.df = choice2.trial.data %>%
  select(IntransTripleCounted, Image_right, Image_left, f.id) %>%
  gather(key, value, -IntransTripleCounted, -f.id) %>%
  group_by(f.id,value) %>%
  summarise(sum_Intrans = sum(IntransTripleCounted))

summary(aov(sum_Intrans ~ factor(value) + Error(f.id/value), data = id.df))

Error: f.id
          Df Sum Sq Mean Sq F value Pr(>F)
Residuals 90 251369    2793               

Error: f.id:value
              Df Sum Sq Mean Sq F value Pr(>F)
factor(value)  1      1     0.5       0   0.96
Residuals     90  15112   167.9               

Error: Within
                Df Sum Sq Mean Sq F value Pr(>F)   
factor(value)   18   3572   198.5    2.12 0.0039 **
Residuals     1620 151328    93.4                  
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

#For exact values
# str(summary(aov(sum_Intrans ~ factor(value) + Error(f.id/value), data = id.df)))

Appendix

Correlations between response times for each trial in a triplet

both.triplet.data %>%
  select(Task, RT_A, RT_B, RT_C, Group, f.id) %>%
  gather(key, value, -Task, -RT_C, -Group, -f.id) %>%
  ggplot(aes(RT_C, value, group=Group, color=Group))+
  geom_smooth(method="lm")+
  theme_bw()+
  facet_grid(key ~ Task)+
  xlab("Response time for the third trial in the triplet")+
  ylab("Response times for the first \n and second trials in triplet")