<jupyter_start><jupyter_text>Experiment 2Generates all figures, tables, and statistical comparisons for the __Results / Experiment 2__ section<jupyter_code>rm(list=ls())
library(dplyr)
library(gtools)
library(ggplot2)
library(binom)
library(BayesFactor)
library(lme4)
library(lmerTest)
library(Hmisc)<jupyter_output><empty_output><jupyter_text>1. Combining files from individual observers into a single table.While loading, also marking out individual trials and removing duplicate key presses<jupyter_code>all.csv.files <- list.files(path= 'Experiment 2', pattern= '*.csv')
results <- data.frame()

observers.with.extremely.strong.bias <- c('BPM1997w')

for(current.filename in all.csv.files){
    # excluding observer with extreme bias, see spatially specific bias analysis for details
    if (substr(current.filename, 1, nchar(current.filename)-4) %in% observers.with.extremely.strong.bias){
        cat(sprintf('Excluding %s due to extreme perceptual bias\n', substr(current.filename, 1, nchar(current.filename)-4)))      
        next;
    }
    
    current.observer <- read.csv2(normalizePath(paste('Experiment 2', current.filename, sep='/')))
    
    # correcting any rounding errors on the RelativeLocation
    current.observer$RelativeLocation[which(current.observer$RelativeLocation< -1)]<- -1
    current.observer$RelativeLocation[current.observer$RelativeLocation> 1]<-  1
    
   
    # adding trial index (time between two extremeties)
    current.observer$Trial <- NA
    for(current.block in unique(current.observer$Block)){
        iTrialStartEnd <- which((current.observer$Block==current.block) &  
                                 current.observer$Event %in% c('Block start', 'Trajectory limit', 'Block end'))
        
        for(iCurrentTrial in 1:(length(iTrialStartEnd)-1)){
            current.observer$Trial[iTrialStartEnd[iCurrentTrial]:iTrialStartEnd[iCurrentTrial+1]]<- iCurrentTrial
        }
    }
    

    # removing duplicate percepts, i.e. a participant reporting same percept again after the blink
    i.duplicates <- c()
    for(current.block in unique(current.observer$Block)){
        current.percept <- 100 # percept code that is definitely different from anything being reported
        for(i.row in which(current.observer$Block==current.block)){
            if (current.observer$Event[i.row]=='Percept'){
                if (current.percept == current.observer$Value[i.row]){
                    i.duplicates <- c(i.duplicates, i.row)
                }
                current.percept <- current.observer$Value[i.row]
            }
        }
    }
    if (length(i.duplicates)>0){
        current.observer <- current.observer[-i.duplicates, ]
    }        
    
    # binding tables
    results <- rbind(results, current.observer)
}

# being tidy
rm(list=setdiff(ls(), "results"))


# replacing "b" and "y" with "background" and "yellow"
levels(results$Color)[levels(results$Color)=='b'] <- 'background'
levels(results$Color)[levels(results$Color)=='y'] <- 'yellow'

# making yellow to be the fist, so that legend matches color order
results$Color <- relevel(results$Color, "yellow")<jupyter_output>Excluding BPM1997w due to extreme perceptual bias<jupyter_text>2. Computing persistence, as a function of the aperture area (occluder) and visibility (color)First, we compute number of switches reported when object was __past__ the midpoint, _i.e._ it re-emerged and should either persist (few trials with switches) or not (many trials with switches) and convert it into the Psurvival for each combination of factors.<jupyter_code># computing whether an event occured before or after the mid-point
results$normalized.location <- results$RelativeLocation*results$Direction
results$after.midpoint <- results$normalized.location>0

# counting post midpoint percept onsets (e.i., switches) for each trial
percepts.per.trial <- results %>%
    dplyr::group_by(Observer, Trial, Occluder, Color) %>%
    dplyr::summarize(percepts.after.midpoint= sum(Event=='Percept' & normalized.location>0.1 & normalized.location<0.4))
percepts.per.trial <- data.frame(percepts.per.trial)

# computing number of trials with at least on post-mid-point switch and total number of trials
percepts.per.observer <- percepts.per.trial %>%
                          dplyr::group_by(Observer, Occluder, Color) %>%
                          dplyr::summarise(
                            trials.with.at.least.one.percept= sum(percepts.after.midpoint>0),
                            trials.with.one.percept= sum(percepts.after.midpoint==1),
                            trial.count= n(),
                            P.switch= trials.with.at.least.one.percept/trial.count)

# logit conversion (for ANOVA)
tiny.correction <- 1e-10
percepts.per.observer$logit.P.switch <- logit(percepts.per.observer$P.switch*(1-2*tiny.correction)+tiny.correction)
percepts.per.observer <- data.frame(percepts.per.observer)<jupyter_output><empty_output><jupyter_text>Computing group averages for each condition, plus a 95% binomial confidence interval from a beta posterior (see `binom.bayes` for details).<jupyter_code>switches.per.condition <- percepts.per.observer %>%
    dplyr::group_by(Occluder, Color) %>%
    dplyr::summarise(total.trial.count= sum(trial.count),
                     P.mean= mean(P.switch),
                     P.ci.lower= binom.bayes(round(total.trial.count*mean(P.switch)), total.trial.count)$lower,
                     P.ci.upper= binom.bayes(round(total.trial.count*mean(P.switch)), total.trial.count)$upper)<jupyter_output><empty_output><jupyter_text>Plotting group averages__Figure 3__ in the manuscript.<jupyter_code>switches.per.condition$i.Occluder <- as.numeric(as.factor(switches.per.condition$Occluder))
survival.plot <- ggplot(switches.per.condition, aes(x= i.Occluder, y= P.mean, ymin=P.ci.lower, ymax= P.ci.upper, 
                                                    color= Color))+
    geom_line()+
    geom_errorbar()+
    geom_point(size= 5)+ 
    xlab("Area occupied by apertures [%]") + ylab(expression(P[switch]))+
    scale_x_continuous(breaks= unique(switches.per.condition$i.Occluder), labels= unique(switches.per.condition$Occluder))
print(survival.plot)
ggsave(file="Generated figures/Experiment 2 - group averages.pdf", plot= survival.plot, width= 12, height= 8, units = 'cm')<jupyter_output><empty_output><jupyter_text>Using `anovaBF` from `BayesFactor` to identify the best model for the data.Column __Bayes Factor__ of the __Table 2__ in the manuscript.<jupyter_code>percepts.per.observer$Occluder <- as.factor(percepts.per.observer$Occluder)
percepts.per.observer$Color <- as.factor(percepts.per.observer$Color)
percepts.per.observer$Observer <- as.factor(percepts.per.observer$Observer)
switches.aov.bf <-   anovaBF(logit.P.switch ~ Occluder * Color + Observer, 
    data= percepts.per.observer,
    whichRandom= 'Observer')
sort(switches.aov.bf, decreasing = TRUE)<jupyter_output><empty_output><jupyter_text>Same analysis but using linear mixed models.__Table 2__ in the manuscript.<jupyter_code>percepts.per.observer$Occluder <- as.numeric(as.character(percepts.per.observer$Occluder))

lm.null <- lme4::lmer(P.switch ~ 1 + (1|Observer), data=percepts.per.observer, REML= FALSE)
lm.aperture.area <- update(lm.null, .~. + Occluder)
lm.color <- update(lm.aperture.area, .~. + Color)
lm.interaction <- update(lm.color, .~. + Occluder:Color)
anova(lm.null, lm.aperture.area, lm.color, lm.interaction)
# x<- anova(lm.null, lm.aperture.area, lm.color, lm.interaction)

# x$'Pr(>Chisq)'<jupyter_output><empty_output><jupyter_text>Computing effect size for the best model that uses only width as a fixed factor<jupyter_code>r2.corr.mer <- function(m) {
  lmfit <-  lm(model.response(model.frame(m)) ~ fitted(m))
  summary(lmfit)$r.squared
}
cat(sprintf('Effect size for width-only model R2: %.3f\n', r2.corr.mer(lm.aperture.area)))<jupyter_output>Effect size for width-only model R2: 0.397<jupyter_text>3. Using aperture area as the only fixed factor, comparing each aperture area condition to the maximal visibility condition (100% aperture area, occluder==100). This is for the __the visible occluder__. In the manuscript, this is the upper table above the plot in __Figure 3__.<jupyter_code>percepts.per.observer$occluder.as.factor <- factor(percepts.per.observer$Occluder, levels=sort(unique(percepts.per.observer$Occluder), decreasing = TRUE))

aperture.levels <- sort(unique(percepts.per.observer$Occluder), decreasing = TRUE)
aperture.labels <- sprintf('%.2f', aperture.levels)
percepts.per.observer$aperture.as.factor <- factor(percepts.per.observer$Occluder, levels= aperture.levels, labels= aperture.labels)

rcontrast<-function(t, df) {
    return (sqrt(t^2/(t^2 + df)))
}


lm.aperture.as.factor <- lmerTest::lmer(P.switch ~ aperture.as.factor + (1|Observer), 
                                        data=percepts.per.observer[percepts.per.observer$Color=='yellow', ])
summary(lm.aperture.as.factor)
coefficients.only <- data.frame(summary(lm.aperture.as.factor)$coefficients)
colnames(coefficients.only) <- c('Estimate', 'Std.Error', 'df', 't.value', 'p.value')
dplyr::mutate(coefficients.only, R.sqr= rcontrast(t.value, df))<jupyter_output><empty_output><jupyter_text>4. Same analysis, but for __invisible occluder__.In the manuscript, this is the lower table above the plot in __Figure 3__.<jupyter_code>lm.aperture.as.factor <- lmerTest::lmer(P.switch ~ aperture.as.factor + (1|Observer), 
                                        data=percepts.per.observer[percepts.per.observer$Color=='background', ])
summary(lm.aperture.as.factor)
coefficients.only <- data.frame(summary(lm.aperture.as.factor)$coefficients)
colnames(coefficients.only) <- c('Estimate', 'Std.Error', 'df', 't.value', 'p.value')
dplyr::mutate(coefficients.only, R.sqr= rcontrast(t.value, df))<jupyter_output><empty_output><jupyter_text>5. Linear mixed models analysis but only for 5% and 10% aperture area to see whether the visibility of the occluding object matters<jupyter_code>percepts.per.observer$Occluder <- as.numeric(as.character(percepts.per.observer$Occluder))
stability.5.10 <- percepts.per.observer %>% dplyr::filter(Occluder %in% c(5, 10))

lm.null <- lme4::lmer(P.switch ~ 1 + (1|Observer), data=stability.5.10, REML= FALSE)
lm.aperture.area <- update(lm.null, .~. + Occluder)
lm.color <- update(lm.aperture.area, .~. + Color)
anova(lm.null, lm.aperture.area, lm.color)


r2.corr.mer <- function(m) {
  lmfit <-  lm(model.response(model.frame(m)) ~ fitted(m))
  summary(lmfit)$r.squared
}
lm.color <- update(lm.null, .~. + Color)
cat(sprintf('Effect size for 5%% and 10%% for color only model R2: %.3f\n', r2.corr.mer(lm.color)))<jupyter_output><empty_output><jupyter_text>6. Analysis of exogenous trigger effectiveness. This information is reported in the __Methods / Displays__ section in the manuscript.<jupyter_code>trigger.effectiveness <- results %>%
    dplyr::group_by(Observer, Block, Trial) %>%
    dplyr::summarize(
        has.trigger= sum(Event=='Reverse rotation')>0,
        report.before.midpoint= sum((Event=='Percept')>0 & (after.midpoint==FALSE))>0
    ) %>%
    dplyr::filter(has.trigger==TRUE) %>%
    dplyr::group_by(Observer) %>%
    dplyr::summarize(
        triggered.fraction= sum(report.before.midpoint)/n()
    )

cat(sprintf('Percentage of trials with an effective exogenous trigger: %.2f±%.2f [%.2f..%.2f] (mean±sd [min..max])\n', 
        mean(trigger.effectiveness$triggered.fraction),
        sd(trigger.effectiveness$triggered.fraction),
        min(trigger.effectiveness$triggered.fraction),
        max(trigger.effectiveness$triggered.fraction)))

trigger.effectiveness<jupyter_output>Percentage of trials with an effective exogenous trigger: 0.90±0.10 [0.70..0.97] (mean±sd [min..max])
{'fraction_non_alphanumeric': 0.13661480914732635, 'fraction_numerical': 0.010541528521635637, 'mean_word_length': 4.161913119789381, 'pattern_counts': {'":': 0, '<': 109, '<?xml version=': 0, '>': 67, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 5, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '12', 'dir': 'jupyter-structured-clean-dedup', 'id': '30779', 'license': 'CC-BY-4.0', 'n_tokens_mistral': 4170, 'n_tokens_neox': 3779, 'n_words': 776, 'path': 'Experiment 2.ipynb', 'repo_name': 'alexander-pastukhov/persistence-of-moving-SFM', 'size': '11763'}