code_examples.R

df = read.csv('data/study_1_data_ready.csv')
df$category_level = df$morph_level_p
df$choice = as.numeric(df$choice=="Black")

# This function estimate the PSE
find_PSE <- function(data) {
  # Get predictor and choice columns as they are named in the specific data
  predictor <- names(data)[1]
  choice <- names(data)[2]
  # Fit the logistic model
  formula <- as.formula(paste(choice, "~", predictor))
  model <- glm(formula, data = data, family = binomial)
  # Extract the model coefficients
  intercept_coef <-coef(model)[1]
  slope_coef <- coef(model)[2]
  #Find PSE (where choice probability is equal to 50%)
  PSE <- (log(0.5/(1-0.5))-intercept_coef)/slope_coef #this finds x when y=0.5
  # Return a data frame with the PSE
  data.frame(PSE = PSE)
}
#Apply the find_PSE function for each subject using plyr (change variable names
#according to your data, note that "choice" should be coded as 0 and 1)
library(plyr)
PSE_df <- ddply(df[, c("category_level", "choice","subject_id")],
                            .(subject_id), find_PSE)
#merge the resulted PSE data frame with the original data frame
df <- merge(df,PSE_df,by=c("subject_id"),all = FALSE)

#Define P0 as the middle point of the category continuum
P0 <- min(df$category_level)+(max(df$category_level)-min(df$category_level))/2
#Create a variable of adjusted category levels
df$adjusted_category_level <- df$category_level-(df$PSE-P0)
#Create a variable of ambiguity level
df$ambiguity <- abs(df$adjusted_category_level-P0)

library(lmerTest)
#Polynomial model of RT as a function of (adjusted) category level
model_poly = glmer(rt ~ poly(adjusted_category_level,2) + (1|subject_id),
                  data = df, family = inverse.gaussian(link = "log"))
summary(model_poly)

# (Generalized) linear mixed model of RT as a function of ambiguity
#Create a variable of ambiguity level
df$ambiguity <- abs(df$adjusted_category_level-P0)
model_ambiguity = glmer(rt ~ ambiguity + (1|subject_id),
                       data = df, family = inverse.gaussian(link = "log"))
#Show model results
summary(model_ambiguity)