Markdown script for finding the Lethal concentration of B.b spores

#library(ggthemes)
#library(ecotox)
#library(tidyverse)
#library(dplyr)
#library(readxl)
#library(ggplot2)

LC for Beauveria

Import dataset
Probit and Logit models both produce predictions of probabilities that lie inside the interval [0,1]
1) Doses need to be in spores/ml format
2) Data is coming directly from an uploaded excel file with columns: colony, dose (spores/ml), total (flies/per cage), and response (# dead)
3) Each cage should be its own row in the excel file.
dip <- read_excel("dip_11_24_24_R.xlsx", 
    sheet = "Sheet2")

##remove dose "zero" from data #This function only accepts non-zero variables to determine proportion response (probit link)
dip <- dip %>%
  filter(dose != "0")

##using ecotox data_analysis tutorial https://search.r-project.org/CRAN/refmans/ecotox/html/LC_probit.html

#This is for testing many different colonies at once. You can use this with one colony though without modification. 
colony_names <- unique(dip$colony) #pull out colony names into a new list 
all_models <- list() #create an empty list to store all of the models 

#Loop to find LC50 and LC90 
for (colony_name in colony_names) {
  colony_data <- subset(dip, colony == colony_name)  # Subset data for the current colony
  
  # Fit the LC_probit model
  modeldip <- LC_probit(response / total ~ log10(dose), #total = total flies per cage   #response = # dead after 14 days 
                     p = c(50, 90), #LC 50 and LC 90 
                     weights = total, 
                     data = colony_data)
 modeldip$"dose(spores/ml)" <- modeldip$dose
 modeldip$colony <- colony_name 
 modeldip <- modeldip[, c(21, 1, 20, 4:8, 11)] #Tells what columns to save
 model_name <- paste0(colony_name, ".modeldip")
 assign(model_name, modeldip)
 all_models[[colony_name]] <- modeldip
}
compiled_modelsdip <- do.call(rbind, all_models)

write.csv(compiled_modelsdip, "dipfall2024LC") #Create a new excel file with the results. Reminder that the dose results will be in spores/ml

Make a lethal concentration curve

ggplot(data = dip, aes(x=log10(dose),y=(response/total), color = colony)) + #You can include "colony" for one colony and this still works 
  geom_point() + 
  geom_smooth(method = "glm", method.args = list(family = binomial(link = "probit")), 
              aes(weight = total), se = TRUE) +
  labs(x = "Log (dose(spores/ml))", y = "Proportion mortality", title = "Lethal Concentration Curve", subtitle = "desired subtitle")+ scale_y_continuous(labels = scales::percent_format(scale = 100)) +
  theme_clean()
## `geom_smooth()` using formula = 'y ~ x'