This is an R Markdown Notebook. When you execute code within the notebook, the results appear beneath the code.

Try executing this chunk by clicking the Run button within the chunk or by placing your cursor inside it and pressing Ctrl+Shift+Enter.

getwd()
## [1] "C:/Users/J Lab/Box/Plant Architecture Lab/Hayley/Cowpea/Cowpea_screen_01_March2022/Cowpea_Phenotype_exp1_analysis"
list.files(pattern=".csv")
## [1] "EVT.csv"                               
## [2] "EVT_exp1.csv"                          
## [3] "FW_DW_watercontent.csv"                
## [4] "multispeq_data_exp1.csv"               
## [5] "Pot_Geno.csv"                          
## [6] "watering 3_12-3_17 - Watering Data.csv"
pot_geno <- read.csv("Pot_Geno.csv")
colnames(pot_geno)[1] <- "Pot.Number"
pot_geno
photo_1exp <- read.csv("multispeq_data_exp1.csv")
photo_1exp
colnames(photo_1exp)
##  [1] "ID"                            "Series"                       
##  [3] "Repeat"                        "Pot.Number"                   
##  [5] "Treatment"                     "air_temp_kinetics"            
##  [7] "Ambient.Humidity"              "Ambient.Pressure"             
##  [9] "Ambient.Temperature"           "data_raw_PAM"                 
## [11] "ECS_averaged_trace"            "ECS_tau"                      
## [13] "ECSt.mAU"                      "fitinput"                     
## [15] "FmPrime"                       "FoPrime"                      
## [17] "Fs"                            "FvP_over_FmP"                 
## [19] "gH."                           "humidity2_K"                  
## [21] "humidity_K"                    "kP700"                        
## [23] "leaf.angle"                    "Leaf.Temperature"             
## [25] "Leaf.Temperature.Differenial"  "Leaf.Temperature.Differential"
## [27] "LEAF_temp"                     "leaf_thickness"               
## [29] "LEF"                           "LEFd_trace"                   
## [31] "Light.Intensity..PAR."         "NPQt"                         
## [33] "outdata"                       "P700_DIRK_ampl"               
## [35] "P700_DIRK_averaged_trace"      "P700_fitinput"                
## [37] "P700_outdata"                  "Phi2"                         
## [39] "PhiNO"                         "PhiNPQ"                       
## [41] "PS1.Active.Centers"            "PS1.Open.Centers"             
## [43] "PS1.Over.Reduced.Centers"      "PS1.Oxidized.Centers"         
## [45] "PSI_data_absorbance"           "pump"                         
## [47] "qL"                            "SPAD"                         
## [49] "test_data_raw_PAM"             "time"                         
## [51] "Time.of.Day"                   "tP700"                        
## [53] "v_initial_P700"                "vH."                          
## [55] "User"                          "Device.ID"                    
## [57] "Latitude"                      "Longitude"                    
## [59] "Issues"
photo_interest <- photo_1exp[,c(4, 5, 50, 15:16, 18, 24, 28, 40:44, 48)]
photo_interest
photo_geno_1exp <- merge(photo_interest, pot_geno, all = TRUE)
photo_geno_1exp <- na.omit(photo_geno_1exp)
photo_geno_1exp[order(photo_geno_1exp$time),]
photo_geno_1exp$time <- as.factor(photo_geno_1exp$time)
photo_geno_1exp$day.of.stress <- photo_geno_1exp$time
photo_geno_1exp$day.of.stress <- gsub("3/16/2022 0:00", "14", photo_geno_1exp$day.of.stress)
photo_geno_1exp$day.of.stress <- gsub("3/9/2022 0:00", "7", photo_geno_1exp$day.of.stress)
colnames(photo_geno_1exp)
##  [1] "Pot.Number"               "Treatment"               
##  [3] "time"                     "FmPrime"                 
##  [5] "FoPrime"                  "FvP_over_FmP"            
##  [7] "Leaf.Temperature"         "leaf_thickness"          
##  [9] "PhiNPQ"                   "PS1.Active.Centers"      
## [11] "PS1.Open.Centers"         "PS1.Over.Reduced.Centers"
## [13] "PS1.Oxidized.Centers"     "SPAD"                    
## [15] "Genotype"                 "day.of.stress"
photo_geno_1exp <- photo_geno_1exp[,c(1:2, 16, 4:15)]
#When examining the data, pot 59 on day 14 measurements for PS1 open, over reduced and oxidized centers is "null" and was not recorded successfully. This pot on this day will be removed. 
photo_geno_1expb <- subset(photo_geno_1exp, photo_geno_1exp$PS1.Open.Centers != "null")
photo_geno_1expb$PS1.Open.Centers <- as.numeric(photo_geno_1expb$PS1.Open.Centers)
photo_geno_1expb$PS1.Over.Reduced.Centers <- as.numeric(photo_geno_1expb$PS1.Over.Reduced.Centers)
photo_geno_1expb$PS1.Oxidized.Centers <- as.numeric(photo_geno_1expb$PS1.Oxidized.Centers)
photo_geno_1expb
#write.csv(photo_geno_1expb, "photosynQ_geno_exp1.csv", row.names = TRUE)
library("ggplot2")
library("cowplot")
library("ggpubr")
## 
## Attaching package: 'ggpubr'
## The following object is masked from 'package:cowplot':
## 
##     get_legend
library("reshape2")

###Histograms

FoPrime_exp1_h <- gghistogram(photo_geno_1expb, x = "FoPrime", bins = 30, binwidth = 20,
                            add = "mean", rug = TRUE,
                            color = "Treatment", fill = "Treatment", facet.by = "day.of.stress",
                            palette = c("blue", "orange")) + xlab("FoPrime (a.u)")
FoPrime_exp1_h

FmPrime_exp1_h <- gghistogram(photo_geno_1expb, x = "FmPrime", bins = 30, binwidth = 100,
                            add = "mean", rug = TRUE,
                            color = "Treatment", fill = "Treatment", facet.by = "day.of.stress",
                            palette = c("blue", "orange")) + xlab("FmPrime (a.u)")
FmPrime_exp1_h

FvP_over_FmP_exp1_h <- gghistogram(photo_geno_1expb, x = "FvP_over_FmP", binwidth = 0.01,
                            add = "mean", rug = TRUE,
                            color = "Treatment", fill = "Treatment", facet.by = "day.of.stress",
                            palette = c("blue", "orange")) + xlab("FvP over FmP (a.u)")
FvP_over_FmP_exp1_h

leaftemp_exp1_h <- gghistogram(photo_geno_1expb, x = "Leaf.Temperature", binwidth = 0.5, 
                            add = "mean", rug = TRUE,
                            color = "Treatment", fill = "Treatment", facet.by = "day.of.stress",
                            palette = c("blue", "orange")) + xlab("Leaf Temperature (C)")
leaftemp_exp1_h

leafthickness_exp1_h <- gghistogram(photo_geno_1expb, x = "leaf_thickness", binwidth = 0.05,
                            add = "mean", rug = TRUE,
                            color = "Treatment", fill = "Treatment", facet.by = "day.of.stress",
                            palette = c("blue", "orange")) + xlab("Leaf Thickness (mm)")
leafthickness_exp1_h

PhiNPQ_exp1_h <- gghistogram(photo_geno_1expb, x = "PhiNPQ", binwidth = 0.01,
                            add = "mean", rug = TRUE,
                            color = "Treatment", fill = "Treatment", facet.by = "day.of.stress",
                            palette = c("blue", "orange")) + xlab("PhiNPQ (a.u.)")
PhiNPQ_exp1_h

PS1.Active.Centers_exp1_h <- gghistogram(photo_geno_1expb, x = "PS1.Active.Centers", binwidth = 1,
                            add = "mean", rug = TRUE,
                            color = "Treatment", fill = "Treatment", facet.by = "day.of.stress",
                            palette = c("blue", "orange")) + xlab("PS1.Active.Centers (a.u.)")
PS1.Active.Centers_exp1_h

PS1.Open.Centers_exp1_h <- gghistogram(photo_geno_1expb, x = "PS1.Open.Centers",
                            add = "mean", rug = TRUE,
                            color = "Treatment", fill = "Treatment", facet.by = "day.of.stress",
                            palette = c("blue", "orange")) + xlab("PS1.Open.Centers (a.u.)")
## Warning: Using `bins = 30` by default. Pick better value with the argument
## `bins`.
PS1.Open.Centers_exp1_h

PS1.Over.Reduced.Centers_exp1_h <- gghistogram(photo_geno_1expb, x = "PS1.Over.Reduced.Centers",
                            add = "mean", rug = TRUE,
                            color = "Treatment", fill = "Treatment", facet.by = "day.of.stress",
                            palette = c("blue", "orange")) + xlab("PS1.Over.Reduced.Centers (a.u.)")
## Warning: Using `bins = 30` by default. Pick better value with the argument
## `bins`.
PS1.Over.Reduced.Centers_exp1_h

PS1.Oxidized.Centers_exp1_h <- gghistogram(photo_geno_1expb, x = "PS1.Oxidized.Centers",
                            add = "mean", rug = TRUE,
                            color = "Treatment", fill = "Treatment", facet.by = "day.of.stress",
                            palette = c("blue", "orange")) + xlab("PS1.Oxidized.Centers (a.u.)")
## Warning: Using `bins = 30` by default. Pick better value with the argument
## `bins`.
PS1.Oxidized.Centers_exp1_h

SPAD_exp1_h <- gghistogram(photo_geno_1expb, x = "SPAD", binwidth = 2.5,
                            add = "mean", rug = TRUE,
                            color = "Treatment", fill = "Treatment", facet.by = "day.of.stress",
                            palette = c("blue", "orange")) + xlab("SPAD (a.u.)")
SPAD_exp1_h

### Look closer for potential outliers

off_Foprime <- subset(photo_geno_1expb, photo_geno_1expb$FoPrime < 200)
off_Foprime
off_leaf_thickness <- subset(photo_geno_1expb, photo_geno_1expb$leaf_thickness > 2.5)
off_leaf_thickness
off_PS1.Active.Centers <- subset(photo_geno_1expb, photo_geno_1expb$PS1.Active.Centers > 40)
off_PS1.Active.Centers
off_PS1.Open.Centers <- subset(photo_geno_1expb, photo_geno_1expb$PS1.Open.Centers < -50)
off_PS1.Open.Centers
off_PS1.Over.Reduced.Centers <- subset(photo_geno_1expb, photo_geno_1expb$PS1.Over.Reduced.Centers < -100)
off_PS1.Over.Reduced.Centers
off_PS1.Oxidized.Centers <- subset(photo_geno_1expb, photo_geno_1expb$PS1.Oxidized.Centers > 600)
off_PS1.Oxidized.Centers