SW_Respirometry_standard_markdown_DEBUGGED.Rmd

---
title: "SW (NO E. coli) Standard  Respirometry Data Analysis"
author: "Caroline Fleming Ianniello (adapted from S. Speroff with help from Brian Kennedy and Ethan Deyle)"
date: "last updated 12/12/2023"
output: html_document
---


# Some important notes for working with this code: 

  ##Remember that negative control with sodium nitrate (Oxycal, https://www.pyroscience.com/en/products/accessories/calibration-capsules/oxcal)  always has to be A1 vial, and that order goes A1-D1, etc.
  ## LINES TO ADAPT FOR EACH RUN: 
  ### Input data file straight from PreSens software
  ### Number of vials filled (SEE NOTE ABOVE)
  ### CSV output name
  
  ## R hates 0s, so we may want to replace with something small in the weights column


```{r set working directory}
# #For some reason knitr is wanting me to set the WD this way
# require("knitr")
# opts_knit$set(root.dir = "C:/Users/carol/Documents/PhD Boston University 2019-/Astrangia_nutrient_stressors_Rotjan_lab/Surface_area_resp_reruns")
```

```{r install packages}

#installing new respR
#install.packages("respR")
# library necessary packages
#install.packages("knitr")
library("knitr")
library(respR)
library(tidyverse)
library(tibble)
#readr is installed in tidyverse
```


```{r EDIT ME - IMPORT DATA}
# raw data import
#the following is a command from readr
data <- read_csv("EXP1A_1_07_26_22_NO_e_coli_Oxygen.csv", skip=12, na=c("No Sensor","NA"))


```


```{r Import & clean CSV}
# rename headers
names(data)[names(data) == "Time/Min."] <- "Minutes" 

# remove unnecessary columns (keep only data/time, minutes, and data)
data <- data[,1:26]

# new time elapsed (sec) interval column
for (i in 1:nrow(data)){ 
  x <- nrow(data)
  s <- seq(from = 0, by = 15, length.out = x)
  data$Time.Elapsed.Sec <- (rep(s))
}

# new time elapsed (min) interval column 
for (i in 1:nrow(data)){
  data$Time.Elapsed.Min <- data$Time.Elapsed.Sec/60
}

#Rename Date/Time, Date.Time
names(data)[1] <-"Date.Time"
# final data check if you want
#str(data)


#make sure all are numeric

data$A1 <- as.numeric(data$A1)
data$A2 <- as.numeric(data$A2)
data$A3 <- as.numeric(data$A3)
data$A4 <- as.numeric(data$A4)
data$A5 <- as.numeric(data$A5)
data$A6 <- as.numeric(data$A6)
data$B1 <- as.numeric(data$B1)
data$B2 <- as.numeric(data$B2)
data$B3 <- as.numeric(data$B3)
data$B4 <- as.numeric(data$B4)
data$B5 <- as.numeric(data$B5)
data$B6 <- as.numeric(data$B6)
data$C1 <- as.numeric(data$C1)
data$C2 <- as.numeric(data$C2)
data$C3 <- as.numeric(data$C3)
data$C4 <- as.numeric(data$C4)
data$C5 <- as.numeric(data$C5)
data$C6 <- as.numeric(data$C6)
data$D1 <- as.numeric(data$D1)
data$D2 <- as.numeric(data$D2)
data$D3 <- as.numeric(data$D3)
data$D4 <- as.numeric(data$D4)
data$D5 <- as.numeric(data$D5)
data$D6 <- as.numeric(data$D6)


```

```{r cut off everything past four and a half hours}
data<-subset(data,data$Minutes<=270)
#270 minutes is four and a half hours. We are manually cutting off all experiments after four and a half hours.
```

```{r Scale negative O2 control and re-zero baseline}

# calculate the mean value of our negative control in a new df
#selecting columns with data and then subtract the average respiration rate from column with negative control
#note: na.rm=TRUE tells R to ignore NAs
neg <- data[,3:26] - mean(data$A1, na.rm=TRUE)

# create new columns in neg to match the time values in data
neg <- neg %>%
  add_column(Minutes = data$Minutes,
             .before = "A1")
neg <- neg %>%
  add_column(Date.Time = data$Date.Time,
             .before = "Minutes")
neg <- neg %>%
  add_column(Time.Elapsed.Sec = data$Time.Elapsed.Sec,
             .after = "D6") 
neg <- neg %>%
  add_column(Time.Elapsed.Min = data$Time.Elapsed.Min,
             .after = "Time.Elapsed.Sec")

# change name back to data
data <- neg


```

```{r EDIT ME - TEMPLATE UPLOAD FOR SURFACE AREA}

resp_info <- read.csv("SA_template_resp_EXP1A_1_07_26_22_NO_e_coli_Oxygen.csv", header=TRUE) # NUMBER ORDER, NOT LETTER ORDER
  volume= as.numeric(resp_info$volume) #this will always be 0.0019L so we can just hard code it
  t= as.numeric(resp_info$temp)
  S= as.numeric(resp_info$salinity)
  P= as.numeric(resp_info$pressure)
```


```{r num vials filled}
#RUN BASED ON NUMBERS NOT LETTERS! (NEED IT IN ORDER)
#enter number of vials filled in experiment
num_vials_filled <- 24

```

```{r generate background respiration rate, SW}

#background SW rate of SW only
# Make a SW control dataframe using columns 4 (our SW positive control) & the time
SW.data <- na.omit(data[c(28, 4)]) #MAKE SURE THIS IS THE MINUTES COLUMN


# SW inspect object, with all oxygen values pooled
SW_insp <- inspect(SW.data,
                   time = 1, 
                   oxygen = 2, 
                   plot = TRUE)

# check assumptions of inspect object 
print(SW_insp)

# Calc background rate of sw subset, this is the big calculation from the RespR package
SW.rate <- calc_rate.bg(SW.data, 
                        time = 1,
                        oxygen = 2)

# print calculated rates & average rate
print(SW.rate)



```

```{r Negative oxygen control}
# Make a -O2 control data_frame using column 3, which will ALWAYS be the negative control
NEG.O2.data<- na.omit(data[c(28, 3)])

# Inspect 
NEG_O2_insp <- inspect(NEG.O2.data,
                       time = 1, 
                       oxygen = 2,
                       plot = TRUE)
# Print
print(NEG_O2_insp)

# Calc background rate of negative O2 control
NEG.O2.rate <- calc_rate.bg(NEG.O2.data, 
                            time = 1,
                            oxygen = 2)

# print calculated rates & average rate
print(NEG.O2.rate)

#these are really just checks
```

```{r For loop}

# create empty output vector to fill in in for loop
output_area_specific_O2 <- rep("NA", num_vials_filled)
# selecting columns with data based on number of vials filled
data2 <- data[, c(1:2, 3:(3 + num_vials_filled), 27:28)] #update last columns for a full run

#BIG FOR LOOP
#NOTE: this for loop starts at B1, what should be the positive oxygen control, because we don't care about A1 the negative control
i=4
for (i in 2:num_vials_filled) { 
  temp.data <- data2[, c(ncol(data2), i+3)]
  #  temp.data <- data2[c(ncol(data2), i+3)]
 # temp.data <- na.omit(temp.data) #get rid of NAs
  #checking inspect file for any flags 
  #I want the time.elapsed.min rather than minutes column to be selected
  temp_inspect <- inspect(temp.data,
                   time = 1, #time column 1
                   oxygen = 2, #oxygen column 2
                   plot = FALSE)
  print(paste("vial", i, "NA Time Elapsed", temp_inspect$checks[1])) # FALSE is good
  print(paste("vial", i, "sequential", temp_inspect$checks[2]))
  print(paste("vial", i, "duplicated", temp_inspect$checks[3]))
  print(paste("vial", i, "evenly-spaced", temp_inspect$checks[4]))
  print(paste("vial", i, "NA Oxygen", temp_inspect$checks[5]))
  if (temp_inspect$checks[1]== TRUE) {
    print(paste("ALERT", "vial", i, "NA Time Elapsed FAILED"))
  }
  if (temp_inspect$checks[2]== TRUE)  {
    print(paste("ALERT", "vial", i, "sequential FAILED"))
    
  }
  if (temp_inspect$checks[3]==TRUE) {
    print(paste("ALERT", "vial", i, "duplicated FAILED"))
    
  }
  if (temp_inspect$checks[4]==TRUE) {
    print(paste("ALERT", "vial", i, "evenly-spaced FAILED"))
    
  }
  if (temp_inspect$checks[5]==TRUE) {
    print(paste("ALERT", "vial", i, "NA Oxygen FAILED"))
    
  }
  
  temp_rate <- calc_rate(temp_inspect, 
                     from = 1,
                     to = nrow(temp.data),
                     by = "row")
#now adjusting the rate to SW
  temp_adjusted_rate <- adjust_rate(temp_rate, SW.rate) #adjust to background seawater respiration
  
  area_animal = as.numeric(resp_info[1, i + 8]) #selecting 8 because that is the first with real data -- first vial column will always be negative control
  
  MO2_temp <- convert_rate(temp_adjusted_rate,
                       oxy.unit = "%Air",            # original O2 units
                       time.unit = "min",        # time units 
                       output.unit = "umol/hr/cm2",     # desired rate units 
                       volume = 0.0019,             # chamber liquid volume, MUST BE IN LITERS, will always be 1.9mL
                       area = convert_val(area_animal, from="cm2"),           # surface area MUST BE in meters squared (take cm squared / 10,000) using convert_val function to get m2
                       t = t,                 # temperature
                       S = S,                   # salinity 
                       P = P)               # atmospheric pressure

  #now we extract the area specific rates and put them in the empty vector we made
  output_area_specific_O2[i] <- as.numeric(MO2_temp$summary$rate.a.spec)

}
# t
# S
# P
# area_animal
# 
#  output_area_specific_O2[6]
```



```{r EDIT ME - OUTPUT FILE}
#taking output from output_area_specific_O2 and putting it into a CSV by vial name
names <- c("B1", "C1", "D1", "A2", "B2", "C2", "D2", "A3", "B3", "C3", "D3" ,"A4" ,"B4", "C4" ,"D4", "A5", "B5", "C5" ,"D5" ,"A6", "B6", "C6", "D6")
names_temp <- names[1: num_vials_filled]
output_temp<-rbind(names_temp,output_area_specific_O2)
write.csv(output_temp, "OUTPUT_SA_EXP1A_1_07_26_22_NO_e_coli_Oxygen.csv", row.names=FALSE)


```