### Visualizing observational mouse data:

#### Load in packages ####
library(chron)
library(hms)
library(ggplot2)
library(reshape2)
library(tidyr)
library(dplyr)
library(purrr)
library(reshape)
library(ggpubr)
library(ggrepel)
library(mgcv)
library(scales)
library(ggridges)
library(cowplot)

#### Bring in and clean the data ####
mousedata<-read.csv("mousedata.csv", na.strings = "NA")
str(mousedata)

#Create a new column for date and time
mousedata$date_time <- as.POSIXct(paste(mousedata$date_of_obs, mousedata$time_of_obs), tz="America/Chicago",
                                  format = "%m/%d/%Y %H:%M:%OS" )
mousedata$date_time

#tell R to recognize time column as time
mousedata$time_of_obs <- as_hms(mousedata$time_of_obs)
mousedata$time_of_obs
str(mousedata)
#sort to check time reformatting worked
sort(mousedata$time_of_obs)

#clarify format for date column:
mousedata$date_of_obs
mousedata$date_of_obs <- as.Date(mousedata$date_of_obs, format = "%m/%d/%Y")
mousedata$date_of_obs

#make replicate number categorical
mousedata$replicate <- as.factor(mousedata$replicate)

#correct cant_tell to numeric
mousedata$cant_tell_mean <- as.numeric(mousedata$cant_tell_mean)

#### Next, tidying "focal area" observations, to be able to quantify, also correct any spelling/vocab discrepancies ####

##separate the comma separated observations on focal area
mousedata$foc15_list <- strsplit(mousedata$focal_area_0to15, ",")
mousedata$foc30_list <- strsplit(mousedata$focal_area_15to30, ",")
mousedata$foc45_list <- strsplit(mousedata$focal_area_30to45, ",")
mousedata$foc60_list <- strsplit(mousedata$focal_area_45to60, ",")

## consolidate observations within the lists:
mousedata <- mousedata %>%
  mutate(foc15_list = map(foc15_list, ~ ifelse(.x == "rear legs", "rear leg", .x)))
mousedata <- mousedata %>%
  mutate(foc15_list = map(foc15_list, ~ ifelse(.x == "posterior back", "hip", .x)))
mousedata <- mousedata %>%
  mutate(foc15_list = map(foc15_list, ~ ifelse(.x == "hips", "hip", .x)))
mousedata <- mousedata %>%
  mutate(foc15_list = map(foc15_list, ~ ifelse(.x == "forelegs", "foreleg", .x)))
mousedata <- mousedata %>%
  mutate(foc15_list = map(foc15_list, ~ ifelse(.x == "leg", "rear leg", .x)))
mousedata <- mousedata %>%
  mutate(foc15_list = map(foc15_list, ~ ifelse(.x == "paws", "paw", .x)))
mousedata <- mousedata %>%
  mutate(foc15_list = map(foc15_list, ~ ifelse(.x == "eyes", "eye", .x)))
mousedata <- mousedata %>%
  mutate(foc15_list = map(foc15_list, ~ ifelse(.x == "front paw", "paw", .x)))
mousedata <- mousedata %>%
  mutate(foc15_list = map(foc15_list, ~ ifelse(.x == "front paws", "paw", .x)))
mousedata <- mousedata %>%
  mutate(foc15_list = map(foc15_list, ~ ifelse(.x == "genital", "genitals", .x)))
mousedata <- mousedata %>%
  mutate(foc15_list = map(foc15_list, ~ ifelse(.x == "rear", "hip", .x)))
mousedata <- mousedata %>%
  mutate(foc15_list = map(foc15_list, ~ ifelse(.x == "rear paw", "paw", .x)))


mousedata <- mousedata %>%
  mutate(foc30_list = map(foc30_list, ~ ifelse(.x == "rear legs", "rear leg", .x)))
mousedata <- mousedata %>%
  mutate(foc30_list = map(foc30_list, ~ ifelse(.x == "posterior back", "hip", .x)))
mousedata <- mousedata %>%
  mutate(foc30_list = map(foc30_list, ~ ifelse(.x == "hips", "hip", .x)))
mousedata <- mousedata %>%
  mutate(foc30_list = map(foc30_list, ~ ifelse(.x == "forelegs", "foreleg", .x)))
mousedata <- mousedata %>%
  mutate(foc30_list = map(foc30_list, ~ ifelse(.x == "leg", "rear leg", .x)))
mousedata <- mousedata %>%
  mutate(foc30_list = map(foc30_list, ~ ifelse(.x == "paws", "paw", .x)))
mousedata <- mousedata %>%
  mutate(foc30_list = map(foc30_list, ~ ifelse(.x == "eyes", "eye", .x)))
mousedata <- mousedata %>%
  mutate(foc30_list = map(foc30_list, ~ ifelse(.x == "front paw", "paw", .x)))
mousedata <- mousedata %>%
  mutate(foc30_list = map(foc30_list, ~ ifelse(.x == "front paws", "paw", .x)))
mousedata <- mousedata %>%
  mutate(foc30_list = map(foc30_list, ~ ifelse(.x == "genital", "genitals", .x)))
mousedata <- mousedata %>%
  mutate(foc30_list = map(foc30_list, ~ ifelse(.x == "rear", "hip", .x)))
mousedata <- mousedata %>%
  mutate(foc30_list = map(foc30_list, ~ ifelse(.x == "rear paw", "paw", .x)))
mousedata <- mousedata %>%
  mutate(foc30_list = map(foc30_list, ~ ifelse(.x == "abodmen", "abdomen", .x)))



mousedata <- mousedata %>%
  mutate(foc45_list = map(foc45_list, ~ ifelse(.x == "rear legs", "rear leg", .x)))
mousedata <- mousedata %>%
  mutate(foc45_list = map(foc45_list, ~ ifelse(.x == "posterior back", "hip", .x)))
mousedata <- mousedata %>%
  mutate(foc45_list = map(foc45_list, ~ ifelse(.x == "hips", "hip", .x)))
mousedata <- mousedata %>%
  mutate(foc45_list = map(foc45_list, ~ ifelse(.x == "forelegs", "foreleg", .x)))
mousedata <- mousedata %>%
  mutate(foc45_list = map(foc45_list, ~ ifelse(.x == "leg", "rear leg", .x)))
mousedata <- mousedata %>%
  mutate(foc45_list = map(foc45_list, ~ ifelse(.x == "paws", "paw", .x)))
mousedata <- mousedata %>%
  mutate(foc45_list = map(foc45_list, ~ ifelse(.x == "eyes", "eye", .x)))
mousedata <- mousedata %>%
  mutate(foc45_list = map(foc45_list, ~ ifelse(.x == "front paw", "paw", .x)))
mousedata <- mousedata %>%
  mutate(foc45_list = map(foc45_list, ~ ifelse(.x == "front paws", "paw", .x)))
mousedata <- mousedata %>%
  mutate(foc45_list = map(foc45_list, ~ ifelse(.x == "genital", "genitals", .x)))
mousedata <- mousedata %>%
  mutate(foc45_list = map(foc45_list, ~ ifelse(.x == "rear", "hip", .x)))
mousedata <- mousedata %>%
  mutate(foc45_list = map(foc45_list, ~ ifelse(.x == "rear paw", "paw", .x)))
mousedata <- mousedata %>%
  mutate(foc45_list = map(foc45_list, ~ ifelse(.x == "tial", "tail", .x)))

mousedata <- mousedata %>%
  mutate(foc60_list = map(foc60_list, ~ ifelse(.x == "rear legs", "rear leg", .x)))
mousedata <- mousedata %>%
  mutate(foc60_list = map(foc60_list, ~ ifelse(.x == "posterior back", "hip", .x)))
mousedata <- mousedata %>%
  mutate(foc60_list = map(foc60_list, ~ ifelse(.x == "hips", "hip", .x)))
mousedata <- mousedata %>%
  mutate(foc60_list = map(foc60_list, ~ ifelse(.x == "forelegs", "foreleg", .x)))
mousedata <- mousedata %>%
  mutate(foc60_list = map(foc60_list, ~ ifelse(.x == "leg", "rear leg", .x)))
mousedata <- mousedata %>%
  mutate(foc60_list = map(foc60_list, ~ ifelse(.x == "paws", "paw", .x)))
mousedata <- mousedata %>%
  mutate(foc60_list = map(foc60_list, ~ ifelse(.x == "eyes", "eye", .x)))
mousedata <- mousedata %>%
  mutate(foc60_list = map(foc60_list, ~ ifelse(.x == "front paw", "paw", .x)))
mousedata <- mousedata %>%
  mutate(foc60_list = map(foc60_list, ~ ifelse(.x == "front paws", "paw", .x)))
mousedata <- mousedata %>%
  mutate(foc60_list = map(foc60_list, ~ ifelse(.x == "genital", "genitals", .x)))
mousedata <- mousedata %>%
  mutate(foc60_list = map(foc60_list, ~ ifelse(.x == "rear", "hip", .x)))
mousedata <- mousedata %>%
  mutate(foc60_list = map(foc60_list, ~ ifelse(.x == "rear paw", "paw", .x)))


#### This is a function and new dataset to simplify these regions (e.g., eye = head). This is used for figures:####
replace_categories <- function(x) {
  x <- ifelse(x %in% c("abdomen", "chest", "hip", "neck", "rump", "shoulder", "back", "genitals"), "trunk", x)
  x <- ifelse(x %in% c("ear", "eye", "mouth", "nose", "face"), "head", x)
  x <- ifelse(x %in% c("foreleg", "paw", "rear leg", "legs"), "limbs", x)
  return(x)
}

# Apply the function to each list element. simplemousedata is same as mousedata, but the focal areas are simplified.
simplemousedata <- mousedata %>%
  mutate(across(c(foc15_list, foc30_list, foc45_list, foc60_list), ~ map(.x, replace_categories)))


#### Modify focal area data, whole and simplified ####
# Steps for analyzing focal area data:
foc15_data <- table(unlist(mousedata$foc15_list))
foc30_data <- table(unlist(mousedata$foc30_list))
foc45_data <- table(unlist(mousedata$foc45_list))
foc60_data <- table(unlist(mousedata$foc60_list))

foc15_simpdata <- table(unlist(simplemousedata$foc15_list))
foc30_simpdata <- table(unlist(simplemousedata$foc30_list))
foc45_simpdata <- table(unlist(simplemousedata$foc45_list))
foc60_simpdata <- table(unlist(simplemousedata$foc60_list))

foc15_unnested <- unnest(mousedata, foc15_list)
foc15_summ <- foc15_unnested %>%
  group_by(category = foc15_list) %>%
  summarise(count = n())

foc30_unnested <- unnest(mousedata, foc30_list)
foc30_summ <- foc30_unnested %>%
  group_by(category = foc30_list) %>%
  summarise(count = n())

foc45_unnested <- unnest(mousedata, foc45_list)
foc45_summ <- foc45_unnested %>%
  group_by(category = foc45_list) %>%
  summarise(count = n())

foc60_unnested <- unnest(mousedata, foc60_list)
foc60_summ <- foc60_unnested %>%
  group_by(category = foc60_list) %>%
  summarise(count = n())

#now, do the same for the simplified focal data:
  simp_foc15_unnested <- unnest(simplemousedata, foc15_list)
simp_foc15_summ <- simp_foc15_unnested %>%
  group_by(category = foc15_list) %>%
  summarise(count = n())

simp_foc30_unnested <- unnest(simplemousedata, foc30_list)
simp_foc30_summ <- simp_foc30_unnested %>%
  group_by(category = foc30_list) %>%
  summarise(count = n())

simp_foc45_unnested <- unnest(simplemousedata, foc45_list)
simp_foc45_summ <- simp_foc45_unnested %>%
  group_by(category = foc45_list) %>%
  summarise(count = n())

simp_foc60_unnested <- unnest(simplemousedata, foc60_list)
simp_foc60_summ <- simp_foc60_unnested %>%
  group_by(category = foc60_list) %>%
  summarise(count = n())


# Create a column that combines the four observation times to show the focal areas per 1-hour observation:
process_row <- function(row) {
  unique_values <- unique(unlist(row))
  if ("all" %in% unique_values) {
    return(list("all"))
  } else {
    return(list(unique_values))
  }
}

mousedata_processed <- mousedata %>%
  rowwise() %>%
  mutate(focal_area = process_row(c(foc15_list, foc30_list, foc45_list, foc60_list)))

simp_mousedata_processed <- simplemousedata %>%
  rowwise() %>%
  mutate(focal_area = process_row(c(foc15_list, foc30_list, foc45_list, foc60_list)))

mousedata_unnested <- mousedata_processed %>%
  unnest(cols = focal_area)

simp_mousedata_unnested <- simp_mousedata_processed %>%
  unnest(cols = focal_area)

longdata<-mousedata_unnested
simp_longdata<-simp_mousedata_unnested

longdata$replicate <- as.factor(longdata$replicate)
simp_longdata$replicate <- as.factor(simp_longdata$replicate)

longdata <- longdata %>%
  filter(!is.na(focal_area))

simp_longdata <- simp_longdata %>%
  filter(!is.na(focal_area))

longdata_rep1 <- longdata[longdata$replicate==1,]
longdata_rep2 <- longdata[longdata$replicate==2,]
longdata_rep3 <- longdata[longdata$replicate==3,]
longdata_rep4 <- longdata[longdata$replicate==4,]
longdata_rep5 <- longdata[longdata$replicate==5,]
longdata_rep6 <- longdata[longdata$replicate==6,]

simp_longdata_rep1 <- simp_longdata[simp_longdata$replicate==1,]
simp_longdata_rep2 <- simp_longdata[simp_longdata$replicate==2,]
simp_longdata_rep3 <- simp_longdata[simp_longdata$replicate==3,]
simp_longdata_rep4 <- simp_longdata[simp_longdata$replicate==4,]
simp_longdata_rep5 <- simp_longdata[simp_longdata$replicate==5,]
simp_longdata_rep6 <- simp_longdata[simp_longdata$replicate==6,]


###############################################################

#### Plotting breakdown of observations by time of day ####
# Convert hms to numeric seconds
mousedata$time_in_seconds <- as.numeric(mousedata$time_of_obs)

# Now plot observation time breakdown, mice 1-4 only
ggplot(data = subset(mousedata, !as.character(replicate) %in% c("6", "5")), aes(x = time_in_seconds, fill = as.factor(replicate))) +
  geom_histogram(bins = 16, position = "stack") +
  scale_fill_manual(values = c("#D81B60", "#1E88E5", "#FFC107", "#004D40")) +
  scale_x_continuous(
    breaks = seq(0, 79200, by = 3600),  # Every hour (3600 seconds)
    labels = function(x) {
      # Convert seconds to time format without leading zero in the hour
      as.character(format(as.POSIXct(x, origin = "1970-01-01", tz = "UTC"), "%l %p"))
    }
  ) +
  theme_classic() +
  labs(fill = "Mouse ID", x = "Time of Observation Start", y = "Total Observations")+
  theme(axis.text.x = element_text(angle = 45, hjust = 1))
# bins are equivalent to 30.1875 minute intervals


#### Preparing to plot behavioral counts over days ####
## first make a new column:
mousedata <- mousedata %>%
  mutate(
    chewing_mean = as.numeric(chewing_mean),
    tugging_mean = as.numeric(tugging_mean),
    sucking_mean = as.numeric(sucking_mean),
    carry_mean = as.numeric(carry_mean),
    cant_tell_mean = as.numeric(cant_tell_mean),
    mean_behav_total = chewing_mean + tugging_mean + sucking_mean + carry_mean + cant_tell_mean
  )

## re-format into a new dataset
mousedata_longer <- mousedata %>%
  pivot_longer(cols = c(chewing_mean, tugging_mean, sucking_mean, carry_mean, cant_tell_mean),
               names_to = "behavior", 
               values_to = "count")

#Create subsets to easily plot by mouse
mousedata_longer_1 <- subset(mousedata_longer, replicate == 1)
mousedata_longer_2 <- subset(mousedata_longer, replicate == 2)
mousedata_longer_3 <- subset(mousedata_longer, replicate == 3)
mousedata_longer_4 <- subset(mousedata_longer, replicate == 4)
mousedata_longer_5 <- subset(mousedata_longer, replicate == 5)

#### Making kernel density plots of overlapping behaviors, normalized to their own frequency, one per mouse (for Fig 2) ####

## Ensure there are no double-stacked behaviors from days with multiple observations:
mouse1_long_mean <- mousedata_longer_1 %>%
  group_by(days_post_install, behavior) %>%
  summarise(
    mean_count = mean(count, na.rm = TRUE),
    .groups = "drop"
  )

mouse1_kde <- ggplot(mouse1_long_mean,
       aes(x = days_post_install, y = behavior, fill = behavior, weight = mean_count))+
  geom_density_ridges(from = 0, to = 35, bandwidth = 6, scale = 2, alpha = 0.7, panel_scaling = FALSE, rel_min_height = 0.01, color = "white")+
  scale_y_discrete(labels = FALSE, expand = c(0,0))+
  scale_x_continuous(expand = c(0,0))+
  theme_classic()+
  theme(axis.text.y = element_blank(), axis.ticks.y = element_blank())+
  scale_fill_manual(values = c("#D81B60", "#1E88E5", "#FFC107", "#004D40", "#8E24AA"))+
  xlab("Days Post-Installation")+
  ylab("Behavior (frequency distribution)")


## Ensure there are no double-stacked behaviors from days with multiple observations:
mouse2_long_mean <- mousedata_longer_2 %>%
  group_by(days_post_install, behavior) %>%
  summarise(
    mean_count = mean(count, na.rm = TRUE),
    .groups = "drop"
  )

mouse2_kde <- ggplot(mouse2_long_mean,
   aes(x = days_post_install, y = behavior, fill = behavior, weight = mean_count))+
  geom_density_ridges(from = 0, to = 35, bandwidth = 6, scale = 2, alpha = 0.7, panel_scaling = FALSE, rel_min_height = 0.01, color = "white")+
  scale_y_discrete(labels = FALSE, expand = c(0,0))+
  scale_x_continuous(expand = c(0,0))+
  theme_classic()+
  theme(axis.text.y = element_blank(), axis.ticks.y = element_blank())+
  scale_fill_manual(values = c("#D81B60", "#1E88E5", "#FFC107", "#004D40", "#8E24AA"))+
  xlab("Days Post-Installation")+
  ylab("Behavior (frequency distribution")

## Ensure there are no double-stacked behaviors from days with multiple observations:
mouse3_long_mean <- mousedata_longer_3 %>%
  group_by(days_post_install, behavior) %>%
  summarise(
    mean_count = mean(count, na.rm = TRUE),
    .groups = "drop"
  )

mouse3_kde <- ggplot(mouse3_long_mean,
    aes(x = days_post_install, y = behavior, fill = behavior, weight = mean_count))+
  geom_density_ridges(from = 0, to = 35, bandwidth = 6, scale = 2, alpha = 0.7, panel_scaling = FALSE, rel_min_height = 0.01, color = "white")+
  scale_y_discrete(labels = FALSE, expand = c(0,0))+
  scale_x_continuous(expand = c(0,0))+
  theme_classic()+
  theme(axis.text.y = element_blank(), axis.ticks.y = element_blank())+
  scale_fill_manual(values = c("#D81B60", "#1E88E5", "#FFC107", "#004D40", "#8E24AA"))+
  xlab("Days Post-Installation")+
  ylab("Behavior (frequency distribution")

#mouse 4 had some sneaky NAs still there:
mousedata_longer_4 <- mousedata_longer_4 %>%
  filter(!is.na(count))

## Ensure there are no double-stacked behaviors from days with multiple observations:
mouse4_long_mean <- mousedata_longer_4 %>%
  group_by(days_post_install, behavior) %>%
  summarise(
    mean_count = mean(count, na.rm = TRUE),
    .groups = "drop"
  )

mouse4_kde <- ggplot(mouse4_long_mean,
  aes(x = days_post_install, y = behavior, fill = behavior, weight = mean_count))+
  geom_density_ridges(from = 0, to = 35, bandwidth = 6, scale = 2, alpha = 0.7, panel_scaling = FALSE, rel_min_height = 0.01, color = "white")+
  scale_y_discrete(labels = FALSE, expand = c(0,0))+
  scale_x_continuous(expand = c(0,0))+
  theme_classic()+
  theme(axis.text.y = element_blank(), axis.ticks.y = element_blank())+
  scale_fill_manual(values = c("#D81B60", "#1E88E5", "#FFC107", "#004D40", "#8E24AA"))+
  xlab("Days Post-Installation")+
  ylab("Behavior (frequency distribution")


mouse1_kde <- mouse1_kde + theme(axis.title.y = element_blank(),
                                 axis.title.x = element_blank(),
                                 legend.position = "none")

mouse2_kde <- mouse2_kde + theme(axis.title.y = element_blank(),
                                 axis.title.x = element_blank(),
                                 legend.position = "none")

mouse3_kde <- mouse3_kde + theme(axis.title.y = element_blank(),
                                 axis.title.x = element_blank(),
                                 legend.position = "none")

mouse4_kde <- mouse4_kde + theme(axis.title.y = element_blank(),
                                 axis.title.x = element_blank(),
                                 legend.position = "none")

#### bar plots colored by behavior (for Fig 2) ####
## new column:
## mean_behaviors_total <- 
mousedata <- mousedata %>%
  mutate(
    chewing_mean = as.numeric(chewing_mean),
    tugging_mean = as.numeric(tugging_mean),
    sucking_mean = as.numeric(sucking_mean),
    carry_mean = as.numeric(carry_mean),
    cant_tell_mean = as.numeric(cant_tell_mean),
    mean_behav_total = chewing_mean + tugging_mean + sucking_mean + carry_mean + cant_tell_mean
  )
mousedata_longer <- mousedata %>%
  pivot_longer(cols = c(chewing_mean, tugging_mean, sucking_mean, carry_mean, cant_tell_mean),
               names_to = "behavior", 
               values_to = "count")

#Create subsets to easily plot by mouse
mousedata_longer_1 <- subset(mousedata_longer, replicate == 1)
mousedata_longer_2 <- subset(mousedata_longer, replicate == 2)
mousedata_longer_3 <- subset(mousedata_longer, replicate == 3)
mousedata_longer_4 <- subset(mousedata_longer, replicate == 4)
mousedata_longer_5 <- subset(mousedata_longer, replicate == 5)

#include gaps so that x axis is to scale (days)
mouse1gaps <- ggplot(data = mouse1_long_mean %>% filter(days_post_install <= 80), aes(x = days_post_install, y = mean_count, fill = behavior)) +
  geom_col() +  # Use geom_col() instead of geom_bar()
  xlab("Days Post-Installation") +
  ylab("Behavior Count") +
  theme_classic() +
  scale_fill_manual(values = c("#D81B60", "#1E88E5", "#FFC107", "#004D40", "#8E24AA"),  
                    # Custom colors
                    name = "Behavior",  # Change legend title
               labels = c("indiscernible", "carrying", "chewing", "sucking", "tugging")) +
  theme(legend.position = "right",  # Move legend to the bottom
    legend.direction = "vertical",
    plot.title = element_text(vjust = -6, hjust = 0.5, size = 13)) +
  scale_x_continuous(expand = expansion(add = c(0,0)))+
  coord_cartesian(xlim = c(-0.5, 35.5))+
  scale_y_continuous(expand = c(0, 0))# Remove white space below bars
mouse1gaps


mouse2gaps <- ggplot(data = mouse2_long_mean %>% filter(days_post_install <= 80), aes(x = days_post_install, y = mean_count, fill = behavior)) +
  geom_col() +  # Use geom_col() instead of geom_bar()
  xlab("Days Post-Installation") +
  ylab("Behavior Count") +
  theme_classic() +
  scale_fill_manual(values = c("#D81B60", "#1E88E5", "#FFC107", "#004D40", "#8E24AA"),  
                    # Custom colors
                    name = "Behavior",  # Change legend title
                    labels = c("indiscernible", "carrying", "chewing", "sucking", "tugging")) +
  theme(legend.position = "right",  # Move legend to the bottom
    legend.direction = "vertical",
    plot.title = element_text(vjust = -6, hjust = 0.5, size = 13)) +
  scale_x_continuous(expand = expansion(add = c(0,0)))+
  coord_cartesian(xlim = c(-0.5, 35.5))+
  scale_y_continuous(expand = c(0, 0))# Remove white space below bars
mouse2gaps


mouse3gaps <- ggplot(data = mouse3_long_mean %>% filter(days_post_install <= 80), aes(x = days_post_install, y = mean_count, fill = behavior)) +
  geom_col() +  # Use geom_col() instead of geom_bar()
  xlab("Days Post-Installation") +
  ylab("Behavior Count") +
  theme_classic() +
  scale_fill_manual(values = c("#D81B60", "#1E88E5", "#FFC107", "#004D40", "#8E24AA"),  
                    # Custom colors
                    name = "Behavior",  # Change legend title
                    labels = c("indiscernible", "carrying", "chewing", "sucking", "tugging")) +
  theme(legend.position = "right",  # Move legend to the bottom
    legend.direction = "vertical",
    plot.title = element_text(vjust = -6, hjust = 0.5, size = 13)) +
  scale_x_continuous(expand = expansion(add = c(0,0)))+
  coord_cartesian(xlim = c(-0.5, 35.5))+
  scale_y_continuous(expand = c(0, 0))# Remove white space below bars
mouse3gaps

mouse4gaps <- ggplot(data = mouse4_long_mean %>% filter(days_post_install <= 80), aes(x = days_post_install, y = mean_count, fill = behavior)) +
  geom_col() +  # Use geom_col() instead of geom_bar()
  xlab("Days Post-Installation") +
  ylab("Behavior Count") +
  theme_classic() +
  scale_fill_manual(values = c("#D81B60", "#1E88E5", "#FFC107", "#004D40", "#8E24AA"),  
                    # Custom colors
                    name = "Behavior",  # Change legend title
                labels = c("indiscernible", "carrying", "chewing", "sucking", "tugging")) +
  theme(legend.position = "right",  # Move legend to the bottom
    legend.direction = "vertical",
    plot.title = element_text(vjust = -6, hjust = 0.5, size = 13)) +
  scale_x_continuous(expand = expansion(add = c(0,0)))+
  coord_cartesian(xlim = c(-0.5, 35.5))+
  scale_y_continuous(expand = c(0, 0))# Remove white space below bars
mouse4gaps

mouse1gaps <- mouse1gaps + theme(axis.title.y = element_blank(), 
                         axis.title.x = element_blank(),
                         legend.position = "none")

mouse2gaps <- mouse2gaps + theme(axis.title.y = element_blank(),
                         axis.title.x = element_blank(),
                         legend.position = "none")

mouse3gaps <- mouse3gaps + theme(axis.title.y = element_blank(),
                         axis.title.x = element_blank(),
                         legend.position = "none")

mouse4gaps <- mouse4gaps + theme(axis.title.y = element_blank(),
                                 axis.title.x = element_blank(),
                         legend.position = "none")

# Arrange the plots with a common legend at the bottom
combinedfigures_gaps <- ggarrange(mouse1gaps, mouse2gaps, mouse3gaps, mouse4gaps,  
                                  nrow = 4, ncol = 1, 
                                  common.legend = TRUE, legend = "right")

# Add a single y-axis label
combinedfigures_gaps <- annotate_figure(combinedfigures_gaps, 
                                        left = text_grob("Behavior Count (in # bees)", rot = 90, size = 14),
                                        bottom = text_grob("Days Post-Installation", size = 14, hjust = 0.75, 
                                                           vjust = 0.1))

#### Making a combined figure including KDEs and behavior stacks (Fig 2) ####
mouse1gaps <- mouse1gaps + theme(plot.margin = unit(c(14, 0, 14, 21), "pt")) 
# top, right, bottom, left
mouse1_kde <- mouse1_kde + theme(plot.margin = unit(c(14, 21, 14, 28), "pt"))

mouse2gaps <- mouse2gaps + theme(plot.margin = unit(c(14, 0, 14, 21), "pt")) 
# top, right, bottom, left
mouse2_kde <- mouse2_kde + theme(plot.margin = unit(c(14, 21, 14, 28), "pt"))

mouse3gaps <- mouse3gaps + theme(plot.margin = unit(c(14, 0, 14, 21), "pt")) 
# top, right, bottom, left
mouse3_kde <- mouse3_kde + theme(plot.margin = unit(c(14, 21, 14, 28), "pt"))

mouse4gaps <- mouse4gaps + theme(plot.margin = unit(c(14, 0, 14, 21), "pt")) 
# top, right, bottom, left
mouse4_kde <- mouse4_kde + theme(plot.margin = unit(c(14, 21, 14, 28), "pt"))

big_figure_combo <- ggarrange(mouse1gaps, mouse1_kde, mouse2gaps, mouse2_kde,
                              mouse3gaps, mouse3_kde, mouse4gaps, mouse4_kde,
                              nrow = 4, ncol = 2)
#save image at 550:733 ratio


#### Plot bees on mouse by day, with trend line (Fig 4) ####
scatterplot <- ggplot(data = subset(mousedata, !as.character(replicate) %in% c("6", "5")), 
                      aes(y = num_at_start, x = days_post_install)) +
  geom_smooth(method = "loess", color = "black", alpha=0.2, se=TRUE, linewidth = 0.5) +
  geom_point(aes(color = replicate), shape = 16, size = 2.5, alpha = 0.8) +  # Solid color points, no outline
  xlab("Days post-installation") +
  ylab("Bees on mouse at start of observation") +
  theme_classic() +
  theme(axis.text.x = element_text(size = 15), axis.text.y = element_text(size = 15),axis.title.x = element_text(size = 18), axis.title.y = element_text(size = 18), legend.position = "none")+
  scale_color_manual(values = c("#E8000B", "#1AC938", "#FF7C00", "#023EFF"))+
  scale_y_continuous(limits= c(0,50), breaks = c(0,10,20,30,40,50),expand = c(0,.5))+
  scale_x_continuous(breaks = c(5,10,15,20,30,40,50,60,70,80), limits = c(0, 80), expand = expansion(add = c(1, 3)))

scatterplot2 <- scatterplot +  geom_vline(xintercept = 3.5, color = "black", linetype = "solid")

scatterplot3 <- scatterplot2 + geom_vline(xintercept = 10.5, color = "black", linetype = "solid")

scatterplot3 + geom_vline(xintercept=20.5, color = "black", linetype = "solid")
#save at 1200:550 ratio

#### Summaries of observations at various time points, used for in-text stats, see excel file: ####
#### BEHAVIORS DAYS 0-1 ####

# mouse 1:
mousedata %>%
  filter(
    replicate == "1",
    days_post_install >= 0,
    days_post_install <= 1
  ) %>%
  summarise(
    chewing_mean    = sum(chewing_mean, na.rm = TRUE),
    tugging_mean    = sum(tugging_mean, na.rm = TRUE),
    sucking_mean    = sum(sucking_mean, na.rm = TRUE),
    carry_mean      = sum(carry_mean, na.rm = TRUE),
    cant_tell_mean  = sum(cant_tell_mean, na.rm = TRUE),
    cementing_mean  = sum(cementing_mean, na.rm = TRUE)
  )

# Mouse 2:
mousedata %>%
  filter(
    replicate == "2",
    days_post_install >= 0,
    days_post_install <= 1
  ) %>%
  summarise(
    chewing_mean    = sum(chewing_mean, na.rm = TRUE),
    tugging_mean    = sum(tugging_mean, na.rm = TRUE),
    sucking_mean    = sum(sucking_mean, na.rm = TRUE),
    carry_mean      = sum(carry_mean, na.rm = TRUE),
    cant_tell_mean  = sum(cant_tell_mean, na.rm = TRUE),
    cementing_mean  = sum(cementing_mean, na.rm = TRUE)
  )

# Mouse 3:
mousedata %>%
  filter(
    replicate == "3",
    days_post_install >= 0,
    days_post_install <= 1
  ) %>%
  summarise(
    chewing_mean    = sum(chewing_mean, na.rm = TRUE),
    tugging_mean    = sum(tugging_mean, na.rm = TRUE),
    sucking_mean    = sum(sucking_mean, na.rm = TRUE),
    carry_mean      = sum(carry_mean, na.rm = TRUE),
    cant_tell_mean  = sum(cant_tell_mean, na.rm = TRUE),
    cementing_mean  = sum(cementing_mean, na.rm = TRUE)
  )

#Mouse 4:
mousedata %>%
  filter(
    replicate == "4",
    days_post_install >= 0,
    days_post_install <= 1
  ) %>%
  summarise(
    chewing_mean    = sum(chewing_mean, na.rm = TRUE),
    tugging_mean    = sum(tugging_mean, na.rm = TRUE),
    sucking_mean    = sum(sucking_mean, na.rm = TRUE),
    carry_mean      = sum(carry_mean, na.rm = TRUE),
    cant_tell_mean  = sum(cant_tell_mean, na.rm = TRUE),
    cementing_mean  = sum(cementing_mean, na.rm = TRUE)
  )

#### BEHAVIORS DAYS 0-3 ####

# mouse 1:
mousedata %>%
  filter(
    replicate == "1",
    days_post_install >= 0,
    days_post_install <= 3
  ) %>%
  summarise(
    chewing_mean    = sum(chewing_mean, na.rm = TRUE),
    tugging_mean    = sum(tugging_mean, na.rm = TRUE),
    sucking_mean    = sum(sucking_mean, na.rm = TRUE),
    carry_mean      = sum(carry_mean, na.rm = TRUE),
    cant_tell_mean  = sum(cant_tell_mean, na.rm = TRUE),
    cementing_mean  = sum(cementing_mean, na.rm = TRUE)
  )

# Mouse 2:
mousedata %>%
  filter(
    replicate == "2",
    days_post_install >= 0,
    days_post_install <= 3
  ) %>%
  summarise(
    chewing_mean    = sum(chewing_mean, na.rm = TRUE),
    tugging_mean    = sum(tugging_mean, na.rm = TRUE),
    sucking_mean    = sum(sucking_mean, na.rm = TRUE),
    carry_mean      = sum(carry_mean, na.rm = TRUE),
    cant_tell_mean  = sum(cant_tell_mean, na.rm = TRUE),
    cementing_mean  = sum(cementing_mean, na.rm = TRUE)
  )

# Mouse 3:
mousedata %>%
  filter(
    replicate == "3",
    days_post_install >= 0,
    days_post_install <= 3
  ) %>%
  summarise(
    chewing_mean    = sum(chewing_mean, na.rm = TRUE),
    tugging_mean    = sum(tugging_mean, na.rm = TRUE),
    sucking_mean    = sum(sucking_mean, na.rm = TRUE),
    carry_mean      = sum(carry_mean, na.rm = TRUE),
    cant_tell_mean  = sum(cant_tell_mean, na.rm = TRUE),
    cementing_mean  = sum(cementing_mean, na.rm = TRUE)
  )

#Mouse 4:
mousedata %>%
  filter(
    replicate == "4",
    days_post_install >= 0,
    days_post_install <= 3
  ) %>%
  summarise(
    chewing_mean    = sum(chewing_mean, na.rm = TRUE),
    tugging_mean    = sum(tugging_mean, na.rm = TRUE),
    sucking_mean    = sum(sucking_mean, na.rm = TRUE),
    carry_mean      = sum(carry_mean, na.rm = TRUE),
    cant_tell_mean  = sum(cant_tell_mean, na.rm = TRUE),
    cementing_mean  = sum(cementing_mean, na.rm = TRUE)
  )

#### BEHAVIORS Days 0-9 ####
# mouse 1:
mousedata %>%
  filter(
    replicate == "1",
    days_post_install >= 7,
    days_post_install < 10
  ) %>%
  summarise(
    chewing_mean    = sum(chewing_mean, na.rm = TRUE),
    tugging_mean    = sum(tugging_mean, na.rm = TRUE),
    sucking_mean    = sum(sucking_mean, na.rm = TRUE),
    carry_mean      = sum(carry_mean, na.rm = TRUE),
    cant_tell_mean  = sum(cant_tell_mean, na.rm = TRUE),
    cementing_mean  = sum(cementing_mean, na.rm = TRUE)
  )

# mouse 2:
mousedata %>%
  filter(
    replicate == "2",
    days_post_install >= 0,
    days_post_install < 10
  ) %>%
  summarise(
    chewing_mean    = sum(chewing_mean, na.rm = TRUE),
    tugging_mean    = sum(tugging_mean, na.rm = TRUE),
    sucking_mean    = sum(sucking_mean, na.rm = TRUE),
    carry_mean      = sum(carry_mean, na.rm = TRUE),
    cant_tell_mean  = sum(cant_tell_mean, na.rm = TRUE),
    cementing_mean  = sum(cementing_mean, na.rm = TRUE)
  )

# mouse 3:
mousedata %>%
  filter(
    replicate == "3",
    days_post_install >= 0,
    days_post_install < 10
  ) %>%
  summarise(
    chewing_mean    = sum(chewing_mean, na.rm = TRUE),
    tugging_mean    = sum(tugging_mean, na.rm = TRUE),
    sucking_mean    = sum(sucking_mean, na.rm = TRUE),
    carry_mean      = sum(carry_mean, na.rm = TRUE),
    cant_tell_mean  = sum(cant_tell_mean, na.rm = TRUE),
    cementing_mean  = sum(cementing_mean, na.rm = TRUE)
  )

# mouse 4:
mousedata %>%
  filter(
    replicate == "4",
    days_post_install >= 0,
    days_post_install < 10
  ) %>%
  summarise(
    chewing_mean    = sum(chewing_mean, na.rm = TRUE),
    tugging_mean    = sum(tugging_mean, na.rm = TRUE),
    sucking_mean    = sum(sucking_mean, na.rm = TRUE),
    carry_mean      = sum(carry_mean, na.rm = TRUE),
    cant_tell_mean  = sum(cant_tell_mean, na.rm = TRUE),
    cementing_mean  = sum(cementing_mean, na.rm = TRUE)
  )
#### BEHAVIORS DAYS 4-10 ####

# mouse 1:
mousedata %>%
  filter(
    replicate == "1",
    days_post_install >= 4,
    days_post_install <= 10
  ) %>%
  summarise(
    chewing_mean    = sum(chewing_mean, na.rm = TRUE),
    tugging_mean    = sum(tugging_mean, na.rm = TRUE),
    sucking_mean    = sum(sucking_mean, na.rm = TRUE),
    carry_mean      = sum(carry_mean, na.rm = TRUE),
    cant_tell_mean  = sum(cant_tell_mean, na.rm = TRUE),
    cementing_mean  = sum(cementing_mean, na.rm = TRUE)
  )

# Mouse 2:
mousedata %>%
  filter(
    replicate == "2",
    days_post_install >= 4,
    days_post_install <= 10
  ) %>%
  summarise(
    chewing_mean    = sum(chewing_mean, na.rm = TRUE),
    tugging_mean    = sum(tugging_mean, na.rm = TRUE),
    sucking_mean    = sum(sucking_mean, na.rm = TRUE),
    carry_mean      = sum(carry_mean, na.rm = TRUE),
    cant_tell_mean  = sum(cant_tell_mean, na.rm = TRUE),
    cementing_mean  = sum(cementing_mean, na.rm = TRUE)
  )

# Mouse 3:
mousedata %>%
  filter(
    replicate == "3",
    days_post_install >= 4,
    days_post_install <= 10
  ) %>%
  summarise(
    chewing_mean    = sum(chewing_mean, na.rm = TRUE),
    tugging_mean    = sum(tugging_mean, na.rm = TRUE),
    sucking_mean    = sum(sucking_mean, na.rm = TRUE),
    carry_mean      = sum(carry_mean, na.rm = TRUE),
    cant_tell_mean  = sum(cant_tell_mean, na.rm = TRUE),
    cementing_mean  = sum(cementing_mean, na.rm = TRUE)
  )

#Mouse 4:
mousedata %>%
  filter(
    replicate == "4",
    days_post_install >= 4,
    days_post_install <= 10
  ) %>%
  summarise(
    chewing_mean    = sum(chewing_mean, na.rm = TRUE),
    tugging_mean    = sum(tugging_mean, na.rm = TRUE),
    sucking_mean    = sum(sucking_mean, na.rm = TRUE),
    carry_mean      = sum(carry_mean, na.rm = TRUE),
    cant_tell_mean  = sum(cant_tell_mean, na.rm = TRUE),
    cementing_mean  = sum(cementing_mean, na.rm = TRUE)
  )

#### BEHAVIORS DAYS 10+ ####
# mouse 1:
mousedata %>%
  filter(
    replicate == "1",
    days_post_install >= 10
  ) %>%
  summarise(
    chewing_mean    = sum(chewing_mean, na.rm = TRUE),
    tugging_mean    = sum(tugging_mean, na.rm = TRUE),
    sucking_mean    = sum(sucking_mean, na.rm = TRUE),
    carry_mean      = sum(carry_mean, na.rm = TRUE),
    cant_tell_mean  = sum(cant_tell_mean, na.rm = TRUE),
    cementing_mean  = sum(cementing_mean, na.rm = TRUE)
  )

# mouse 2: 
mousedata %>%
  filter(
    replicate == "2",
    days_post_install >= 10
  ) %>%
  summarise(
    chewing_mean    = sum(chewing_mean, na.rm = TRUE),
    tugging_mean    = sum(tugging_mean, na.rm = TRUE),
    sucking_mean    = sum(sucking_mean, na.rm = TRUE),
    carry_mean      = sum(carry_mean, na.rm = TRUE),
    cant_tell_mean  = sum(cant_tell_mean, na.rm = TRUE),
    cementing_mean  = sum(cementing_mean, na.rm = TRUE)
  )

# mouse 3:
mousedata %>%
  filter(
    replicate == "3",
    days_post_install >= 10
  ) %>%
  summarise(
    chewing_mean    = sum(chewing_mean, na.rm = TRUE),
    tugging_mean    = sum(tugging_mean, na.rm = TRUE),
    sucking_mean    = sum(sucking_mean, na.rm = TRUE),
    carry_mean      = sum(carry_mean, na.rm = TRUE),
    cant_tell_mean  = sum(cant_tell_mean, na.rm = TRUE),
    cementing_mean  = sum(cementing_mean, na.rm = TRUE)
  )

#mouse 4: 
mousedata %>%
  filter(
    replicate == "4",
    days_post_install >= 10
  ) %>%
  summarise(
    chewing_mean    = sum(chewing_mean, na.rm = TRUE),
    tugging_mean    = sum(tugging_mean, na.rm = TRUE),
    sucking_mean    = sum(sucking_mean, na.rm = TRUE),
    carry_mean      = sum(carry_mean, na.rm = TRUE),
    cant_tell_mean  = sum(cant_tell_mean, na.rm = TRUE),
    cementing_mean  = sum(cementing_mean, na.rm = TRUE)
  )

#### BEHAVIORS DAYS 10 - 80 ####
# mouse 1:
mousedata %>%
  filter(
    replicate == "1",
    days_post_install >= 10,
    days_post_install <= 80
  ) %>%
  summarise(
    chewing_mean    = sum(chewing_mean, na.rm = TRUE),
    tugging_mean    = sum(tugging_mean, na.rm = TRUE),
    sucking_mean    = sum(sucking_mean, na.rm = TRUE),
    carry_mean      = sum(carry_mean, na.rm = TRUE),
    cant_tell_mean  = sum(cant_tell_mean, na.rm = TRUE),
    cementing_mean  = sum(cementing_mean, na.rm = TRUE)
  )

# mouse 2: 
mousedata %>%
  filter(
    replicate == "2",
    days_post_install >= 10,
    days_post_install <= 80
  ) %>%
  summarise(
    chewing_mean    = sum(chewing_mean, na.rm = TRUE),
    tugging_mean    = sum(tugging_mean, na.rm = TRUE),
    sucking_mean    = sum(sucking_mean, na.rm = TRUE),
    carry_mean      = sum(carry_mean, na.rm = TRUE),
    cant_tell_mean  = sum(cant_tell_mean, na.rm = TRUE),
    cementing_mean  = sum(cementing_mean, na.rm = TRUE)
  )

# mouse 3:
mousedata %>%
  filter(
    replicate == "3",
    days_post_install >= 10,
    days_post_install <= 80
  ) %>%
  summarise(
    chewing_mean    = sum(chewing_mean, na.rm = TRUE),
    tugging_mean    = sum(tugging_mean, na.rm = TRUE),
    sucking_mean    = sum(sucking_mean, na.rm = TRUE),
    carry_mean      = sum(carry_mean, na.rm = TRUE),
    cant_tell_mean  = sum(cant_tell_mean, na.rm = TRUE),
    cementing_mean  = sum(cementing_mean, na.rm = TRUE)
  )

#mouse 4: 
mousedata %>%
  filter(
    replicate == "4",
    days_post_install >= 10,
    days_post_install <= 80
  ) %>%
  summarise(
    chewing_mean    = sum(chewing_mean, na.rm = TRUE),
    tugging_mean    = sum(tugging_mean, na.rm = TRUE),
    sucking_mean    = sum(sucking_mean, na.rm = TRUE),
    carry_mean      = sum(carry_mean, na.rm = TRUE),
    cant_tell_mean  = sum(cant_tell_mean, na.rm = TRUE),
    cementing_mean  = sum(cementing_mean, na.rm = TRUE)
  )

#### BEHAVIORS ALL TIME ####
# mouse 1:
mousedata %>%
  filter(
    replicate == "1",
  ) %>%
  summarise(
    chewing_mean    = sum(chewing_mean, na.rm = TRUE),
    tugging_mean    = sum(tugging_mean, na.rm = TRUE),
    sucking_mean    = sum(sucking_mean, na.rm = TRUE),
    carry_mean      = sum(carry_mean, na.rm = TRUE),
    cant_tell_mean  = sum(cant_tell_mean, na.rm = TRUE),
    cementing_mean  = sum(cementing_mean, na.rm = TRUE)
  )

# Mouse 2:
mousedata %>%
  filter(
    replicate == "2"
  ) %>%
  summarise(
    chewing_mean    = sum(chewing_mean, na.rm = TRUE),
    tugging_mean    = sum(tugging_mean, na.rm = TRUE),
    sucking_mean    = sum(sucking_mean, na.rm = TRUE),
    carry_mean      = sum(carry_mean, na.rm = TRUE),
    cant_tell_mean  = sum(cant_tell_mean, na.rm = TRUE),
    cementing_mean  = sum(cementing_mean, na.rm = TRUE)
  )

# Mouse 3:
mousedata %>%
  filter(
    replicate == "3"
  ) %>%
  summarise(
    chewing_mean    = sum(chewing_mean, na.rm = TRUE),
    tugging_mean    = sum(tugging_mean, na.rm = TRUE),
    sucking_mean    = sum(sucking_mean, na.rm = TRUE),
    carry_mean      = sum(carry_mean, na.rm = TRUE),
    cant_tell_mean  = sum(cant_tell_mean, na.rm = TRUE),
    cementing_mean  = sum(cementing_mean, na.rm = TRUE)
  )

#Mouse 4:
mousedata %>%
  filter(
    replicate == "4"
  ) %>%
  summarise(
    chewing_mean    = sum(chewing_mean, na.rm = TRUE),
    tugging_mean    = sum(tugging_mean, na.rm = TRUE),
    sucking_mean    = sum(sucking_mean, na.rm = TRUE),
    carry_mean      = sum(carry_mean, na.rm = TRUE),
    cant_tell_mean  = sum(cant_tell_mean, na.rm = TRUE),
    cementing_mean  = sum(cementing_mean, na.rm = TRUE)
  )

#### BEHAVIORS DAYS 1-80 ####
# mouse 1:
mousedata %>%
  filter(
    replicate == "1",
    days_post_install <= 80
  ) %>%
  summarise(
    chewing_mean    = sum(chewing_mean, na.rm = TRUE),
    tugging_mean    = sum(tugging_mean, na.rm = TRUE),
    sucking_mean    = sum(sucking_mean, na.rm = TRUE),
    carry_mean      = sum(carry_mean, na.rm = TRUE),
    cant_tell_mean  = sum(cant_tell_mean, na.rm = TRUE),
    cementing_mean  = sum(cementing_mean, na.rm = TRUE)
  )

# Mouse 2:
mousedata %>%
  filter(
    replicate == "2",
    days_post_install <= 80
  ) %>%
  summarise(
    chewing_mean    = sum(chewing_mean, na.rm = TRUE),
    tugging_mean    = sum(tugging_mean, na.rm = TRUE),
    sucking_mean    = sum(sucking_mean, na.rm = TRUE),
    carry_mean      = sum(carry_mean, na.rm = TRUE),
    cant_tell_mean  = sum(cant_tell_mean, na.rm = TRUE),
    cementing_mean  = sum(cementing_mean, na.rm = TRUE)
  )

# Mouse 3:
mousedata %>%
  filter(
    replicate == "3",
    days_post_install <= 80
  ) %>%
  summarise(
    chewing_mean    = sum(chewing_mean, na.rm = TRUE),
    tugging_mean    = sum(tugging_mean, na.rm = TRUE),
    sucking_mean    = sum(sucking_mean, na.rm = TRUE),
    carry_mean      = sum(carry_mean, na.rm = TRUE),
    cant_tell_mean  = sum(cant_tell_mean, na.rm = TRUE),
    cementing_mean  = sum(cementing_mean, na.rm = TRUE)
  )

#Mouse 4:
mousedata %>%
  filter(
    replicate == "4",
    days_post_install <= 80
  ) %>%
  summarise(
    chewing_mean    = sum(chewing_mean, na.rm = TRUE),
    tugging_mean    = sum(tugging_mean, na.rm = TRUE),
    sucking_mean    = sum(sucking_mean, na.rm = TRUE),
    carry_mean      = sum(carry_mean, na.rm = TRUE),
    cant_tell_mean  = sum(cant_tell_mean, na.rm = TRUE),
    cementing_mean  = sum(cementing_mean, na.rm = TRUE)
  )

#### BEES ON MOUSE COUNTS ####
# Days 0-3
mousedata %>%
  filter(
    days_post_install <= 3
  ) %>%
  summarise(
    mean_num_at_start    = mean(num_at_start, na.rm = TRUE),
    sd_num_at_start    = sd(num_at_start, na.rm = TRUE),
    range_num_at_start = diff(range(num_at_start, na.rm = TRUE))
  )

# Days 4-10
mousedata %>%
  filter(
    days_post_install >= 4,
    days_post_install <= 10
  ) %>%
  summarise(
    mean_num_at_start    = mean(num_at_start, na.rm = TRUE),
    sd_num_at_start    = sd(num_at_start, na.rm = TRUE),
    range_num_at_start = diff(range(num_at_start, na.rm = TRUE))
  )

# Days 11-20
mousedata %>%
  filter(
    days_post_install >= 11,
    days_post_install <= 20
  ) %>%
  summarise(
    mean_num_at_start    = mean(num_at_start, na.rm = TRUE),
    sd_num_at_start    = sd(num_at_start, na.rm = TRUE),
    range_num_at_start = diff(range(num_at_start, na.rm = TRUE))
  )

# Days 21-80
mousedata %>%
  filter(
    days_post_install >= 21,
    days_post_install <= 80
  ) %>%
  summarise(
    mean_num_at_start    = mean(num_at_start, na.rm = TRUE),
    sd_num_at_start    = sd(num_at_start, na.rm = TRUE),
    range_num_at_start = diff(range(num_at_start, na.rm = TRUE))
  )

#### Focal area summaries used for figure 4 ####
#make a dataset without mouse replicates 5 & 6
simp_longdata_no56 <- simp_longdata[ !(simp_longdata$replicate %in% c("5", "6")), ]
table(simp_longdata_no56$focal_area)

#summarizing counts for the purpose of figure making outside of R
table(simp_longdata_no56$focal_area[simp_longdata_no56$days_post_install < 7])
table(simp_longdata_no56$focal_area[simp_longdata_no56$days_post_install < 4])
table(simp_longdata_no56$focal_area[simp_longdata_no56$days_post_install < 3])
table(simp_longdata_no56$focal_area[simp_longdata_no56$days_post_install < 2])
table(simp_longdata_no56$focal_area[simp_longdata_no56$days_post_install > 10 & simp_longdata_no56$days_post_install < 21])
table(simp_longdata_no56$focal_area[simp_longdata_no56$days_post_install > 20 & simp_longdata_no56$days_post_install < 81])
table(simp_longdata_no56$focal_area)