Main visual precision descriptive analyses

── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
✔ dplyr     1.2.1          ✔ readr     2.2.0     
✔ forcats   1.0.1          ✔ stringr   1.6.0     
✔ ggplot2   4.0.3.9000     ✔ tibble    3.3.1     
✔ lubridate 1.9.5          ✔ tidyr     1.3.2     
✔ purrr     1.2.2          
── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
✖ dplyr::filter() masks stats::filter()
✖ dplyr::lag()    masks stats::lag()
ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
here() starts at /Users/tspri/Documents/visual-precision

Note: The package "relayer" is highly experimental. Use at your own risk.

Loading required package: viridisLite

Loading required package: lme4

Loading required package: Matrix


Attaching package: 'Matrix'


The following objects are masked from 'package:tidyr':

    expand, pack, unpack



Attaching package: 'lmerTest'


The following object is masked from 'package:lme4':

    lmer


The following object is masked from 'package:stats':

    step



Attaching package: 'cowplot'


The following object is masked from 'package:lubridate':

    stamp



Attaching package: 'scales'


The following object is masked from 'package:viridis':

    viridis_pal


The following object is masked from 'package:purrr':

    discard


The following object is masked from 'package:readr':

    col_factor



Attaching package: 'patchwork'


The following object is masked from 'package:cowplot':

    align_plots



Attaching package: 'rlang'


The following objects are masked from 'package:purrr':

    flatten, flatten_chr, flatten_dbl, flatten_int, flatten_lgl,
    flatten_raw, invoke, splice

Warning in readLines(file): incomplete final line found on
'/Users/tspri/Documents/visual-precision/.env'

Load data

PROCESSED_DATA_PATH = here("data", DATA_FOLDER, "processed_data")
looking_time_resampled_clean <- read.csv(file.path(PROCESSED_DATA_PATH, "level-looks_data.csv"))
trial_metadata <- read.csv(here("data","metadata","level-trialtype_data.csv")) |> distinct(Trials.trialID, .keep_all=TRUE)
trial_summary_data <- read.csv(file.path(PROCESSED_DATA_PATH,"level-trials_data.csv"))
all_similarities <- read.csv(here("data", "embeddings", "similarities-all_data.csv"))
openclip_similarities <- read.csv(here("data", "embeddings", "similarities-openclip_data.csv"))
layerwise <- read.csv(here("data/embeddings/similarities-layerwise_data.csv"))
all_sims <- trial_metadata |> left_join(layerwise, by=c("Trials.targetImage"="text1", "Trials.distractorImage"="text2")) |> left_join(all_similarities, by=c("Trials.targetImage"="text1", "Trials.distractorImage"="text2"))

# Filter by section if PROJECT_VERSION specifies a particular sample
if (!is.null(SECTION_FILTER)) {
  section_trials <- trial_metadata |> filter(section == SECTION_FILTER) |> distinct(Trials.trialID)
  trial_summary_data   <- trial_summary_data   |> semi_join(section_trials, by = "Trials.trialID")
  looking_time_resampled_clean <- looking_time_resampled_clean |> semi_join(section_trials, by = "Trials.trialID")
  all_looking_times    <- all_looking_times    |> semi_join(section_trials, by = "Trials.trialID")
}

looking_data_summarized <- trial_summary_data |>
  filter(trial_exclusion == 0 & exclude_participant == 0 & exclude_participant_insufficient_data == 0 & SubjectInfo.subjID != "VVI128") |> dplyr::select(-section) |>
  left_join(trial_metadata) |>
  arrange(AoA_Est_target) 

Joining with `by = join_by(Trials.trialID, Trials.targetImage,
Trials.distractorImage, Trials.imagePair)`

median_age <- looking_data_summarized %>%
  distinct(SubjectInfo.subjID, .keep_all = TRUE) %>%  # Keep only distinct subjects
  summarize(overall_median_age = median(SubjectInfo.testAge, na.rm = TRUE))  %>%
  pull(overall_median_age) / 30

median_age <- round(median_age, 2)

aoa_summary <- trial_metadata |>
  summarize(median_aoa = median(AoA_Est_target), 
            mean_aoa = mean(AoA_Est_target))
median_aoa <- aoa_summary$median_aoa

# Total counts per section
total_counts <- trial_summary_data %>%
  distinct(SubjectInfo.subjID, section) %>%
  count(section, name = "total")

# Filtered counts per section
filtered_counts <- looking_data_summarized %>%
  distinct(SubjectInfo.subjID, section) %>%
  count(section, name = "filtered") 

# Combine and compute percentage
summary_table <- total_counts %>%
  left_join(filtered_counts, by = "section") %>%
  mutate(
    filtered = coalesce(filtered, 0),
    percentage_filtered = round((filtered / total) * 100, 2)
  )

summary_table

  section total filtered percentage_filtered
1 sample1   124       90               72.58
2 sample2   120       90               75.00

util function to see what the order of trials were for individual participants

order_of_trials <- function(data, subjID) {
  curr_order <- data |>
  filter(SubjectInfo.subjID == subjID) |>
  distinct(Trials.trialID, .keep_all=TRUE) |>
  arrange(Trials.ordinal) |>
  dplyr::select(Trials.ordinal, Trials.trialID, Trials.targetImage, Trials.distractorImage, Trials.trialType, Trials.leftImage, Trials.rightImage, Trials.targetAudio)
  write_csv(curr_order, paste0(subjID,".csv"))
}

#order_of_trials(all_looking_times, "L7Y4Y6")

Overall timecourse plot of proportion target looking

#summarizing within subject for each time point
summarize_subj <- looking_time_resampled_clean %>%
  filter(trial_exclusion == 0 & exclude_participant ==0 & exclude_participant_insufficient_data == 0 & SubjectInfo.subjID != "VVI128") %>%
  filter(!is.na(accuracy_transformed)) %>%
  group_by(SubjectInfo.subjID, time_normalized_corrected, SubjectInfo.testAge) %>%
  summarized_data(., "SubjectInfo.subjID", "accuracy_transformed", c("SubjectInfo.testAge", "time_normalized_corrected")) %>%
  rename(mean_accuracy = mean_value)

Warning: There were 779 warnings in `summarize()`.
The first warning was:
ℹ In argument: `ci = qt(0.975, N - 1) * sd_value/sqrt(N)`.
ℹ In group 1: `SubjectInfo.subjID = "VVI006"`, `SubjectInfo.testAge = 510`,
  `time_normalized_corrected = -3267`.
Caused by warning in `qt()`:
! NaNs produced
ℹ Run `dplyr::last_dplyr_warnings()` to see the 778 remaining warnings.

#summarizing across subjects for each time point
summarize_across_subj <- summarize_subj %>%
  group_by(time_normalized_corrected) %>%
  summarized_data(., "time_normalized_corrected", "mean_accuracy", c("time_normalized_corrected")) |>
  rename(accuracy = mean_value)

Warning: There were 8 warnings in `summarize()`.
The first warning was:
ℹ In argument: `ci = qt(0.975, N - 1) * sd_value/sqrt(N)`.
ℹ In group 303: `time_normalized_corrected = 4667`.
Caused by warning in `qt()`:
! NaNs produced
ℹ Run `dplyr::last_dplyr_warnings()` to see the 7 remaining warnings.

looking_times <- ggplot(summarize_across_subj,aes(time_normalized_corrected,accuracy))+
  xlim(-2000,4000)+
  geom_errorbar(aes(ymin=accuracy-ci,ymax=accuracy+ci),width=0, alpha=0.2)+
  #geom_point(alpha=0.2)+
    geom_smooth(method="gam")+
  geom_vline(xintercept=0,size=1.5)+
  geom_hline(yintercept=0.5,size=1.2,linetype="dashed")+
  geom_vline(xintercept=300,linetype="dotted")+
  ylim(0,1)+
  xlab("Time (normalized to target word onset) in ms")+
  ylab("Proportion Target Looking") 

Warning: Using `size` aesthetic for lines was deprecated in ggplot2 3.4.0.
ℹ Please use `linewidth` instead.

looking_times

`geom_smooth()` using formula = 'y ~ s(x, bs = "cs")'

Warning: Removed 129 rows containing non-finite outside the scale range
(`stat_smooth()`).

ggsave(here("figures",PROJECT_VERSION,"prop_looking_across_time.png"),looking_times,width=9,height=6,bg = "white")

`geom_smooth()` using formula = 'y ~ s(x, bs = "cs")'

Warning: Removed 129 rows containing non-finite outside the scale range
(`stat_smooth()`).

Descriptive plots

Usable trial stats

#Overall baseline-corrected proportion target looking by condition
looking_data_exclusions <- trial_summary_data |>
  mutate(trial_exclusion_reason = as.factor(ifelse(is.na(trial_exclusion_reason), "included", trial_exclusion_reason))) |>
  mutate(trial_exclusion_reason = factor(trial_exclusion_reason, levels=c(setdiff(unique(trial_exclusion_reason), "included"), "included")))

# Compute the number of included trials per subject
included_counts <- looking_data_exclusions |>
  filter(trial_exclusion_reason == "included") |>
  count(SubjectInfo.subjID, name = "included_trials")


excluded_subjects <- setdiff(unique(looking_data_exclusions$SubjectInfo.subjID), 
                             unique(looking_data_summarized$SubjectInfo.subjID))

# Ensure all subjects are included and replace NA with 0
exclusions_by_subject <- looking_data_exclusions |>
  summarize(N = n(), .by = c(SubjectInfo.subjID, trial_exclusion_reason)) |>
  left_join(included_counts, by = "SubjectInfo.subjID") |>
  mutate(included_trials = replace_na(included_trials, 0)) |>  # Replace NA values
  # Ordering subject IDs so that the 'red' colored label to mark excluded subjects is placed on the subjects who are actually excluded
  mutate(SubjectInfo.subjID = factor(SubjectInfo.subjID, 
                                     levels = unique(SubjectInfo.subjID[order(-included_trials)])),
         color = ifelse(SubjectInfo.subjID %in% excluded_subjects, "red", "black")) |>
  arrange(color)

looking_data_exclusions_combined <- looking_data_exclusions |>
  summarize(N = n(),
    .by=c(trial_exclusion_reason))
# Format labels with HTML <span> for colors

exclusions_color_labels <- exclusions_by_subject |> 
  distinct(SubjectInfo.subjID, .keep_all = TRUE) |>
  mutate(colored_labels = paste0("<span style='color:", color, "'>", SubjectInfo.subjID, "</span>")) |>
  arrange(desc(included_trials)) |>
  mutate(colored_labels = factor(colored_labels, levels = colored_labels))

ggplot(exclusions_by_subject, aes(x = SubjectInfo.subjID, y = N, fill = trial_exclusion_reason)) +
  geom_bar(stat = "identity", position = "stack") +  
  labs(x = "Subject ID", y = "Count (N)", title = "Trial Exclusions by Subject", fill = "Trial Exclusion Reason") +
  scale_x_discrete(labels = exclusions_color_labels$colored_labels) +  # Apply colored labels
    theme(
    axis.text.x = ggtext::element_markdown(angle = 45, hjust = 1),  # Enable HTML-styled labels
    strip.text = element_text(size = 5)
  ) +
  scale_fill_brewer(palette = "Set3")

#TODO: we have a participant who participated twice but their first session only had 2 trials. figure out how to deal with this/discard first session trials. Also have to figure out how to deal with participants who are excluded by default

Estimating item and subject-level noise

target_looking_item_subject_level <- summarized_data(looking_data_summarized, "Trials.targetImage", "corrected_target_looking", c("SubjectInfo.subjID", "SubjectInfo.testAge", "AoA_Est_target"))

Warning: There were 2782 warnings in `summarize()`.
The first warning was:
ℹ In argument: `ci = qt(0.975, N - 1) * sd_value/sqrt(N)`.
ℹ In group 2: `Trials.targetImage = "acorn"`, `SubjectInfo.subjID = "VVI008"`,
  `SubjectInfo.testAge = 450`, `AoA_Est_target = 5.95`.
Caused by warning in `qt()`:
! NaNs produced
ℹ Run `dplyr::last_dplyr_warnings()` to see the 2781 remaining warnings.

target_looking_item_level <- summarized_data(target_looking_item_subject_level |> rename(mean_target_looking = mean_value), "Trials.targetImage", "mean_target_looking", "AoA_Est_target")

target_looking_item_age_level <- summarized_data(target_looking_item_subject_level |> rename(mean_target_looking = mean_value) |> add_age_split(), "Trials.targetImage", "mean_target_looking", c("AoA_Est_target", "age_half"))

target_looking_trial_subject_level <- summarized_data(looking_data_summarized, "Trials.trialID", "corrected_target_looking", c("SubjectInfo.subjID", "SubjectInfo.testAge", "AoA_Est_distractor", "AoA_Est_target", "Trials.targetImage", "Trials.distractorImage", "Trials.imagePair"))

Warning: There were 4840 warnings in `summarize()`.
The first warning was:
ℹ In argument: `ci = qt(0.975, N - 1) * sd_value/sqrt(N)`.
ℹ In group 1: `Trials.trialID = "easy-acorn-key"`, `SubjectInfo.subjID =
  "VVI006"`, `SubjectInfo.testAge = 510`, `AoA_Est_distractor = 3.58`,
  `AoA_Est_target = 5.95`, `Trials.targetImage = "acorn"`,
  `Trials.distractorImage = "key"`, `Trials.imagePair = "acorn-key"`.
Caused by warning in `qt()`:
! NaNs produced
ℹ Run `dplyr::last_dplyr_warnings()` to see the 4839 remaining warnings.

target_looking_trial_level <- summarized_data(target_looking_trial_subject_level |> rename(mean_target_looking = mean_value), "Trials.trialID", "mean_target_looking", c("AoA_Est_distractor", "AoA_Est_target", "Trials.targetImage", "Trials.distractorImage", "Trials.imagePair"))

target_looking_subject_level <- summarized_data(looking_data_summarized, "SubjectInfo.subjID", "corrected_target_looking", "SubjectInfo.testAge")

Item-level performance

item_performances <- ggplot(target_looking_item_level, aes(reorder(Trials.targetImage, mean_value), mean_value)) +
   geom_hline(yintercept=0,linetype="dashed")+
  geom_errorbar(aes(ymin=lower_ci,ymax=upper_ci),width=0, alpha=0.6)+
  (geom_point(aes(aoa_size=AoA_Est_target)) |> rename_geom_aes(new_aes = c("size" = "aoa_size")))+
  (geom_jitter(data=target_looking_item_subject_level, aes(x=Trials.targetImage, y=mean_value, color=SubjectInfo.subjID, age_size = SubjectInfo.testAge/30), alpha=0.3, width=0.2) |> rename_geom_aes(new_aes = c("size" = "age_size"))) +
  xlab("Target image")+
  ylab("Proportion of target looking")+
  ggtitle("Mean proportion of target looking across target items") +
  theme(axis.title.x = element_text(face="bold", size=15, vjust=-1),
        axis.text.x  = element_text(size=10,angle=45,vjust=0.5),
        axis.title.y = element_text(face="bold", size=15),
        axis.text.y  = element_text(size=10),
        strip.text.x = element_text(size = 10,face="bold"),
        plot.margin = margin(t = 10, r = 10, b = 30, l = 10)
        ) +
  scale_y_continuous(breaks = seq(-1, 1, by = 0.2)) +
  scale_size_c(aesthetics = "age_size",name = "Age of participant in months", range=c(2,4), guide = guide_legend(order = 2)) +
  scale_size_c( aesthetics = "aoa_size",name = "Est. AoA of target words in years", guide = guide_legend(order = 1)) +
  scale_color_viridis_d(name = "Participant IDs",option="D") +
  guides(
    color = "none"
  ) 

Warning in geom_point(aes(aoa_size = AoA_Est_target)): Ignoring unknown
aesthetics: aoa_size

Warning in geom_jitter(data = target_looking_item_subject_level, aes(x =
Trials.targetImage, : Ignoring unknown aesthetics: age_size

Warning: The `scale_name` argument of `continuous_scale()` is deprecated as of ggplot2
3.5.0.

Warning: The `trans` argument of `continuous_scale()` is deprecated as of ggplot2 3.5.0.
ℹ Please use the `transform` argument instead.

item_performances

ggsave(here("figures",PROJECT_VERSION,"item_performances.png"),item_performances,width=15,height=10,bg = "white")

trial_performances <- ggplot(summarized_data(looking_data_summarized, "Trials.trialID", "corrected_target_looking", c("Trials.trialID", "AoA_Est_target")), aes(reorder(Trials.trialID, mean_value), mean_value)) +
   geom_hline(yintercept=0,linetype="dashed")+
  geom_errorbar(aes(ymin=lower_ci,ymax=upper_ci),width=0, alpha=0.2)+
  (geom_point(aes(aoa_size=AoA_Est_target)) |> rename_geom_aes(new_aes = c("size" = "aoa_size")))+
  (geom_jitter(data=looking_data_summarized,aes(x=Trials.trialID, y=corrected_target_looking, color=SubjectInfo.subjID, age_size = SubjectInfo.testAge/30), alpha=0.3, width=0.2) |> rename_geom_aes(new_aes = c("size" = "age_size"))) +
  xlab("Trial type")+
  ylab("Proportion of target looking")+
  ggtitle("Mean proportion of target looking across trial types") +
  scale_size_c(aesthetics = "age_size",name = "Age of participant in months", range=c(2,4), guide = guide_legend(order = 2)) +
  scale_size_c( aesthetics = "aoa_size",name = "Est. AoA of target words in years", guide = guide_legend(order = 1)) +
  scale_color_viridis_d(name = "Participant IDs",option="D") +
  coord_cartesian(ylim = c(-1, 1)) +
    theme(axis.title.x = element_text(face="bold", size=15, vjust=-1),
        axis.text.x  = element_text(size=8,angle=45,hjust=1),
        axis.title.y = element_text(face="bold", size=15),
        axis.text.y  = element_text(size=15),
        strip.text.x = element_text(size = 10,face="bold"),
        aspect.ratio = 1,
        plot.margin = margin(t = 10, r = 10, b = 30, l = 10)
        ) +
  guides(
    color = "none"
  ) 

Warning in geom_point(aes(aoa_size = AoA_Est_target)): Ignoring unknown
aesthetics: aoa_size

Warning in geom_jitter(data = looking_data_summarized, aes(x = Trials.trialID,
: Ignoring unknown aesthetics: age_size

trial_performances

####Item-level timecourses

# First, let's make sure we have the correct data
summarize_target <- looking_time_resampled_clean %>%
  filter(trial_exclusion == 0 & exclude_participant == 0 & exclude_participant_insufficient_data == 0) %>%
  #dplyr::select(-section) |>
  left_join(trial_metadata) %>%
  filter(!is.na(accuracy_transformed)) %>%
  group_by(Trials.targetImage, time_normalized_corrected) %>%
  #filter(Trials.targetImage %in% c("acorn", "bulldozer", "cheese", "potato", "squirrel", "snail", "turkey", "turtle")) %>%
  summarized_data(., "Trials.targetImage", "accuracy_transformed", c("Trials.targetImage", "AoA_Est_target", "time_normalized_corrected")) %>%
  rename(mean_accuracy = mean_value)

Joining with `by = join_by(Trials.trialID, Trials.targetImage,
Trials.distractorImage, Trials.imagePair)`

Warning: There were 705 warnings in `summarize()`.
The first warning was:
ℹ In argument: `ci = qt(0.975, N - 1) * sd_value/sqrt(N)`.
ℹ In group 1: `Trials.targetImage = "acorn"`, `AoA_Est_target = 5.95`,
  `time_normalized_corrected = -3800`.
Caused by warning in `qt()`:
! NaNs produced
ℹ Run `dplyr::last_dplyr_warnings()` to see the 704 remaining warnings.

# Create a mapping of target images to their AoA values
target_aoa_mapping <- summarize_target %>%
  dplyr::select(Trials.targetImage, AoA_Est_target) %>%
  distinct() %>%
  arrange(AoA_Est_target)  # Sort by AoA

# Convert Trials.targetImage to a factor with levels ordered by AoA (descending)
summarize_target <- summarize_target %>%
  mutate(Trials.targetImage = factor(Trials.targetImage, 
                                     levels = target_aoa_mapping$Trials.targetImage))
aoa_values <- target_aoa_mapping$AoA_Est_target
# Normalize aoa_values to range between 0 and 1
aoa_scaled <- rescale(aoa_values, to = c(0, 1))

# Generate a high-resolution color gradient (e.g., 1000 colors)
num_colors <- 1000
color_palette <- viridis(num_colors, option = "E",alpha=0.8)

# Map scaled values to colors using interpolation
color_values <- color_palette[ceiling((1-aoa_scaled) * (num_colors - 1)) + 1]


# Create named color vector
target_colors <- setNames(color_values, target_aoa_mapping$Trials.targetImage)


# Create the plot with a discrete legend but gradient colors
looking_times_target <- ggplot(summarize_target, 
                              aes(time_normalized_corrected, 
                                  mean_accuracy, 
                                  color = AoA_Est_target, group=Trials.targetImage)) +  # Group by target image
  annotate("rect", xmin = 300, xmax = 3500, ymin = 0, ymax = 1, 
          fill = "gray95", alpha = 0.5) +  
 annotate("rect", xmin = -2000, xmax = 0, ymin = 0, ymax = 1, 
           fill = "gray95", alpha = 0.5) + 

  coord_cartesian(xlim = c(-2000, 3500), ylim = c(0, 1)) +
     geom_line (alpha=0.2, size=0.8) +
   geom_smooth(data=summarize_across_subj,
              aes(x=time_normalized_corrected,
                  y=accuracy,
                  ymin=lower_ci,  # Assuming 'se' is your standard error column
                  ymax=upper_ci,  # Adjust these to match your error metric
                  group=NA),
              stat="identity",  # Use identity instead of a smoothing method
              color="black",
              fill="gray40",
              size=1.5) +
   scale_color_viridis_c(option="B",direction=-1,name="Estimated target\nword AoA") +
  guides(color = guide_colorbar(reverse = TRUE))+
  geom_vline(xintercept = 0, size = 1.5) +
  geom_hline(yintercept = 0.5, size = 1.2, linetype = "dashed") +
  geom_vline(xintercept = 300, linetype = "dotted", size=1.2) +
  xlab("Time (normalized to target word onset) in ms") +
  ylab("Proportion\nTarget Looking") +
    # Add arrows
  annotate("segment", x = -2000, xend = -50, y = 0.85, yend = 0.85, 
           arrow = arrow(length = unit(0.2, "cm"), type="closed", ends="both"), size = 1, color = "black") + 
  annotate("segment", x = 340, xend = 3490, y = 0.85, yend = 0.85, 
           arrow = arrow(length = unit(0.2, "cm"),type="closed", ends="both"), size = 1, color = "black") + 

  # Add text above arrows
  annotate("text", x = -1000, y = 0.93, label = "Baseline window", 
           size = 6, fontface = "bold") +
  annotate("text", x = 2000, y = 0.93, label = "Critical window", 
           size = 6, fontface = "bold") +
  theme_minimal() +
  scale_x_continuous(breaks=seq(-2000,3500,by=1000)) +
  scale_y_continuous(breaks=seq(0,1,by=0.5))+
  theme(
    text = element_text(size=10,face = "bold"),              # All text bold
    # Increase space between x-axis and its title
    axis.title.x = element_text(
      face = "bold", 
      size = 18,
      margin = margin(t = 10, r = 0, b = 0, l = 0)
    ),
    
    # Increase space between y-axis and its title
    axis.title.y = element_text(
      face = "bold", 
      size= 18,
      margin = margin(t=0,r = 5, b = 0, l = 0)
    ),      
    axis.text = element_text(size=14,face = "bold"),         # Axis text bold
    legend.title = element_text(size = 14, face = "bold"),  # Small bold legend title
    legend.text = element_text(size = 14, face = "bold")    # Small bold legend text
  )

# Display the plot
looking_times_target

ggsave(here("figures",PROJECT_VERSION,"timecourse_plot.svg"),looking_times_target, device="pdf",width=8.5,height=3.3,bg = "white")
ggsave(here("figures",PROJECT_VERSION,"timecourse_plot.png"),looking_times_target,width=9.5,height=5.3,bg = "white")

Age plots

Age-split timecourses

summarize_across_age_halves <- summarize_subj |> add_age_split() |>
  group_by(time_normalized_corrected, age_half) |>
  dplyr::summarize(n=n(),
            non_na_n = sum(!is.na(mean_accuracy)),
            mean_value=mean(mean_accuracy,na.rm=TRUE),
            sd_accuracy=sd(mean_accuracy,na.rm=TRUE),
            se_accuracy=sd_accuracy/sqrt(n),
            ci = ifelse(n > 1, qt(0.975, n - 1) * sd_accuracy / sqrt(n), NA)) %>%
  filter(n > 5)

`summarise()` has regrouped the output.
ℹ Summaries were computed grouped by time_normalized_corrected and age_half.
ℹ Output is grouped by time_normalized_corrected.
ℹ Use `summarise(.groups = "drop_last")` to silence this message.
ℹ Use `summarise(.by = c(time_normalized_corrected, age_half))` for
  per-operation grouping (`?dplyr::dplyr_by`) instead.

looking_times_age <- ggplot(summarize_across_age_halves,aes(time_normalized_corrected,mean_value,color=age_half))+
  xlim(-2000,4000)+
  geom_errorbar(aes(ymin=mean_value-ci,ymax=mean_value+ci),width=0, alpha=0.2)+
  #geom_point(alpha=0.2)+
    geom_smooth(method="gam")+
  geom_vline(xintercept=0,size=1.5)+
  geom_hline(yintercept=0.5,size=1.2,linetype="dashed")+
  geom_vline(xintercept=300,linetype="dotted")+
  ylim(0,1)+
  xlab("Time (normalized to target word onset) in ms")+
  ylab("Proportion Target Looking") +
  labs(title="Proportion of looking time across age and time", caption = (paste0("14-24 month olds age split at median=",median_age," months"))) +
  scale_color_brewer(palette = "Set1", name="Age half")

looking_times_age

`geom_smooth()` using formula = 'y ~ s(x, bs = "cs")'

Warning: Removed 157 rows containing non-finite outside the scale range
(`stat_smooth()`).

ggsave(here("figures",PROJECT_VERSION,"prop_looking_across_time_age.png"),looking_times,width=9,height=6,bg = "white")

`geom_smooth()` using formula = 'y ~ s(x, bs = "cs")'

Warning: Removed 129 rows containing non-finite outside the scale range
(`stat_smooth()`).

Age-split performance

# note: using shorter window just because there's clearly no age-related effect with longer window 
participant_age_halves = summarized_data(looking_data_summarized |> left_join(target_looking_item_level) |> rename(mean_looking = mean_value) |>  add_age_split() |> mutate(acc = corrected_target_looking), "age_half", "acc", 
                                         c("age_half", "SubjectInfo.subjID", "SubjectInfo.testAge")) 

Joining with `by = join_by(Trials.targetImage, AoA_Est_target)`

participant_age_halves$age_half <- factor(participant_age_halves$age_half, levels = rev(levels(factor(participant_age_halves$age_half))))

participant_age_halves_summarized = summarized_data(participant_age_halves |> rename(looking_time = mean_value), "age_half", "looking_time", c("age_half"))
accuracy_age <- ggplot(participant_age_halves,
       aes(x = age_half, y = mean_value, color = age_half)) +
    geom_violin(aes(color=age_half), alpha=0.3) +
  geom_jitter(size = 3, alpha = 0.3, width=0.2) +
  geom_hline(yintercept=0,size=1.2,linetype="dashed")+ 
  geom_point(data = participant_age_halves_summarized, size = 3) +
  geom_linerange(data = participant_age_halves_summarized, aes(ymin = mean_value - ci, ymax = mean_value + ci), alpha = 0.9) + 
  scale_color_brewer(palette = "Set1", name="Age half", direction=-1) +  # Using RColorBrewer for colors
  ylab("Baseline-corrected proportion of target looking") +
  xlab("Age half of child") +
  ggpubr::stat_cor(method = "pearson") +
  labs(caption = (paste0("14-24 month olds age split at median=",median_age, " months")), title="Word recognition acuracy by age")

Registered S3 method overwritten by 'car':
  method           from
  na.action.merMod lme4

accuracy_age

ggsave(here("figures",PROJECT_VERSION,"accuracy_age.png"),accuracy_age,width=9,height=6,bg = "white")

Age continuous

participant_age_jittered <- participant_age_halves %>%
  rowwise() %>%
  mutate(x_jittered = SubjectInfo.testAge/30 + runif(n(), -0.25, 0.25))
jitterer <- position_jitter(width = .5,seed=1)
age_continuous <- ggplot(participant_age_jittered, aes(x = x_jittered, y = mean_value)) +
  geom_hline(yintercept = 0, linetype = "dashed") +
  geom_point(size = 4, alpha = 0.4, color="#2171B5") +
  geom_smooth(method="lm",color="#2171B5") +
  geom_linerange(aes(ymin = mean_value - ci, ymax = mean_value + ci), alpha = 0.1) + 
  scale_y_continuous(breaks = seq(-0.4, 0.3, by = 0.1), limits = c(-0.3, 0.3)) +
  ggpubr::stat_cor()+
  labs(x = "Age in months", y = "Baseline-corrected proportion target looking") +
  #ggtitle("Corrected proportion of target looking by age") +
  theme_minimal()
age_continuous

`geom_smooth()` using formula = 'y ~ x'

Warning: Removed 31 rows containing missing values or values outside the scale range
(`geom_segment()`).

participant_age_jittered <- participant_age_halves %>%
  rowwise() %>%
  mutate(x_jittered = SubjectInfo.testAge/30 + runif(n(), -0.25, 0.25))
jitterer <- position_jitter(width = .5,seed=1)
age_continuous_formatted <- ggplot(participant_age_jittered, aes(x = x_jittered, y = mean_value)) +
  geom_hline(yintercept = 0, linetype = "dashed") +
  geom_point(size = 4, alpha = 0.4, color="#2171B5") +
  geom_smooth(method="lm",color="#2171B5") +
  geom_linerange(aes(ymin = mean_value - ci, ymax = mean_value + ci), alpha = 0.1) + 
  scale_y_continuous(breaks = seq(-0.4, 0.3, by = 0.1), limits = c(-0.3, 0.3)) +
  labs(x = "Age in months", y = "Baseline-corrected\nproportion target looking") +
  #ggtitle("Corrected proportion of target looking by age") +
  theme(
    text = element_text(size=10,face = "bold"),              # All text bold
    # Increase space between x-axis and its title
    axis.title.x = element_text(
      face = "bold", 
      size = 18,
      margin = margin(t = 10, r = 0, b = 0, l = 0)
    ),
    
    # Increase space between y-axis and its title
    axis.title.y = element_text(
      face = "bold", 
      size= 18,
      margin = margin(t=0,r = 5, b = 0, l = 0)
    ),      
    axis.text = element_text(size=14,face = "bold"),         # Axis text bold
    legend.title = element_text(size = 14, face = "bold"),  # Small bold legend title
    legend.text = element_text(size = 14, face = "bold")    # Small bold legend text
  )
age_continuous_formatted

`geom_smooth()` using formula = 'y ~ x'

Warning: Removed 31 rows containing missing values or values outside the scale range
(`geom_segment()`).

ggsave(here("figures",PROJECT_VERSION,"age_continuous_formatted.png"),age_continuous_formatted,width=9,height=6,bg = "white")

`geom_smooth()` using formula = 'y ~ x'

Warning: Removed 31 rows containing missing values or values outside the scale range
(`geom_segment()`).

younger group mean

younger_mean <- participant_age_halves |>
  filter(age_half == "younger") |>
  summarize(mean_value = mean(SubjectInfo.testAge)/30) |>
  pull(mean_value)

Order effects

First vs second instances for primary targets

first_instance_target <- looking_data_summarized |>
  group_by(SubjectInfo.subjID, Trials.targetImage) |>
  arrange(Trials.ordinal, .by_group = TRUE) |>
  filter(Trials.trialType %in% c("easy", "hard")) |>
  slice(1) |>
  ungroup() |>
  mutate(order_id = "one")
  
second_instance_target <- looking_data_summarized |>
  group_by(SubjectInfo.subjID, Trials.targetImage) |>
  arrange(Trials.ordinal, .by_group = TRUE) |>
  filter(Trials.trialType %in% c("easy", "hard")) |>
  slice(2) |>
  ungroup() |>
  mutate(order_id = "two")

ordered_conditions <- rbind(first_instance_target, second_instance_target)

order_plot <- half_violins_plot(ordered_conditions, "order_id", "corrected_target_looking", "SubjectInfo.subjID", c("one", "two"), input_xlab = "Target image order")

Warning: There were 4 warnings in `summarize()`.
The first warning was:
ℹ In argument: `ci = qt(0.975, N - 1) * sd_value/sqrt(N)`.
ℹ In group 188: `order_id = "two"`, `SubjectInfo.subjID = "VVI015"`.
Caused by warning in `qt()`:
! NaNs produced
ℹ Run `dplyr::last_dplyr_warnings()` to see the 3 remaining warnings.

# A tibble: 2 × 8
  order_id mean_value sd_value     N     se     ci lower_ci upper_ci
  <chr>         <dbl>    <dbl> <int>  <dbl>  <dbl>    <dbl>    <dbl>
1 one          0.0641    0.142   180 0.0106 0.0209   0.0431   0.0850
2 two          0.0714    0.157   167 0.0122 0.0240   0.0473   0.0954

Warning in geom_path(aes(group = SubjectInfo.subjID), color = "black", fill =
NA, : Ignoring unknown parameters: `fill`

order_plot

Sections

violins

section_plot <- half_violins_plot(looking_data_summarized, "section", "corrected_target_looking", "SubjectInfo.subjID", c("sample1", "sample2"), input_xlab = "Data sample")

# A tibble: 2 × 8
  section mean_value sd_value     N      se     ci lower_ci upper_ci
  <chr>        <dbl>    <dbl> <int>   <dbl>  <dbl>    <dbl>    <dbl>
1 sample1     0.0500   0.0852    90 0.00898 0.0178   0.0322   0.0679
2 sample2     0.0600   0.0964    90 0.0102  0.0202   0.0399   0.0802

Warning in geom_path(aes(group = SubjectInfo.subjID), color = "black", fill =
NA, : Ignoring unknown parameters: `fill`

section_plot

`geom_path()`: Each group consists of only one observation.
ℹ Do you need to adjust the group aesthetic?

timecourse

summarize_across_sections <- summarize_subj |>
  left_join(looking_data_summarized |> distinct(section, SubjectInfo.subjID)) |>
  group_by(time_normalized_corrected, section) |>
  dplyr::summarize(n=n(),
            non_na_n = sum(!is.na(mean_accuracy)),
            mean_value=mean(mean_accuracy,na.rm=TRUE),
            sd_accuracy=sd(mean_accuracy,na.rm=TRUE),
            se_accuracy=sd_accuracy/sqrt(n),
            ci = ifelse(n > 1, qt(0.975, n - 1) * sd_accuracy / sqrt(n), NA)) %>%
  filter(n > 5)

Joining with `by = join_by(SubjectInfo.subjID)`
`summarise()` has regrouped the output.

looking_times_section <- ggplot(summarize_across_sections,aes(time_normalized_corrected,mean_value,color=section))+
  xlim(-2000,4000)+
  geom_errorbar(aes(ymin=mean_value-ci,ymax=mean_value+ci),width=0, alpha=0.2)+
  #geom_point(alpha=0.2)+
    geom_smooth(method="gam")+
  geom_vline(xintercept=0,size=1.5)+
  geom_hline(yintercept=0.5,size=1.2,linetype="dashed")+
  geom_vline(xintercept=300,linetype="dotted")+
  ylim(0,1)+
  xlab("Time (normalized to target word onset) in ms")+
  ylab("Proportion Target Looking") +
  labs(title="Proportion of looking time by section and time") +
  scale_color_brewer(palette = "Set2", name="Section")

looking_times_section

`geom_smooth()` using formula = 'y ~ s(x, bs = "cs")'

Warning: Removed 149 rows containing non-finite outside the scale range
(`stat_smooth()`).

Easy vs hard trial plots

Proportion of target looking for easy vs hard trials

avg_corrected_target_looking_by_condition <- summarized_data(looking_data_summarized,
                                                             "Trials.trialType", "corrected_target_looking", "SubjectInfo.subjID")
high_accuracy_subjects <- avg_corrected_target_looking_by_condition |>
  filter(lower_ci > 0 & N > 5) |>
  distinct(SubjectInfo.subjID)

overall_condition_plot <- half_violins_plot(looking_data_summarized |> group_by(SubjectInfo.subjID, Trials.targetImage) |>
  arrange(Trials.ordinal, .by_group = TRUE) |>
  #slice(2) |>
  ungroup() |> mutate(
  Trials.trialType = ifelse(Trials.trialType == "easy-distractor", "easy", Trials.trialType),
  Trials.trialType = ifelse(Trials.trialType == "hard-distractor", "hard", Trials.trialType)), "Trials.trialType", "corrected_target_looking", "SubjectInfo.subjID", c("easy", "hard"))

# A tibble: 2 × 8
  Trials.trialType mean_value sd_value     N      se     ci lower_ci upper_ci
  <chr>                 <dbl>    <dbl> <int>   <dbl>  <dbl>    <dbl>    <dbl>
1 easy                 0.0695    0.110   180 0.00818 0.0161   0.0533   0.0856
2 hard                 0.0408    0.114   180 0.00851 0.0168   0.0240   0.0576

Warning in geom_path(aes(group = SubjectInfo.subjID), color = "black", fill =
NA, : Ignoring unknown parameters: `fill`

overall_condition_plot

timecourse

summarize_trial_difficulty <- looking_time_resampled_clean %>%
  mutate(
  Trials.trialType = ifelse(Trials.trialType == "easy-distractor", "easy", Trials.trialType),
  Trials.trialType = ifelse(Trials.trialType == "hard-distractor", "hard", Trials.trialType)) |>
  filter(trial_exclusion == 0 & exclude_participant ==0 & exclude_participant_insufficient_data == 0) %>%
  filter(!is.na(accuracy_transformed)) %>%
  group_by(SubjectInfo.subjID, time_normalized_corrected, SubjectInfo.testAge) %>%
  summarized_data(., "Trials.trialID", "accuracy_transformed", c("Trials.trialType", "time_normalized_corrected")) %>%
  rename(mean_accuracy = mean_value) 

Warning: There were 904 warnings in `summarize()`.
The first warning was:
ℹ In argument: `ci = qt(0.975, N - 1) * sd_value/sqrt(N)`.
ℹ In group 1: `Trials.trialID = "easy-acorn-key"`, `Trials.trialType = "easy"`,
  `time_normalized_corrected = -3800`.
Caused by warning in `qt()`:
! NaNs produced
ℹ Run `dplyr::last_dplyr_warnings()` to see the 903 remaining warnings.

summarize_across_trial_type <- summarize_trial_difficulty |>
  group_by(time_normalized_corrected, Trials.trialType) |>
  dplyr::summarize(n=n(),
            non_na_n = sum(!is.na(mean_accuracy)),
            mean_value=mean(mean_accuracy,na.rm=TRUE),
            sd_accuracy=sd(mean_accuracy,na.rm=TRUE),
            se_accuracy=sd_accuracy/sqrt(n),
            ci = ifelse(n > 1, qt(0.975, n - 1) * sd_accuracy / sqrt(n), NA)) 

`summarise()` has regrouped the output.
ℹ Summaries were computed grouped by time_normalized_corrected and
  Trials.trialType.
ℹ Output is grouped by time_normalized_corrected.
ℹ Use `summarise(.groups = "drop_last")` to silence this message.
ℹ Use `summarise(.by = c(time_normalized_corrected, Trials.trialType))` for
  per-operation grouping (`?dplyr::dplyr_by`) instead.

looking_times_trial_type <- ggplot(summarize_across_trial_type,aes(time_normalized_corrected,mean_value,color=Trials.trialType))+
  xlim(-2000,4000)+
  geom_errorbar(aes(ymin=mean_value-ci,ymax=mean_value+ci),width=0, alpha=0.2)+
  #geom_point(alpha=0.2)+
    geom_smooth(method="gam")+
  geom_vline(xintercept=0,size=1.5)+
  geom_hline(yintercept=0.5,size=1.2,linetype="dashed")+
  geom_vline(xintercept=300,linetype="dotted")+
  ylim(0,1)+
  xlab("Time (normalized to target word onset) in ms")+
  ylab("Proportion Target Looking") +
 # labs(title="Proportion of looking time across time, split at mean age-of-acquisition") +
  scale_color_manual(values=c("forestgreen", "#4292C6"), name="Trial type", labels=c( "easy", "hard"))+
   guides(color = guide_legend(reverse = TRUE))

looking_times_trial_type

`geom_smooth()` using formula = 'y ~ s(x, bs = "cs")'

Warning: Removed 251 rows containing non-finite outside the scale range
(`stat_smooth()`).

Target image difficulty

# summarize average accuracy within participant (by word alone)
condition_based_looking <- looking_data_summarized  |>
  filter(!grepl("distractor", Trials.trialType)) |>
  distinct(SubjectInfo.subjID, Trials.trialID, Trials.targetImage, corrected_target_looking, Trials.trialType, AoA_Est_target)

target_looking_by_target_word <- looking_data_summarized  |>
  filter(!grepl("distractor", Trials.trialType)) |>
  group_by(SubjectInfo.subjID, Trials.targetImage) |>
  summarize(target_looking_diff = corrected_target_looking[match("easy", Trials.trialType)] - corrected_target_looking[match("hard", Trials.trialType)],
            baseline_window_looking_diff = mean_target_looking_baseline_window[match("easy", Trials.trialType)] - mean_target_looking_baseline_window[match("hard", Trials.trialType)], 
            critical_window_looking_diff = mean_target_looking_critical_window[match("easy", Trials.trialType)] - mean_target_looking_critical_window[match("hard", Trials.trialType)]) |>
  filter(!is.na(target_looking_diff))

`summarise()` has regrouped the output.
ℹ Summaries were computed grouped by SubjectInfo.subjID and Trials.targetImage.
ℹ Output is grouped by SubjectInfo.subjID.
ℹ Use `summarise(.groups = "drop_last")` to silence this message.
ℹ Use `summarise(.by = c(SubjectInfo.subjID, Trials.targetImage))` for
  per-operation grouping (`?dplyr::dplyr_by`) instead.

#clean names for individual images for plot
overall_target_looking_by_word <- target_looking_by_target_word %>%
  filter(!is.na(target_looking_diff)) %>%
  group_by(Trials.targetImage) %>%
  summarize(N=n(),
            corrected_target_looking=mean(target_looking_diff,na.rm=TRUE),
            ci=qt(0.975, N-1)*sd(target_looking_diff,na.rm=T)/sqrt(N),
            lower_ci=corrected_target_looking-ci,
            upper_ci=corrected_target_looking+ci
          ) 


word_prefs <- ggplot(overall_target_looking_by_word,aes(reorder(Trials.targetImage,corrected_target_looking),corrected_target_looking))+
  geom_hline(yintercept=0,linetype="dashed")+
  geom_errorbar(aes(ymin=lower_ci,ymax=upper_ci),width=0, alpha=0.3)+
  geom_jitter(data=target_looking_by_target_word, aes(x=Trials.targetImage, y=target_looking_diff, color=SubjectInfo.subjID), size = 2.5, alpha = 0.3, width=0.1) +
  geom_point(size = 2.4)+
  xlab("Target Word")+
  ylab("Easy trial prefered target looking")+
  theme(axis.title.x = element_text(face="bold", size=15),
        axis.text.x  = element_text(size=10,angle=90,vjust=0.5,hjust=1),
        axis.title.y = element_text(face="bold", size=15),
        axis.text.y  = element_text(size=10),
        strip.text.x = element_text(size = 10,face="bold")
        ) +
  scale_y_continuous(breaks = seq(-1, 1, by = 0.2)) +
  scale_size_continuous(name = "Number of participants") + 
  scale_color_viridis_d(name = "Participant IDs",option="D") +
  guides(color = "none")
  
word_prefs

Baseline image-pair preferences

saliency_effects <- looking_data_summarized |>
  # Calculating the proportion of time looking at the target word even if it isn't the target word in that particular study
 mutate(original_target_looking_baseline_window = ifelse(grepl("distractor", Trials.trialType), 1 - mean_target_looking_baseline_window,   mean_target_looking_baseline_window),
        original_target_looking_critical_window = ifelse(grepl("distractor", Trials.trialType), 1 - mean_target_looking_critical_window, mean_target_looking_critical_window)) |>
  add_age_split() |>
 summarize(mean_baseline_looking = mean(original_target_looking_baseline_window),
           mean_critical_looking = mean(original_target_looking_critical_window),
           .by = c(Trials.imagePair, age_half))
saliency_effects

    Trials.imagePair age_half mean_baseline_looking mean_critical_looking
1       bathtub-sock    older             0.4045799             0.4437562
2       bathtub-sock  younger             0.5055230             0.4108362
3     bathtub-shower    older             0.4214784             0.4647545
4     bathtub-shower  younger             0.4527236             0.4582721
5   bulldozer-orange  younger             0.4343250             0.3825479
6   bulldozer-orange    older             0.5748909             0.4443989
7       duck-chicken  younger             0.4709926             0.4598276
8       duck-chicken    older             0.4341072             0.4726071
9          truck-car  younger             0.4603808             0.4945374
10         truck-car    older             0.5259696             0.4820767
11    duck-butterfly    older             0.4700572             0.4677344
12    duck-butterfly  younger             0.4626382             0.4664256
13        cloud-tree  younger             0.4225677             0.4191716
14        cloud-tree    older             0.5086454             0.4811971
15         acorn-key  younger             0.4469638             0.4186809
16         acorn-key    older             0.4788331             0.4480183
17      cloud-zipper    older             0.3922510             0.3500712
18      cloud-zipper  younger             0.3537901             0.3767473
19    truck-airplane  younger             0.4644456             0.5409015
20    truck-airplane    older             0.3575887             0.5309715
21        snail-worm    older             0.3996053             0.4479406
22        snail-worm  younger             0.4054677             0.4260101
23         snail-cow    older             0.4532581             0.4402747
24         snail-cow  younger             0.3983266             0.3683577
25       turkey-swan    older             0.4403477             0.4877313
26       turkey-goat  younger             0.4106288             0.4254979
27       turkey-swan  younger             0.4929389             0.4957185
28       turkey-goat    older             0.4253035             0.4760688
29        cheese-mud  younger             0.5015125             0.5914350
30        cheese-mud    older             0.4545625             0.5276830
31      turtle-horse  younger             0.4527153             0.4231412
32      turtle-horse    older             0.4166157             0.3823349
33       turtle-frog    older             0.4551658             0.5246384
34       turtle-frog  younger             0.4927704             0.5100310
35   squirrel-monkey  younger             0.4322144             0.4753925
36   squirrel-monkey    older             0.3600953             0.3757038
37   balloon-dresser    older             0.5991749             0.5560835
38   balloon-dresser  younger             0.6903595             0.6788103
39     cheese-butter    older             0.5058167             0.5619523
40     cheese-butter  younger             0.6411010             0.5965848
41    potato-glasses    older             0.5362731             0.5233058
42    potato-glasses  younger             0.5746260             0.6033294
43      balloon-kite    older             0.5899602             0.5873214
44      balloon-kite  younger             0.6230225             0.6552561
45    squirrel-eagle    older             0.4370153             0.5333088
46    squirrel-eagle  younger             0.5074317             0.5251362
47    camel-elephant  younger             0.4310795             0.4060744
48    camel-elephant    older             0.3300443             0.3892423
49        potato-pot    older             0.5062483             0.5873245
50        potato-pot  younger             0.5675511             0.6441189
51     camel-ladybug    older             0.4104379             0.4337154
52     camel-ladybug  younger             0.4732745             0.5157344
53     lizard-walrus  younger             0.4015119             0.3941967
54     lizard-walrus    older             0.3891069             0.4015273
55      lizard-snake    older             0.3776323             0.3423290
56      lizard-snake  younger             0.3395491             0.3828282
57 bulldozer-tractor    older             0.4571464             0.4083064
58 bulldozer-tractor  younger             0.4087415             0.4220546
59      crow-penguin    older             0.4429395             0.4766056
60      crow-penguin  younger             0.4716791             0.5302641
61     crow-seahorse    older             0.3901576             0.4119352
62     crow-seahorse  younger             0.4873987             0.4789381
63     acorn-coconut  younger             0.5080370             0.5027286
64     acorn-coconut    older             0.5163894             0.5122532

CLIP similarity plots

Comparing proportion of target looking time to CLIP and other model cosine similarities

# CLIP
clip_data_summarized <- summarize_similarity_data(looking_data_summarized)
clip_data_summarized_low_aoa <- summarize_similarity_data(looking_data_summarized |> filter(AoA_Est_target < 4.78))
# |>filter(SubjectInfo.subjID %in% high_accuracy_subjects$SubjectInfo.subjID)))
clip_age_half_summarized <- looking_data_summarized |>
  add_age_split() |>
  summarize_similarity_data(extra_fields = c("age_half", "AoA_Est_target"))

clip_plots <- generate_multimodal_plots(clip_data_summarized, "CLIP")

`geom_smooth()` using formula = 'y ~ x'
`geom_smooth()` using formula = 'y ~ x'
`geom_smooth()` using formula = 'y ~ x'

clip_plots_low_aoa <- generate_multimodal_plots(clip_data_summarized_low_aoa, "CLIP AoA < 4.78")

`geom_smooth()` using formula = 'y ~ x'
`geom_smooth()` using formula = 'y ~ x'
`geom_smooth()` using formula = 'y ~ x'

input_median_age = median_age
clip_age_plots <- generate_multimodal_age_effect_plots(clip_age_half_summarized, model_type="CLIP", median_age=input_median_age)

`geom_smooth()` using formula = 'y ~ x'
`geom_smooth()` using formula = 'y ~ x'
`geom_smooth()` using formula = 'y ~ x'

clip_plots

clip_age_plots

# CVCL
create_model_plots(all_similarities |> transmute(image_similarity=cvcl_image_similarity, text_similarity=cvcl_text_similarity, multimodal_similarity=cvcl_multimodal_similarity, text1, text2), name="CVCL", median_age=input_median_age)

`geom_smooth()` using formula = 'y ~ x'
`geom_smooth()` using formula = 'y ~ x'
`geom_smooth()` using formula = 'y ~ x'
`geom_smooth()` using formula = 'y ~ x'
`geom_smooth()` using formula = 'y ~ x'
`geom_smooth()` using formula = 'y ~ x'

create_model_plots(all_similarities |> transmute(image_similarity=dino_imagenet100_vitb14_image_similarity, text1, text2), name="ImageNetVIT", median_age=input_median_age)

`geom_smooth()` using formula = 'y ~ x'
`geom_smooth()` using formula = 'y ~ x'

create_model_plots(all_similarities |> transmute(image_similarity=dino_say_vitb14_image_similarity, text1, text2), name="SayCamVIT", median_age=input_median_age)

`geom_smooth()` using formula = 'y ~ x'
`geom_smooth()` using formula = 'y ~ x'

create_model_plots(all_similarities |> transmute(image_similarity=dinov3.babyview_image_similarity, text1, text2), name="BV-DINO", median_age=input_median_age)

`geom_smooth()` using formula = 'y ~ x'
`geom_smooth()` using formula = 'y ~ x'

create_model_plots(all_similarities |> transmute(image_similarity=clip_image_similarity, text_similarity=clip_text_similarity, multimodal_similarity=clip_multimodal_similarity,text1, text2), name="CLIP", median_age=input_median_age)

`geom_smooth()` using formula = 'y ~ x'
`geom_smooth()` using formula = 'y ~ x'
`geom_smooth()` using formula = 'y ~ x'
`geom_smooth()` using formula = 'y ~ x'
`geom_smooth()` using formula = 'y ~ x'
`geom_smooth()` using formula = 'y ~ x'

combined

bv_dino_data <- looking_data_summarized |> left_join(all_similarities, by=c("Trials.targetImage"="text1", "Trials.distractorImage"="text2")) |> rename("BabyView Similarity"="dinov3.babyview_image_similarity")
group_vars = c("Trials.trialID", "Trials.targetImage", "Trials.distractorImage",  "image_similarity", "text_similarity", "multimodal_similarity", "BabyView Similarity", "section")
dino_data_summarized <- (summarized_data(
  bv_dino_data,
  "Trials.trialID", 
  "corrected_target_looking", 
  group_vars
))
dino_data_long <- dino_data_summarized |>
  rename(
    `Text Similarity` = text_similarity,
    `Image Similarity` = image_similarity,
    `Multimodal Similarity` = multimodal_similarity
  ) |>
  pivot_longer(
    cols = c("Text Similarity", "Image Similarity", "Multimodal Similarity", "BabyView Similarity"),
    names_to = "sim_type",
    values_to = "similarity"
  ) |>
  mutate(
    sim_type = factor(
      sim_type,
      levels = c("Text Similarity", "Image Similarity", "Multimodal Similarity", "BabyView Similarity")
    )
  )

all_sections_plot <- multiple_similarity_effects_plot_refactored(
  dino_data_long,
  x_var = "similarity",
  group_var = "sim_type",
  input_title = "Looking time and embedding correlations",
  label_filter = "Trials.targetImage == 'acorn'",
  facet_scales_x = list(
    `Text Similarity` = scale_x_continuous(
      breaks = seq(0.7, 0.9, by = 0.1),
      limits = c(0.62, 0.93)
    ),
    `Image Similarity` = scale_x_continuous(
      breaks = seq(0.4, 0.8, by = 0.1),
      limits = c(0.38, 0.86)
    ),
    `Multimodal Similarity` = scale_x_continuous(
      breaks = seq(0.1, 0.4, by=0.1),
      limits = c(0.05, 0.45)
    ),
        `Babyview Similarity` = scale_x_continuous(breaks = seq(0, 0.4, by = 0.2), limits = c(-0.1, 0.5))
  ),
  color_map = c("Text Similarity" = "#C76D48", "Multimodal Similarity" = "#749B69"),
  n_col = 4
)
all_sections_plot

`geom_smooth()` using formula = 'y ~ x'
`geom_smooth()` using formula = 'y ~ x'

ggsave(here("figures",PROJECT_VERSION,"clip_plot_all.svg"),all_sections_plot,width=20, height=10,bg = "white",device="pdf")

`geom_smooth()` using formula = 'y ~ x'
`geom_smooth()` using formula = 'y ~ x'

ggsave(here("figures",PROJECT_VERSION,"clip_plot_all.png"),all_sections_plot,width=20, height=10,bg = "white")

`geom_smooth()` using formula = 'y ~ x'
`geom_smooth()` using formula = 'y ~ x'

scaled axes.

dino_data_long_scaled <- dino_data_summarized |>
  transmute(
    `Text Similarity` = scale(text_similarity),
    `Image Similarity` = scale(image_similarity),
    `Multimodal Similarity` = scale(multimodal_similarity),
    `BabyView Similarity` = scale(`BabyView Similarity`),
    mean_value = mean_value,
    Trials.targetImage = Trials.targetImage,
    Trials.distractorImage = Trials.distractorImage
  ) |>
  pivot_longer(
    cols = c("Text Similarity", "Image Similarity", "Multimodal Similarity", "BabyView Similarity"),
    names_to = "sim_type",
    values_to = "similarity"
  ) |>
  mutate(
    sim_type = factor(
      sim_type,
      levels = c("Text Similarity", "Image Similarity", "Multimodal Similarity", "BabyView Similarity")
    )
  )

all_sections_plot <- multiple_similarity_effects_plot_refactored(
  dino_data_long_scaled,
  x_var = "similarity",
  group_var = "sim_type",
  input_title = "Looking time and embedding correlations",
  label_filter = "Trials.targetImage == 'acorn'",
  facet_scales_x = list(
    `Text Similarity` = scale_x_continuous(
      breaks = seq(-2, 2, by = 0.5),
      limits = c(-2, 2)
    ),
    `Image Similarity` = scale_x_continuous(
      breaks = seq(-2, 2, by = 0.5),
      limits = c(-2, 2)
    ),
    `Multimodal Similarity` = scale_x_continuous(
      breaks = seq(-2, 2, by = 0.5),
      limits = c(-2, 2)
    ),
        `Babyview Similarity` = scale_x_continuous(
      breaks = seq(-2, 2, by = 0.5),
      limits = c(-2, 2)
          )
  ),
  color_map = c("Text Similarity" = "#C76D48", "Multimodal Similarity" = "#749B69"),
  n_col = 4
)
all_sections_plot

`geom_smooth()` using formula = 'y ~ x'
`geom_smooth()` using formula = 'y ~ x'

section 1 plot

section_1_plot <- multiple_similarity_effects_plot_refactored(
  dino_data_long |> filter(section == "sample1" & sim_type %in% c("Image Similarity", "Text Similarity")) |> mutate(
    sim_type = factor(
      sim_type,
      levels = c("Text Similarity", "Image Similarity")
    )
  ),
  x_var = "similarity",
  group_var = "sim_type",
  input_title = "Looking time and embedding correlations",
  label_filter = "Trials.targetImage == 'acorn' | Trials.distractorImage == 'acorn'",
  facet_scales_x = list(
    `Text Similarity` = scale_x_continuous(
      breaks = seq(0.75, 0.9, by = 0.05),
      limits = c(0.7, 0.91)
    ),
    `Image Similarity` = scale_x_continuous(
      breaks = seq(0.5, 0.9, by = 0.1),
      limits = c(0.45, 0.9)
    )
  ),
  color_map = c("Text Similarity" = "#C76D48")
)
section_1_plot

`geom_smooth()` using formula = 'y ~ x'
`geom_smooth()` using formula = 'y ~ x'

ggsave(here("figures",PROJECT_VERSION,"clip_plot_section11.svg"),section_1_plot,width=15, height=10,bg = "white",device="pdf")

`geom_smooth()` using formula = 'y ~ x'
`geom_smooth()` using formula = 'y ~ x'

ggsave(here("figures",PROJECT_VERSION,"clip_plot_section1.png"),section_1_plot,width=15, height=10,bg = "white")

`geom_smooth()` using formula = 'y ~ x'
`geom_smooth()` using formula = 'y ~ x'

correlation between bv and CLIP models

clip_bv_df <- all_sims |>
  mutate(
    fitted = predict(lm(image_similarity ~ dinov3.babyview_image_similarity, data = all_sims)),
    residual = abs(image_similarity - fitted),
    is_outlier = residual > (mean(residual) + 0.2*sd(residual))
  ) 

clip_bv <- ggplot(clip_bv_df, aes(x=dinov3.babyview_image_similarity, y=image_similarity)) +
  geom_smooth(method="lm", color="#215D89", alpha=0.5, se=FALSE) +
  geom_smooth(method="lm", color="#215D89", se=TRUE, alpha=0.2) +
  geom_point(size=8, alpha=0.8, color="#215D89") +
  scale_y_continuous(breaks = seq(0.4, 0.8, by = 0.2)) +
  geom_label_repel(
    data = clip_bv_df |> filter(is_outlier),
    aes(label=Trials.imagePair),
    color="black",
    fill="white",
    alpha=0.8,
    segment.alpha=0.7,
    nudge_y=0.02,
    force=10,
    force_pull=0.1,
    segment.size=1.2,
    size=6
  ) +
  ylab("CLIP image similarity") +
  xlab("BabyView model similarity") +
  theme_classic() +
  theme(
    text = element_text(size=14, face="bold"),
    axis.title.x = element_text(face="bold", size=20, margin=margin(t=15, r=0, b=0, l=0)),
    axis.title.y = element_text(face="bold", size=20, margin=margin(t=0, r=10, b=0, l=0)),
    axis.text = element_text(size=16, face="bold"),
    legend.key = element_blank(),
    legend.title = element_text(size=22, face="bold"),
    legend.text = element_text(size=22, face="bold"),
    legend.position = "bottom",
    panel.spacing = unit(0.5, "cm")
  )
clip_bv

`geom_smooth()` using formula = 'y ~ x'
`geom_smooth()` using formula = 'y ~ x'

ggsave(here("figures", PROJECT_VERSION,"clip_bv.svg"),clip_bv,width=6.67, height=6.67,units="in", bg = "white",device="pdf")

`geom_smooth()` using formula = 'y ~ x'
`geom_smooth()` using formula = 'y ~ x'

ggsave(here("figures", PROJECT_VERSION,"clip_bv.png"),clip_bv,width=6.67, height=6.67,units="in", bg = "white")

`geom_smooth()` using formula = 'y ~ x'
`geom_smooth()` using formula = 'y ~ x'

Some exploratory age effects with three age bins

clip_age_bins_summarized <- looking_data_summarized |>
  mutate(age_bin = case_when(
    SubjectInfo.testAge < 19*30 ~ "14-18",
    SubjectInfo.testAge < 22*30 ~ "18-21",
    TRUE ~ "21-24")) |>
  summarize_similarity_data(extra_fields = c("age_bin", "AoA_Est_target"))

similarity_age_plot(clip_age_bins_summarized, "image_similarity", age_grouping="age_bin", model_type="CLIP", median_age=median_age)

`geom_smooth()` using formula = 'y ~ x'

ggplot(clip_age_half_summarized |> left_join(all_similarities |> dplyr::select(text1, text2, dinov3.babyview_image_similarity), by=c("Trials.targetImage"="text1", "Trials.distractorImage"="text2")), aes(x = dinov3.babyview_image_similarity, y = mean_value, 
                   color = age_half)) +
    geom_hline(yintercept = 0, linetype = "dashed", alpha = 0.7) +
    geom_point(size = 3, alpha = 0.7) +
    geom_smooth(method = "lm", se = TRUE, alpha = 0.3) +
    
    # Improved styling
    scale_color_brewer(palette = "Set1", name="Age half") +
    theme_minimal() +
    theme(
      legend.position = "bottom",
      axis.title = element_text(size = 11),
      legend.title = element_text(size = 10),
      plot.title = element_text(size = 12)
    ) +
    
    # Clear labels
    labs(
      title = "BabyView DINO correlation across age",
      x = paste("Target-Distractor", "Image Similarity"),
      y = "Baseline-corrected Proportion\nTarget Looking",
      caption = paste0("Median age: ", median_age, " months")
    ) 

`geom_smooth()` using formula = 'y ~ x'

AoA effects

input_median_age = median_age
facet_plot <- generate_aoa_facet_plots(
  data = clip_age_half_summarized, 
  model_type = "CLIP",
  median_age = input_median_age,  
  age_grouping = "age_half",
  include_multimodal = FALSE
)

`geom_smooth()` using formula = 'y ~ x'

# For separate plots by AoA group:
#low_aoa_plots <- generate_multimodal_age_effect_plots(
#  data = clip_age_half_summarized %>% filter(AoA_Est_target <= 4.44),
#  model_type = "CLIP",
#  median_age = input_median_age,  
#  age_grouping = "age_half",
#  include_multimodal = FALSE
#)
facet_plot

`geom_smooth()` using formula = 'y ~ x'

aoa_target_median <- median(
  trial_metadata$AoA_Est_target,
  na.rm = TRUE
)

aoa_diff_median <- median(
  trial_metadata |>
    mutate(aoa_diff = AoA_Est_target - AoA_Est_distractor) |>
    pull(aoa_diff),
  na.rm = TRUE
)

# Split by target AoA
aoa_halved <- looking_data_summarized |>
  mutate(
    aoa_half = case_when(
      AoA_Est_target <= aoa_target_median ~ "lower_aoa",
      TRUE ~ "higher_aoa"
    )
  ) |>
  summarize_similarity_data(extra_fields = "aoa_half")

plot_grid(
  age_half_plot(
    aoa_halved,
    group_var = "aoa_half",
    x_var = "text_similarity",
    x_label = "text similarity"
  ),
  age_half_plot(
    aoa_halved,
    group_var = "aoa_half",
    x_var = "image_similarity",
    x_label = "image similarity"
  )
)

`geom_smooth()` using formula = 'y ~ x'
`geom_smooth()` using formula = 'y ~ x'

check for interactions here. hypothesis: higher-level text similarity effect for lower aoa items for older kids and lower-level image similarity effect for higher aoa items for younger kids.

AoA Difference

# Split by target-distractor AoA difference
aoa_diff_halved <- looking_data_summarized |>
  mutate(
    aoa_diff = AoA_Est_target - AoA_Est_distractor,
    aoa_diff_half = case_when(
      aoa_diff <= aoa_diff_median ~ "lower_aoa_diff",
      TRUE ~ "higher_aoa_diff"
    )
  ) |>
  summarize_similarity_data(extra_fields = "aoa_diff_half")

plot_grid(
  age_half_plot(
    aoa_diff_halved,
    group_var = "aoa_diff_half",
    x_var = "text_similarity",
    x_label = "text similarity"
  ),
  age_half_plot(
    aoa_diff_halved,
    group_var = "aoa_diff_half",
    x_var = "image_similarity",
    x_label = "image similarity"
  )
)

`geom_smooth()` using formula = 'y ~ x'
`geom_smooth()` using formula = 'y ~ x'

was there a difference between samples? we are seeing a stronger effect in sample 1 compared to sample 2

library(dplyr)
library(ggplot2)

plot_data <- trial_metadata |>
  mutate(
    aoa_diff = AoA_Est_target - AoA_Est_distractor
  )

summary_data <- plot_data |>
  group_by(section) |>
  summarise(
    mean_aoa_diff = mean(aoa_diff, na.rm = TRUE),
    sd_aoa_diff = sd(aoa_diff, na.rm = TRUE),
    n = sum(!is.na(aoa_diff)),
    se = sd_aoa_diff / sqrt(n),
    ci_lower = mean_aoa_diff - qt(0.975, df = n - 1) * se,
    ci_upper = mean_aoa_diff + qt(0.975, df = n - 1) * se,
    .groups = "drop"
  )

ggplot(plot_data, aes(x = aoa_diff, y = section)) +

  # Individual observations
  geom_jitter(
    height = 0.12,
    width = 0,
    alpha = 0.35,
    size = 2
  ) +

  # 95% CI
  geom_errorbar(
    data = summary_data,
    aes(
      xmin = ci_lower,
      xmax = ci_upper,
      y = section
    ),
    width = 0.12,
    linewidth = 1,
    inherit.aes = FALSE
  ) +

  # Mean
  geom_point(
    data = summary_data,
    aes(
      x = mean_aoa_diff,
      y = section
    ),
    size = 4,
    inherit.aes = FALSE
  ) +

  labs(
    x = "Target-distractor AoA difference",
    y = NULL
  ) +

  theme_classic(base_size = 14)

so the issue might not be higher or greater but about difference from centroid. let me divide by that convention instead.

aoa_diff_halved <- looking_data_summarized |>
  mutate(
    aoa_diff = AoA_Est_target - AoA_Est_distractor,
    aoa_diff_half = case_when(
      aoa_diff <= -0.32 ~ "high_aoa_diff",
      aoa_diff >= 0.32 ~ "high_aoa_diff",
      TRUE ~ "lower_aoa_diff"
    )
  ) |>
  summarize_similarity_data(extra_fields = "aoa_diff_half")

plot_grid(
  age_half_plot(
    aoa_diff_halved,
    group_var = "aoa_diff_half",
    x_var = "text_similarity",
    x_label = "text similarity"
  ),
  age_half_plot(
    aoa_diff_halved,
    group_var = "aoa_diff_half",
    x_var = "image_similarity",
    x_label = "image similarity"
  )
)

`geom_smooth()` using formula = 'y ~ x'
`geom_smooth()` using formula = 'y ~ x'

aoa_diff_halved |> count(aoa_diff_half)

# A tibble: 2 × 2
  aoa_diff_half      n
  <chr>          <int>
1 high_aoa_diff     32
2 lower_aoa_diff    32

High AoA difference may make image similarity effects appear. However, to be noted that there are many confounding variables. Could just be random effect of sample 1.

section effects

clip_section_summarized <- looking_data_summarized |>
  summarize_similarity_data(extra_fields = c("section", "AoA_Est_target"))

similarity_age_plot(clip_section_summarized, "image_similarity",
                    age_grouping = "section", model_type = "CLIP",
                    median_age = median_age)

`geom_smooth()` using formula = 'y ~ x'

OpenCLIP training

looking_data_w_openclip <- looking_data_summarized |>
  dplyr::select(-text_similarity, -multimodal_similarity, -image_similarity) |>
  left_join(openclip_similarities, by = c("Trials.distractorImage"="text2", "Trials.targetImage"="text1"))

Warning in left_join(dplyr::select(looking_data_summarized, -text_similarity, : Detected an unexpected many-to-many relationship between `x` and `y`.
ℹ Row 1 of `x` matches multiple rows in `y`.
ℹ Row 51 of `y` matches multiple rows in `x`.
ℹ If a many-to-many relationship is expected, set `relationship =
  "many-to-many"` to silence this warning.

openclip_data_summarized <- summarize_similarity_data(looking_data_w_openclip, extra_fields=c("epoch")) |> filter(!is.na(image_similarity))
openclip_age_half_summarized <- looking_data_w_openclip |>
  add_age_split() |>
  summarize_similarity_data(extra_fields = c("age_half", "epoch")) |> filter(!is.na(image_similarity))

# Function to calculate Pearson's correlation per epoch
calculate_correlations <- function(data, x_var, y_var, group_var = c("epoch"), conf_level = 0.95, method="pearson") {
  data |>
    group_by(across(all_of(group_var))) |>
    summarize(
      {
        cor_test <- cor.test(.data[[x_var]], .data[[y_var]], method = method, conf.level = conf_level)
        tibble(
          pearson_cor = cor_test$estimate,
          p_value = cor_test$p.value,
          ci_lower = cor_test$conf.int[1],
          ci_upper = cor_test$conf.int[2]
        )
      },
      .groups = "drop"
    )
}

image_correlation_results_age_split <- calculate_correlations(openclip_age_half_summarized, "image_similarity", "mean_value", group_var=c("epoch", "age_half"))

text_correlation_results <- calculate_correlations(openclip_data_summarized, "text_similarity", "mean_value")
image_correlation_results <- calculate_correlations(openclip_data_summarized, "image_similarity", "mean_value")
multimodal_correlation_results <- calculate_correlations(openclip_data_summarized, "multimodal_similarity", "mean_value")

ggplot(image_correlation_results, aes(x = log(epoch), y = pearson_cor)) +
  geom_point(aes(color = p_value < 0.05), size = 3) + 
  geom_smooth(span = 2) +
  labs(title = "Image similarity correlation across Open-CLIP training",
       x = "log(Epoch)",
       y = "Pearson Correlation") +
  scale_color_manual(values = c("TRUE" = "black", "FALSE" = "gray")) +  # Set color for significance
  theme_minimal() +
  theme(legend.position = "none")

`geom_smooth()` using method = 'loess' and formula = 'y ~ x'

ggplot(text_correlation_results, aes(x = log(epoch), y = pearson_cor)) +
  geom_point(aes(color = p_value < 0.05), size = 3) + 
  geom_smooth(span = 2) +
  labs(title = "Word similarity correlation across Open-CLIP training",
       x = "log (Epoch)",
       y = "Pearson Correlation") +
  scale_color_manual(values = c("TRUE" = "black", "FALSE" = "gray")) +  # Set color for significance
  theme_minimal() +
  theme(legend.position = "none")

`geom_smooth()` using method = 'loess' and formula = 'y ~ x'

ggplot(multimodal_correlation_results, aes(x = log(epoch), y = pearson_cor)) +
  geom_point(aes(color = p_value < 0.05), size = 3) + 
  geom_smooth(span = 2) +
  labs(title = "Multimodal similarity correlation across Open-CLIP training",
       x = "log (Epoch)",
       y = "Pearson Correlation") +
  scale_color_manual(values = c("TRUE" = "black", "FALSE" = "gray")) +  # Set color for significance
  theme_minimal() +
  theme(legend.position = "none")

`geom_smooth()` using method = 'loess' and formula = 'y ~ x'

image_openclip_age <- epoch_age_half_plot(calculate_correlations(openclip_age_half_summarized , "image_similarity", "mean_value", group_var=c("epoch", "age_half")), "image_similarity")
text_openclip_age <- epoch_age_half_plot(calculate_correlations(openclip_age_half_summarized, "text_similarity", "mean_value", group_var=c("epoch", "age_half")), "text_similarity")
image_openclip_age

`geom_smooth()` using method = 'loess' and formula = 'y ~ x'

text_openclip_age

`geom_smooth()` using method = 'loess' and formula = 'y ~ x'

ggsave(here("figures",PROJECT_VERSION,"openclip_image_embeddings_age.png"),image_openclip_age,width=9,height=6,bg = "white")

`geom_smooth()` using method = 'loess' and formula = 'y ~ x'

# --- set to FALSE to drop multimodal ---
include_multimodal <- TRUE

sim_colors <- c(
  "Image"      = "#4E79A7",  
  "Text"       = "#C0724A",  
  "Multimodal" = "#749B69"   
)

combined <- bind_rows(
  image_correlation_results %>% mutate(similarity = "Image"),
  text_correlation_results  %>% mutate(similarity = "Text"),
  if (include_multimodal)
    multimodal_correlation_results %>% mutate(similarity = "Multimodal")
) %>%
  mutate(similarity = factor(similarity, levels = names(sim_colors)))

ggplot(combined, aes(log(epoch), pearson_cor, color = similarity)) +
  geom_point(size = 3, alpha = 0.9) +
  geom_smooth(span = 2, se = TRUE) +
  scale_color_manual(values = sim_colors, drop = TRUE) +
  labs(x = "log(OpenCLIP epoch)",
       y = "Pearson's correlation (r)",
       color = NULL) +
  theme_minimal(base_size = 14) +
  theme(legend.position = c(0.02, 0.98),
        legend.justification = c(0, 1),
        legend.background = element_rect(fill = "white", color = "grey80"))

`geom_smooth()` using method = 'loess' and formula = 'y ~ x'

ggsave(here("figures",PROJECT_VERSION,"openclip_all_measures.png"),width=9,height=6,bg = "white")

`geom_smooth()` using method = 'loess' and formula = 'y ~ x'

mixed-models

# Prepare looking data with all similarity types and epochs
looking_data_w_all_sims <- looking_data_summarized |>
  dplyr::select(-image_similarity, -text_similarity, -multimodal_similarity) |>
  left_join(
    openclip_similarities, 
    by = c("Trials.distractorImage" = "text2", "Trials.targetImage" = "text1"),
    relationship = "many-to-many"
  ) |>
  rename(image_similarity_openclip = image_similarity, text_similarity_openclip = text_similarity, multimodal_similarity_openclip = multimodal_similarity)

# If you want to include text and multimodal as well in the same plot:
# You'll need to fit separate models for each similarity type and combine them

unique_epochs <- looking_data_w_all_sims |>
  filter(!is.na(image_similarity_openclip)) |>
  distinct(epoch) |>
  arrange(epoch) |>
  pull(epoch)

sim_types_to_fit <- c("image_similarity_openclip", "text_similarity_openclip", "multimodal_similarity_openclip")

epoch_models_all_types <- expand_grid(
  epoch = unique_epochs,
  sim_type = sim_types_to_fit
) |>
  pmap_df(function(epoch, sim_type) {
    d <- looking_data_w_all_sims |> 
      filter(.data$epoch == !!epoch) |>
      mutate(age_in_months = SubjectInfo.testAge) |>
      drop_na(all_of(c(sim_type, "corrected_target_looking", "age_in_months", 
                        "AoA_Est_target", "MeanSaliencyDiff")))
    
    if (nrow(d) < 10) return(tibble())
    
    formula_str <- paste0(
      "scale(corrected_target_looking) ~ scale(",
      sim_type,
      ") * scale(age_in_months) + scale(AoA_Est_target) + scale(MeanSaliencyDiff) + ",
      "(scale(",
      sim_type,
      ") | SubjectInfo.subjID) + (1 | Trials.targetImage) + (1 | Trials.imagePair)"
    )
    
    tryCatch({
      m <- lmer(as.formula(formula_str), data = d, 
                control = lmerControl(optimizer = "bobyqa", optCtrl = list(maxfun = 2e4)))
      
      coefs <- summary(m)$coefficients
      var_name <- paste0("scale(", sim_type, ")")
      int_name <- paste0("scale(", sim_type, "):scale(age_in_months)")
      
      tibble(
        epoch = epoch,
        sim_type = case_when(
          sim_type == "image_similarity_openclip" ~ "Image",
          sim_type == "text_similarity_openclip" ~ "Text",
          sim_type == "multimodal_similarity_openclip" ~ "Multimodal"
        ),
        main_estimate = coefs[var_name, "Estimate"],
        main_se = coefs[var_name, "Std. Error"],
        main_p = coefs[var_name, "Pr(>|t|)"],
        interaction_estimate = coefs[int_name, "Estimate"],
        interaction_se = coefs[int_name, "Std. Error"],
        interaction_p = coefs[int_name, "Pr(>|t|)"]
      )
    }, error = function(e) {
      warning(paste("Error fitting epoch", epoch, "sim_type", sim_type))
      return(tibble())
    })
  })

boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')

if (nrow(epoch_models_all_types) == 0) {
  warning("No models fitted for multi-type analysis. Check data availability.")
}

# Color palette matching your combined plot
sim_colors <- c(
  "Image"      = "#4E79A7",  
  "Text"       = "#C0724A",  
  "Multimodal" = "#749B69"   
)

# Diagnostic for multi-type
cat("\nEpoch models (all types) summary:\n")

Epoch models (all types) summary:

if (nrow(epoch_models_all_types) > 0) {
  print(epoch_models_all_types |> count(sim_type))
} else {
  cat("No multi-type models fitted.\n")
}

# A tibble: 3 × 2
  sim_type       n
  <chr>      <int>
1 Image         32
2 Multimodal    32
3 Text          32

# Plot 1: All similarity types - main effects
if (nrow(epoch_models_all_types) > 0) {
  p_all_main <- ggplot(epoch_models_all_types, 
                       aes(x = log(epoch), y = main_estimate, color = sim_type)) +
    geom_hline(yintercept = 0, linetype = "dashed", color = "gray50") +
      geom_ribbon(
      aes(x = log(epoch),                          # ← explicitly added
          ymin = main_estimate - 1.96 * main_se,
          ymax = main_estimate + 1.96 * main_se,
          fill = sim_type),
      alpha = 0.15,
      color = NA,
      data = epoch_models_all_types
    )+
    geom_line(size = 1) +
    geom_point(aes(fill = main_p < 0.05), size = 2.5, shape = 21, stroke = 1.5) +
    scale_color_manual(values = sim_colors, guide = "none") +
    scale_fill_manual(
      values = c("TRUE" = "black", "FALSE" = "white"),
      guide = "none"
    ) +
    facet_wrap(~sim_type) +
    labs(
      title = "Main Effects of Similarity Measures Across OpenCLIP Training",
      x = "log(Epoch)",
      y = "β (scaled)"
    ) +
    theme_minimal(base_size = 12) +
    theme(
      panel.border = element_rect(fill = NA, color = "grey80")
    )

  # Plot 2: All similarity types - interactions with age
  p_all_int <- ggplot(epoch_models_all_types, 
                      aes(x = log(epoch), y = interaction_estimate, color = sim_type)) +
    geom_hline(yintercept = 0, linetype = "dashed", color = "gray50") +
        geom_ribbon(
      aes(x = log(epoch),               
          ymin = interaction_estimate - 1.96 * interaction_se,
          ymax = interaction_estimate + 1.96 * interaction_se,
          fill = sim_type),
      alpha = 0.15,
      color = NA,
      data = epoch_models_all_types
    )+
    geom_line(size = 1) +
    geom_point(aes(fill = interaction_p < 0.05), size = 2.5, shape = 21, stroke = 1.5) +
    scale_color_manual(values = sim_colors, guide = "none") +
    scale_fill_manual(
      values = c("TRUE" = "black", "FALSE" = "white"),
      guide = "none"
    ) +
    facet_wrap(~sim_type) +
    labs(
      title = "Similarity × Age Interactions Across OpenCLIP Training",
      x = "log(Epoch)",
      y = "β (scaled)"
    ) +
    theme_minimal(base_size = 12) +
    theme(
      panel.border = element_rect(fill = NA, color = "grey80")
    )

  p_all_main / p_all_int

  ggsave(
    here("figures", PROJECT_VERSION, "openclip_epoch_models_all_types.png"),
    width = 12,
    height = 10,
    bg = "white"
  )
} else {
  cat("Skipping multi-type plots - no data available\n")
}

layer-wise

layerwise_long <- layerwise |>
  pivot_longer(
    cols = starts_with("layer"),
    names_to = "layer",
    values_to = "image_similarity"
  ) |>
  mutate(layer = as.integer(str_extract(layer, "\\d+")))


looking_data_w_layerwise <- looking_data_summarized |>
  dplyr::select(-image_similarity) |>
  left_join(layerwise_long, by = c("Trials.distractorImage" = "text1", "Trials.targetImage" = "text2"), relationship="many-to-many")

layerwise_data_summarized <- summarize_similarity_data(looking_data_w_layerwise, extra_fields = c("layer"))

image_correlation_results_layerwise <- calculate_correlations(
  layerwise_data_summarized, "image_similarity", "mean_value", group_var = "layer", method="pearson"
)

ggplot(image_correlation_results_layerwise, aes(x = layer, y = pearson_cor)) +
  geom_point(aes(color = p_value < 0.05), size = 3) +
  geom_line() +
  geom_ribbon(aes(ymin = ci_lower, ymax = ci_upper), alpha = 0.15) +
  scale_x_continuous(breaks = 0:12) +
  scale_color_manual(values = c("TRUE" = "black", "FALSE" = "gray")) +
  geom_hline(yintercept=0, linetype="dotted") +
  labs(
    title = "Image similarity correlation across CLIP ViT layers",
    x = "Layer",
    y = "Pearson Correlation"
  ) +
  theme_classic(base_size=16) +
  theme(legend.position = "none")

ggsave(here("figures", PROJECT_VERSION, "clip_layerwise.png"), width=12, height=9, bg="white")

layerwise_models <- layerwise_long |>
  group_by(layer) |>
  group_split() |>
  map(function(layer_data) {
    layer_num <- unique(layer_data$layer)
    
  d <- looking_data_summarized |> mutate(age_in_months = SubjectInfo.testAge) |>
  dplyr::select(-any_of("image_similarity")) |>
  left_join(layer_data, by = c("Trials.distractorImage" = "text2", "Trials.targetImage" = "text1"),
            relationship = "many-to-many")
    
    m <- lmer(scale(corrected_target_looking) ~ scale(image_similarity) * scale(age_in_months) 
              + scale(AoA_Est_target)
              + scale(MeanSaliencyDiff)
              + (scale(image_similarity) | SubjectInfo.subjID)
              + (1 | Trials.targetImage)
              + (1 | Trials.imagePair),
              data = d)
    
    coefs <- summary(m)$coefficients
    
    tibble(
      layer = layer_num,
      image_estimate = coefs["scale(image_similarity)", "Estimate"],
      image_se = coefs["scale(image_similarity)", "Std. Error"],
      image_p = coefs["scale(image_similarity)", "Pr(>|t|)"],
      interaction_estimate = coefs["scale(image_similarity):scale(age_in_months)", "Estimate"],
      interaction_se = coefs["scale(image_similarity):scale(age_in_months)", "Std. Error"],
      interaction_p = coefs["scale(image_similarity):scale(age_in_months)", "Pr(>|t|)"]
    )
  }) |>
  bind_rows()

boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')

Warning in checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv, : Model failed to converge with max|grad| = 0.00308501 (tol = 0.002, component 1)
  See ?lme4::convergence and ?lme4::troubleshooting.

boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')
boundary (singular) fit: see help('isSingular')

# Plot: image_similarity main effect
p1 <- ggplot(layerwise_models, aes(x = layer, y = image_estimate)) +
  geom_hline(yintercept = 0, linetype = "dashed", color = "gray50") +
  geom_ribbon(aes(ymin = image_estimate - 1.96 * image_se,
                  ymax = image_estimate + 1.96 * image_se), alpha = 0.15) +
  geom_line() +
  geom_point(aes(color = image_p < 0.05), size = 3) +
  scale_color_manual(values = c("TRUE" = "black", "FALSE" = "gray")) +
  scale_x_continuous(breaks = 0:12) +
  labs(title = "Image similarity effect across ViT layers",
       x = "Layer", y = "β (scaled)") +
  theme_minimal() + theme(legend.position = "none")

# Plot: interaction with age
p2 <- ggplot(layerwise_models, aes(x = layer, y = interaction_estimate)) +
  geom_hline(yintercept = 0, linetype = "dashed", color = "gray50") +
  geom_ribbon(aes(ymin = interaction_estimate - 1.96 * interaction_se,
                  ymax = interaction_estimate + 1.96 * interaction_se), alpha = 0.15) +
  geom_line() +
  geom_point(aes(color = interaction_p < 0.05), size = 3) +
  scale_color_manual(values = c("TRUE" = "black", "FALSE" = "gray")) +
  scale_x_continuous(breaks = 0:12) +
  labs(title = "Image similarity × age interaction across ViT layers",
       x = "Layer", y = "β (scaled)") +
  theme_minimal() + theme(legend.position = "none")

p1 / p2

cor.test(all_sims$layer12_image_similarity, all_sims$image_similarity)

    Pearson's product-moment correlation

data:  all_sims$layer12_image_similarity and all_sims$image_similarity
t = 4.4909, df = 62, p-value = 3.151e-05
alternative hypothesis: true correlation is not equal to 0
95 percent confidence interval:
 0.2842322 0.6607680
sample estimates:
      cor 
0.4954275 

cor.test(all_sims$clip.hf_image_similarity, all_sims$layer12_image_similarity)

    Pearson's product-moment correlation

data:  all_sims$clip.hf_image_similarity and all_sims$layer12_image_similarity
t = 4.4885, df = 62, p-value = 3.178e-05
alternative hypothesis: true correlation is not equal to 0
95 percent confidence interval:
 0.2839927 0.6606212
sample estimates:
      cor 
0.4952309 

cor.test(all_sims$clip.hf_image_similarity, all_sims$image_similarity)

    Pearson's product-moment correlation

data:  all_sims$clip.hf_image_similarity and all_sims$image_similarity
t = 4419.2, df = 62, p-value < 2.2e-16
alternative hypothesis: true correlation is not equal to 0
95 percent confidence interval:
 0.9999974 0.9999990
sample estimates:
      cor 
0.9999984 

Examining linear correlation between similarity types

cor(data.frame(trial_metadata$multimodal_similarity, trial_metadata$image_similarity, trial_metadata$text_similarity, trial_metadata$cor))

                                     trial_metadata.multimodal_similarity
trial_metadata.multimodal_similarity                            1.0000000
trial_metadata.image_similarity                                 0.4456994
trial_metadata.text_similarity                                  0.6434561
trial_metadata.cor                                                     NA
                                     trial_metadata.image_similarity
trial_metadata.multimodal_similarity                       0.4456994
trial_metadata.image_similarity                            1.0000000
trial_metadata.text_similarity                             0.4070992
trial_metadata.cor                                                NA
                                     trial_metadata.text_similarity
trial_metadata.multimodal_similarity                      0.6434561
trial_metadata.image_similarity                           0.4070992
trial_metadata.text_similarity                            1.0000000
trial_metadata.cor                                               NA
                                     trial_metadata.cor
trial_metadata.multimodal_similarity                 NA
trial_metadata.image_similarity                      NA
trial_metadata.text_similarity                       NA
trial_metadata.cor                                    1

cor.test(trial_metadata$image_similarity, trial_metadata$text_similarity)

    Pearson's product-moment correlation

data:  trial_metadata$image_similarity and trial_metadata$text_similarity
t = 3.5095, df = 62, p-value = 0.0008424
alternative hypothesis: true correlation is not equal to 0
95 percent confidence interval:
 0.1792247 0.5935161
sample estimates:
      cor 
0.4070992 

cor.test(all_sims$dinov3.babyview_image_similarity, all_sims$image_similarity)

    Pearson's product-moment correlation

data:  all_sims$dinov3.babyview_image_similarity and all_sims$image_similarity
t = 3.4883, df = 62, p-value = 9e-04
alternative hypothesis: true correlation is not equal to 0
95 percent confidence interval:
 0.1768460 0.5919226
sample estimates:
      cor 
0.4050477 

cor.test(all_sims$dino_imagenet100_vitb14_image_similarity, all_sims$image_similarity, method="spearman")

Warning in cor.test.default(all_sims$dino_imagenet100_vitb14_image_similarity,
: cannot compute exact p-value with ties

    Spearman's rank correlation rho

data:  all_sims$dino_imagenet100_vitb14_image_similarity and all_sims$image_similarity
S = 17869, p-value = 2.742e-07
alternative hypothesis: true rho is not equal to 0
sample estimates:
      rho 
0.5909091 

library(GGally)

similarity_cols <- all_sims |>
  dplyr::select(matches("similarity$")) |>
  dplyr::select(-matches("^layer\\d+_"))

ggpairs(
  similarity_cols,
  upper = list(continuous = wrap("cor", method = "spearman", size = 3)),
  lower = list(continuous = wrap("points", alpha = 0.3, size = 0.5)),
  diag = list(continuous = "densityDiag")
)

Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties
Warning in cor.test.default(x, y, method = method): cannot compute exact
p-value with ties

AoA plots

aoa_jittered <- target_looking_item_level %>%
  rowwise() %>%
  mutate(aoa_x_jittered = AoA_Est_target + runif(n(), -0.01, 0.01))

aoa_across_time <- ggplot(aoa_jittered, aes(x=aoa_x_jittered, y=mean_value))+
  geom_hline(yintercept=0,linetype="dashed")+
  geom_point(size=8, alpha=0.5, color="#CF784B")+
  #geom_linerange(aes(ymin = mean_value - ci, ymax = mean_value + ci), alpha = 0.1) + 
  geom_smooth(method="lm", color="#CF784B", alpha=0.5)+
  geom_label_repel(size=7,aes(label=Trials.targetImage),force_pull=10,nudge_y=-0.03) +
  xlab("Estimated target word age-of-acquisition (AoA)")+
  ylab("Baseline-corrected\nprop. target looking")+
  #ggtitle("Prop. of target looking by estimated AoA of target words")+
  theme_classic()+
 theme(
    text = element_text(size=10,face = "bold"),              # All text bold
    # Increase space between x-axis and its title
    axis.title.x = element_text(
      face = "bold", 
      size = 18,
      margin = margin(t = 10, r = 0, b = 0, l = 0)
    ),
    
    # Increase space between y-axis and its title
    axis.title.y = element_text(
      face = "bold", 
      size= 18,
      margin = margin(t = 10, r = 5, b = 0, l = 0)
    ),      
    axis.text = element_text(size=16,face = "bold"),         # Axis text bold
    legend.title = element_text(size = 12, face = "bold"),  # Small bold legend title
    legend.text = element_text(size = 12, face = "bold")    # Small bold legend text
  )+
 scale_x_continuous(breaks = seq(3, 6, by = 0.5)) +
   coord_cartesian(xlim = c(3.2,6.5), ylim = c(-0.12, 0.18))
aoa_across_time

`geom_smooth()` using formula = 'y ~ x'

ggsave(here("figures",PROJECT_VERSION,"aoa_plot_new.svg"),aoa_across_time,width=7,height=7,bg = "white", device="pdf")

`geom_smooth()` using formula = 'y ~ x'

  #ggpubr::stat_cor(method = "spearman", label.y = max(target_looking_item_level$mean_value) + 0.05)


age_half_plot(target_looking_item_age_level, "AoA_Est_target", x_label="AoA of target word", title="Target looking by estimated AoA of target word across age half")

`geom_smooth()` using formula = 'y ~ x'

ggplot(target_looking_trial_level, aes(x=AoA_Est_distractor, y=mean_value)) +
  geom_point()+
  geom_label_repel(aes(label=Trials.imagePair)) +
  geom_smooth(method="gam")+
  geom_hline(yintercept=0,linetype="dashed")+
  xlab("Estimated AoA of distractor words")+
  ylab("Baseline-corrected prop. of target looking")+
  ggtitle("Prop. of target looking by estimated AoA of distractor words")+
  theme_minimal()+
  theme(axis.title.x = element_text(face="bold", size=15, vjust=-1),
        axis.text.x  = element_text(size=10,angle=0,vjust=0.5),
        axis.title.y = element_text(face="bold", size=15),
        axis.text.y  = element_text(size=10),
        strip.text.x = element_text(size = 10,face="bold"),
        
        ) +
  ggpubr::stat_cor(method = "spearman")

`geom_smooth()` using formula = 'y ~ s(x, bs = "cs")'

hist(target_looking_item_level$AoA_Est_target)

aoa_jittered <- target_looking_item_level %>%
  rowwise() %>%
  mutate(aoa_x_jittered = AoA_Est_target + runif(n(), -0.01, 0.01))

aoa_across_time <- ggplot(aoa_jittered, aes(x=aoa_x_jittered, y=mean_value))+
  geom_hline(yintercept=0,linetype="dashed")+
  geom_point(size=8, color="#E95454", alpha=0.6)+ #, aes(alpha=section)
  scale_alpha_manual(values = c(sample1 = 0.8, sample2 = 0.3)) +
  #geom_linerange(aes(ymin = mean_value - ci, ymax = mean_value + ci), alpha = 0.1) + 
  geom_smooth(method="lm", color="#E95454", alpha=0.5, se=FALSE) +
geom_smooth(method="lm", color="#E95454", se=TRUE, alpha=0.2) +
  geom_label_repel(data=aoa_jittered |> filter(AoA_Est_target > 5.7 | (AoA_Est_target < 5 & AoA_Est_target > 4.5) | Trials.targetImage %in% c("airplane", "mud", "worm", "car", "butterfly", "duck", "acorn", "coconut", "key", "truck") | (mean_value > 0.1 & AoA_Est_target > 4)), size=7,aes(label=Trials.targetImage),force_pull=15,nudge_y=-0.03) +
  xlab("Estimated target word age-of-acquisition (AoA)")+
  ylab("Baseline-corrected\nproportion target looking")+
  #ggtitle("Prop. of target looking by estimated AoA of target words")+
  theme_classic()+
  guides(shape = "none", alpha="none") +
 theme(
    text = element_text(size=10,face = "bold"),              # All text bold
    # Increase space between x-axis and its title
    axis.title.x = element_text(
      face = "bold", 
      size = 24,
      margin = margin(t = 10, r = 0, b = 0, l = 0)
    ),
    
    # Increase space between y-axis and its title
    axis.title.y = element_text(
      face = "bold", 
      size= 24,
      margin = margin(t = 10, r = 5, b = 0, l = 0)
    ),      
    axis.text = element_text(size=20),         # Axis text bold
    #legend.title = element_text(size = 12, face = "bold"),  # Small bold legend title
    #legend.text = element_text(size = 12, face = "bold")    # Small bold legend text
  )+
 scale_x_continuous(breaks = seq(3, 6, by = 0.5)) +
   coord_cartesian(xlim = c(3.2,6.5), ylim = c(-0.16, 0.18))
aoa_across_time

`geom_smooth()` using formula = 'y ~ x'
`geom_smooth()` using formula = 'y ~ x'

ggsave(here("figures",PROJECT_VERSION,"aoa_plot.svg"),aoa_across_time,height=6.936,width=9.437,units="in",bg = "white", device="pdf")

`geom_smooth()` using formula = 'y ~ x'
`geom_smooth()` using formula = 'y ~ x'

ggsave(here("figures",PROJECT_VERSION,"aoa_plot.png"),aoa_across_time,height=6.936,width=9.437,units="in",bg = "white")

`geom_smooth()` using formula = 'y ~ x'
`geom_smooth()` using formula = 'y ~ x'

aoa_data <- looking_data_summarized |> add_aoa_split() 
aoa_data$aoa_half <- factor(aoa_data$aoa_half, levels = c("lower", "higher"))

aoa_plot <- half_violins_plot(aoa_data, "aoa_half", "corrected_target_looking", "SubjectInfo.subjID", c("lower", "higher"), input_xlab = "Age-of-acquisition of target word", input_caption=paste0("AoA split at median=", median_aoa))

# A tibble: 2 × 8
  aoa_half mean_value sd_value     N      se     ci lower_ci upper_ci
  <fct>         <dbl>    <dbl> <int>   <dbl>  <dbl>    <dbl>    <dbl>
1 lower        0.0760    0.122   180 0.00912 0.0180   0.0580   0.0940
2 higher       0.0310    0.122   180 0.00913 0.0180   0.0130   0.0491

Warning in geom_path(aes(group = SubjectInfo.subjID), color = "black", fill =
NA, : Ignoring unknown parameters: `fill`

aoa_plot

summarize_trial <- looking_time_resampled_clean %>%
  filter(trial_exclusion == 0 & exclude_participant ==0 & exclude_participant_insufficient_data == 0) %>%
  filter(!is.na(accuracy_transformed)) %>%
  group_by(SubjectInfo.subjID, time_normalized_corrected, SubjectInfo.testAge) %>%
  summarized_data(., "Trials.trialID", "accuracy_transformed", c("time_normalized_corrected")) %>%
  rename(mean_accuracy = mean_value)

Warning: There were 904 warnings in `summarize()`.
The first warning was:
ℹ In argument: `ci = qt(0.975, N - 1) * sd_value/sqrt(N)`.
ℹ In group 1: `Trials.trialID = "easy-acorn-key"`, `time_normalized_corrected =
  -3800`.
Caused by warning in `qt()`:
! NaNs produced
ℹ Run `dplyr::last_dplyr_warnings()` to see the 903 remaining warnings.

 summarize_across_aoa <- summarize_trial |> left_join(trial_metadata) |> add_aoa_split() |>
  group_by(time_normalized_corrected, aoa_half) |>
  dplyr::summarize(n=n(),
            non_na_n = sum(!is.na(mean_accuracy)),
            mean_value=mean(mean_accuracy,na.rm=TRUE),
            sd_accuracy=sd(mean_accuracy,na.rm=TRUE),
            se_accuracy=sd_accuracy/sqrt(n),
            ci = ifelse(n > 1, qt(0.975, n - 1) * sd_accuracy / sqrt(n), NA)) 

Joining with `by = join_by(Trials.trialID)`
`summarise()` has regrouped the output.

looking_times_aoa <- ggplot(summarize_across_aoa,aes(time_normalized_corrected,mean_value,color=aoa_half))+
  xlim(-2000,4000)+
  geom_errorbar(aes(ymin=mean_value-ci,ymax=mean_value+ci),width=0, alpha=0.2)+
  #geom_point(alpha=0.2)+
    geom_smooth(method="gam")+
  geom_vline(xintercept=0,size=1.5)+
  geom_hline(yintercept=0.5,size=1.2,linetype="dashed")+
  geom_vline(xintercept=300,linetype="dotted")+
  ylim(0,1)+
  xlab("Time (normalized to target word onset) in ms")+
  ylab("Proportion Target Looking") +
 # labs(title="Proportion of looking time across time, split at mean age-of-acquisition") +
  scale_color_manual(values=c("#6BAED6", "#4292C6"), name="Target word \nage-of-acquisition (AoA)", labels=c( "AoA >= 4.78", "AoA < 4.78"))+
   guides(color = guide_legend(reverse = TRUE))

looking_times_aoa

`geom_smooth()` using formula = 'y ~ s(x, bs = "cs")'

Warning: Removed 254 rows containing non-finite outside the scale range
(`stat_smooth()`).

AoA correlations

aoa_wb_trials <- trial_metadata |> filter(!is.na(AoA_target_WB))
cor.test(aoa_wb_trials$AoA_target_WB, aoa_wb_trials$AoA_Est_target, method="spearman")

Warning in cor.test.default(aoa_wb_trials$AoA_target_WB,
aoa_wb_trials$AoA_Est_target, : cannot compute exact p-value with ties

    Spearman's rank correlation rho

data:  aoa_wb_trials$AoA_target_WB and aoa_wb_trials$AoA_Est_target
S = 603.22, p-value = 0.3121
alternative hypothesis: true rho is not equal to 0
sample estimates:
      rho 
0.2607651 

cor.test(trial_metadata$AoA_Est_target, trial_metadata$image_similarity, method="pearson")

    Pearson's product-moment correlation

data:  trial_metadata$AoA_Est_target and trial_metadata$image_similarity
t = 0.87156, df = 62, p-value = 0.3868
alternative hypothesis: true correlation is not equal to 0
95 percent confidence interval:
 -0.1395672  0.3464564
sample estimates:
      cor 
0.1100164 

subjects_data <- read.csv(file.path(PROCESSED_DATA_PATH, "level-participant_added-trials_data.csv")) |> filter(exclude_participant_insufficient_data == 0)
trial_data <- looking_data_summarized |> distinct(SubjectInfo.subjID, SubjectInfo.testAge)
subjects_data <- subjects_data |> 
  left_join(trial_data)

Joining with `by = join_by(SubjectInfo.subjID)`

hist(subjects_data$SubjectInfo.testAge/30, breaks=10)

Saliency correlations

target_looking_trial_age_level <- summarized_data(target_looking_trial_subject_level |> rename(mean_target_looking = mean_value) |> add_age_split(), "Trials.trialID", "mean_target_looking", c("AoA_Est_target", "age_half"))
trialid_age_level_summary <- target_looking_trial_age_level |> left_join(trial_metadata)

Joining with `by = join_by(Trials.trialID, AoA_Est_target)`

trialid_level_summary <- target_looking_trial_level |> left_join(trial_metadata)

Joining with `by = join_by(Trials.trialID, AoA_Est_distractor, AoA_Est_target,
Trials.targetImage, Trials.distractorImage, Trials.imagePair)`

age_half_plot(trialid_age_level_summary, "MeanTargetSaliency", x_label="Mean target saliency")

`geom_smooth()` using formula = 'y ~ x'

age_half_plot(trialid_age_level_summary, "MeanSaliencyDiff", x_label="Mean saliency difference", title="Target looking by mean saliency difference across age and trial types")

`geom_smooth()` using formula = 'y ~ x'

target_looking_image_pair <- summarized_data(looking_data_summarized, "Trials.imagePair", "corrected_target_looking", c("SubjectInfo.subjID", "SubjectInfo.testAge"))

Warning: There were 709 warnings in `summarize()`.
The first warning was:
ℹ In argument: `ci = qt(0.975, N - 1) * sd_value/sqrt(N)`.
ℹ In group 2: `Trials.imagePair = "acorn-coconut"`, `SubjectInfo.subjID =
  "VVI008"`, `SubjectInfo.testAge = 450`.
Caused by warning in `qt()`:
! NaNs produced
ℹ Run `dplyr::last_dplyr_warnings()` to see the 708 remaining warnings.

target_looking_imagepair_level <- summarized_data(target_looking_image_pair |> rename("mean_looking"="mean_value") |> add_age_split(), "Trials.imagePair", "mean_looking", c("age_half")) |> left_join(trial_metadata |> distinct(Trials.imagePair, .keep_all=TRUE)) |> mutate(MeanSaliencyDiff = abs(MeanSaliencyDiff))

Joining with `by = join_by(Trials.imagePair)`

age_half_plot(target_looking_imagepair_level, "MeanSaliencyDiff", x_label="Mean saliency difference", title="Target looking by mean saliency difference across age half and image pairs")

`geom_smooth()` using formula = 'y ~ x'

Saliency vs baseline

baseline_looking_trial_age_level <- summarized_data(looking_data_summarized |> add_age_split(), "Trials.trialID", "mean_target_looking_baseline_window", c("age_half")) |> left_join(trial_metadata)

Joining with `by = join_by(Trials.trialID)`

baseline_looking_trial_level <- summarized_data(looking_data_summarized |> rename("mean_looking"="mean_target_looking_critical_window"), "Trials.trialID", "mean_looking", c("Trials.trialID")) |> left_join(trial_metadata)

Joining with `by = join_by(Trials.trialID)`

age_half_plot(baseline_looking_trial_age_level, "MeanSaliencyDiff", x_label="Mean saliency difference", y_label = "Target looking in baseline window")

`geom_smooth()` using formula = 'y ~ x'

ggplot(baseline_looking_trial_level, aes(x=MeanSaliencyDiff, y=mean_value)) +
  geom_smooth(method="lm") +
  geom_point() +
  ylab("Baseline looking") +
  xlab("Mean target-distractor saliency difference") +
  ggpubr::stat_cor()

`geom_smooth()` using formula = 'y ~ x'

# Saliency plots
# Prepare the data by combining both datasets with a grouping variable
baseline_looking_trial_level <- summarized_data(looking_data_summarized |> rename("mean_looking"="mean_target_looking_baseline_window"), "Trials.trialID", "mean_looking", c("Trials.trialID", "section")) |> left_join(trial_metadata) 

Joining with `by = join_by(Trials.trialID, section)`

critical_looking_trial_level <- summarized_data(looking_data_summarized |> rename("mean_looking"="mean_target_looking_critical_window"), "Trials.trialID", "mean_looking", c("Trials.trialID", "section")) |> left_join(trial_metadata)

Joining with `by = join_by(Trials.trialID, section)`

trialid_level_summary$plot_type <- "Baseline-corrected"
baseline_looking_trial_level$plot_type <- "Baseline window"
critical_looking_trial_level$plot_type <- "Critical window"

combined_data <- bind_rows(trialid_level_summary, baseline_looking_trial_level, critical_looking_trial_level)

# Create labels for trials where target is acorn and distractor is coconut
combined_data$label_text <- ifelse(
  combined_data$Trials.targetImage == "acorn",
  paste("Target:", combined_data$Trials.targetImage, "\nDistractor:", combined_data$Trials.distractorImage),
  ""
)

label_data <- subset(combined_data, label_text != "")

# Per-facet dashed reference lines
hline_data <- data.frame(
  plot_type = c("Baseline window", "Critical window", "Baseline-corrected"),
  yint      = c(0.5, 0.5, 0.0)
)

combined_saliency_plot <- ggplot(combined_data, aes(x = MeanSaliencyDiff, y = mean_value)) +
  geom_smooth(method = "lm", color = "#182B49", fill = "#182B49", alpha = 0.15) +
  geom_hline(data = hline_data, aes(yintercept = yint),
             linetype = "dashed", color = "gray40", linewidth = 1) +
  geom_point(size = 8, color = "#182B49", alpha = 0.8) +
  #scale_alpha_manual(values = c(sample1 = 0.8, sample2 = 0.3)) +
  geom_label_repel(data = label_data,
                   aes(label = label_text),
                   color = "black",
                   fill = "white",
                   alpha = 0.8,
                   segment.alpha = 0.7,
                   nudge_y = 0.02,
                   force = 10,
                   force_pull = 0.1,
                   segment.size = 1.2,
                   max.overlaps = Inf,
                   size = 8) +
  facet_wrap(~ plot_type, scales = "free_y", ncol = 3) +
  xlab("Target-distractor GBVS mean saliency difference") +
  ylab("Proportion of target looking") +
  theme_minimal() +
  theme(
    text = element_text(size = 16, face = "bold"),
    axis.title.x = element_text(face = "bold", size = 29,
                                margin = margin(t = 15, r = 0, b = 0, l = 0)),
    legend.key = element_blank(),
    axis.title.y = element_text(face = "bold", size = 29,
                                margin = margin(t = 0, r = 10, b = 0, l = 0)),
    axis.text = element_text(size = 24, face = "bold"),
    legend.title = element_text(size = 22, face = "bold"),
    legend.text = element_text(size = 22, face = "bold"),
    legend.position = "bottom",
    strip.text = element_text(size = 28, face = "bold"),
    strip.background = element_rect(fill = "gray90", color = NA),
    strip.text.x = element_text(margin = margin(t = 8, b = 8)),
    panel.spacing = unit(0.5, "cm")
  )

combined_saliency_plot

`geom_smooth()` using formula = 'y ~ x'

# Save the combined plot
ggsave(here("figures", PROJECT_VERSION, "saliency_plot_combined_full.svg"), 
       combined_saliency_plot, width = 15, height = 12, bg = "white", device = "pdf")

`geom_smooth()` using formula = 'y ~ x'

ggsave(here("figures", PROJECT_VERSION, "saliency_plot_combined_full.png"), 
       combined_saliency_plot, width = 15, height = 12, bg = "white")

`geom_smooth()` using formula = 'y ~ x'

combined_image_similarity_plot <- ggplot(combined_data, aes(x = image_similarity, y = mean_value)) +
  geom_smooth(method = "lm", color = "#182B49", fill = "#182B49", alpha = 0.15) +
  geom_hline(data = hline_data, aes(yintercept = yint),
             linetype = "dashed", color = "gray40", linewidth = 1) +
  geom_point(size = 8, color = "#182B49", alpha = 0.8) +
  #scale_alpha_manual(values = c(sample1 = 0.8, sample2 = 0.3)) +
  geom_label_repel(data = label_data,
                   aes(label = label_text),
                   color = "black",
                   fill = "white",
                   alpha = 0.8,
                   segment.alpha = 0.7,
                   nudge_y = 0.02,
                   force = 10,
                   force_pull = 0.1,
                   segment.size = 1.2,
                   max.overlaps = Inf,
                   size = 8) +
  facet_wrap(~ plot_type, scales = "free_y", ncol = 3) +
  xlab("Target-distractor image similarity difference") +
  ylab("Proportion of target looking") +
  theme_minimal() +
  theme(
    text = element_text(size = 16, face = "bold"),
    axis.title.x = element_text(face = "bold", size = 29,
                                margin = margin(t = 15, r = 0, b = 0, l = 0)),
    legend.key = element_blank(),
    axis.title.y = element_text(face = "bold", size = 29,
                                margin = margin(t = 0, r = 10, b = 0, l = 0)),
    axis.text = element_text(size = 24, face = "bold"),
    legend.title = element_text(size = 22, face = "bold"),
    legend.text = element_text(size = 22, face = "bold"),
    legend.position = "bottom",
    strip.text = element_text(size = 28, face = "bold"),
    strip.background = element_rect(fill = "gray90", color = NA),
    strip.text.x = element_text(margin = margin(t = 8, b = 8)),
    panel.spacing = unit(0.5, "cm")
  )

combined_image_similarity_plot

`geom_smooth()` using formula = 'y ~ x'

# Save the combined plot
ggsave(here("figures", PROJECT_VERSION, "saliency_plot_combined_full.svg"), 
       combined_saliency_plot, width = 15, height = 12, bg = "white", device = "pdf")

`geom_smooth()` using formula = 'y ~ x'

ggsave(here("figures", PROJECT_VERSION, "image_sim_plot_combined_full.png"), 
       combined_image_similarity_plot, width = 15, height = 12, bg = "white")

`geom_smooth()` using formula = 'y ~ x'

Embedding type comparisons

cor.test(trial_metadata$text_similarity, trial_metadata$image_similarity)

    Pearson's product-moment correlation

data:  trial_metadata$text_similarity and trial_metadata$image_similarity
t = 3.5095, df = 62, p-value = 0.0008424
alternative hypothesis: true correlation is not equal to 0
95 percent confidence interval:
 0.1792247 0.5935161
sample estimates:
      cor 
0.4070992 

cor.test(trial_metadata$multimodal_similarity, trial_metadata$image_similarity)

    Pearson's product-moment correlation

data:  trial_metadata$multimodal_similarity and trial_metadata$image_similarity
t = 3.9204, df = 62, p-value = 0.0002236
alternative hypothesis: true correlation is not equal to 0
95 percent confidence interval:
 0.2244840 0.6232295
sample estimates:
      cor 
0.4456994 

cor.test(trial_metadata$AoA_Est_target, trial_metadata$text_similarity)

    Pearson's product-moment correlation

data:  trial_metadata$AoA_Est_target and trial_metadata$text_similarity
t = -0.47989, df = 62, p-value = 0.633
alternative hypothesis: true correlation is not equal to 0
95 percent confidence interval:
 -0.3021248  0.1877840
sample estimates:
        cor 
-0.06083327 

cor.test(trial_metadata$AoA_Est_target, trial_metadata$image_similarity)

    Pearson's product-moment correlation

data:  trial_metadata$AoA_Est_target and trial_metadata$image_similarity
t = 0.87156, df = 62, p-value = 0.3868
alternative hypothesis: true correlation is not equal to 0
95 percent confidence interval:
 -0.1395672  0.3464564
sample estimates:
      cor 
0.1100164 

cor.test(trial_metadata$MeanSaliencyDiff, trial_metadata$AoA_Est_target)

    Pearson's product-moment correlation

data:  trial_metadata$MeanSaliencyDiff and trial_metadata$AoA_Est_target
t = 0.21806, df = 62, p-value = 0.8281
alternative hypothesis: true correlation is not equal to 0
95 percent confidence interval:
 -0.2196208  0.2716437
sample estimates:
       cor 
0.02768289 

CDI

cdi_data <- read.csv(here("data", "main", "data_to_analyze", "level-participant_source-cdi_added-percentile_data.csv"))
looking_data_with_cdi <- looking_data_summarized |> left_join(cdi_data, by=c("SubjectInfo.subjID"="local_id"))
# split-half percentile_rank and plot 
# Median split on CDI percentile_rank
pr_median <- median(looking_data_with_cdi$percentile_rank, na.rm = TRUE)

clip_cdi_bins_summarized <- looking_data_with_cdi |>
  filter(!is.na(percentile_rank)) |>
  mutate(cdi_percentile_bin = if_else(percentile_rank < pr_median,
                           "Low percentile", "High percentile")) |>
  summarize_similarity_data(extra_fields = c("cdi_percentile_bin", "AoA_Est_target"))

similarity_age_plot(clip_cdi_bins_summarized, "image_similarity",
                    age_grouping = "cdi_percentile_bin", model_type = "CLIP",
                    median_age = median_age)

`geom_smooth()` using formula = 'y ~ x'

vocab_score_median <- median(looking_data_with_cdi$vocabulary_score, na.rm = TRUE)
clip_cdi_bins_summarized <- looking_data_with_cdi |>
  filter(!is.na(vocabulary_score)) |>
  mutate(cdi_vocab_bin = if_else(percentile_rank < vocabulary_score,
                           "Low percentile", "High percentile")) |>
  summarize_similarity_data(extra_fields = c("cdi_vocab_bin", "AoA_Est_target"))

similarity_age_plot(clip_cdi_bins_summarized, "image_similarity",
                    age_grouping = "cdi_vocab_bin", model_type = "CLIP",
                    median_age = median_age)

`geom_smooth()` using formula = 'y ~ x'

target_looking_subject_level_cdi <- target_looking_subject_level |> left_join(cdi_data, by=c("SubjectInfo.subjID"="local_id"))
ggplot(data=target_looking_subject_level_cdi, aes(x=mean_value, y=percentile_rank)) +
  geom_point(size=4, alpha=0.6) +
  geom_smooth(method="lm") + theme_classic() + xlab("Baseline-corrected proportion target looking") +
  ylab("Vocabulary percentile rank")

`geom_smooth()` using formula = 'y ~ x'

Warning: Removed 132 rows containing non-finite outside the scale range
(`stat_smooth()`).

Warning: Removed 132 rows containing missing values or values outside the scale range
(`geom_point()`).

cor.test(target_looking_subject_level_cdi$percentile_rank,target_looking_subject_level_cdi$vocabulary_score)

    Pearson's product-moment correlation

data:  target_looking_subject_level_cdi$percentile_rank and target_looking_subject_level_cdi$vocabulary_score
t = 7.9667, df = 46, p-value = 3.312e-10
alternative hypothesis: true correlation is not equal to 0
95 percent confidence interval:
 0.6090766 0.8595973
sample estimates:
      cor 
0.7614371 

all_cdi_data <- read.csv(here("data/main/data_to_analyze/level-words_source-cdi_data.csv"))
looking_data_summarized_with_knowledge <- looking_data_summarized |>
  left_join(all_cdi_data |> filter(type=="produces" & word_type=="extra_words") |> transmute(local_id, word, knows_target=responded), by=c("Trials.targetImage"="word", "SubjectInfo.subjID"="local_id"))

clip_cdi_knowledge_summarized <- looking_data_summarized_with_knowledge |>
  filter(!is.na(knows_target)) |>
  mutate(cdi_knowledge_bin = if_else(knows_target == 1,
                           "Knows", "Does not know")) |>
  summarize_similarity_data(extra_fields = c("cdi_knowledge_bin", "AoA_Est_target"))

Warning: There was 1 warning in `summarize()`.
ℹ In argument: `ci = qt(0.975, N - 1) * sd_value/sqrt(N)`.
ℹ In group 4: `Trials.trialID = "easy-balloon-dresser-distractor"`,
  `Trials.targetImage = "dresser"`, `Trials.distractorImage = "balloon"`,
  `image_similarity = 0.5211465`, `cor = NA`, `text_similarity = 0.6946648`,
  `multimodal_similarity = 0.4043515`, `cdi_knowledge_bin = "Knows"`,
  `AoA_Est_target = 4.28`.
Caused by warning in `qt()`:
! NaNs produced

similarity_age_plot(clip_cdi_knowledge_summarized, "image_similarity",
                    age_grouping = "cdi_knowledge_bin", model_type = "CLIP",
                    median_age = median_age)

`geom_smooth()` using formula = 'y ~ x'

Warning: Removed 1 row containing missing values or values outside the scale range
(`geom_segment()`).