TAC linearity summary
Kelsey Jesser
11/12/22
#clear R env
rm(list = ls())Load libraries
library(ggplot2)
library(ggpmisc)
library(plyr)
library(dplyr)
library(tidyr)Import data file
TAC<-read.csv("sample_data.csv", na.strings=c("","NA"))Plot a standard curve and summarize efficiency for each target
#Format data and add column with logs of quantities
TAC$CT<-as.numeric(gsub(",","",TAC$CT))
TAC$Quantity<-as.numeric(gsub(",","",TAC$Quantity))
TAC$Quantity_log<-log10(TAC$Quantity)
#summarize efficiencies
TAC<-TAC[!is.na(TAC$CT),]
TAC<-TAC[!is.na(TAC$Quantity_log),]
TAC<-TAC[rowSums(is.na(TAC)) != ncol(TAC), ]
efficiency_summ<-ddply(TAC, "Target.Name", function(x) {
model <- lm(CT ~ Quantity_log, data = x)
coef(model)})
colnames(efficiency_summ) <- c("Target.Name", "Y.intercept","Slope")
efficiency_summ$Efficiency<-10^(-1/efficiency_summ$Slope)-1
#create a list of dataframes by target
targets=by(TAC, TAC[,"Target.Name"], function(x) x)
#plot of standard curve points across all targets
p<-ggplot(TAC, aes(x = Quantity_log, y = CT)) +
geom_point()
p
#plot standard curves for all targets
target_plots <- function(targets){
ggplot(targets, aes(x = Quantity_log, y = CT)) +
geom_point()+
stat_poly_line(formula = y~x, se=FALSE) +
stat_poly_eq(aes(label = paste(after_stat(eq.label), after_stat(rr.label), sep = "*\", \"*")), formula= y~x, label.y=0.25)
}
q1 <- lapply(targets, target_plots)
q2 <- lapply(seq_along(q1), function(i) {
q1[[i]] + ggtitle(names(targets)[i])
})
#q2 #careful--this command makes a lot of plots!
#faceted plot with all adjusted standard curves for all targets
r<-ggplot(TAC, aes(x = Quantity_log, y = CT)) +
geom_point() +
facet_wrap(~ Target.Name)+
stat_poly_line(formula = y~x, se=FALSE) +
stat_poly_eq(aes(label = paste(after_stat(eq.label), after_stat(rr.label), sep = "*\", \"*")), formula= y~x, label.y=20)
r
Plot standard curves and summarize efficiencies without last point in the curve (where it flattens out as it approaches the LOD)
#remove last point on standard curve
TAC_adj <- TAC[TAC$Quantity_log != 0, ]
#TAC_adj <- TAC[TAC$Quantity_log != 0, ] #run this to look at standard curves without last two dilutions to improve efficencies of AMR gene assays and other targets with lower LODs
#create a list of dataframes by target
targets_adj=by(TAC_adj, TAC_adj[,"Target.Name"], function(x) x)
#plot adjusted standard curves for all targets
target_plots_adj <- function(targets_adj){
ggplot(targets_adj, aes(x = Quantity_log, y = CT)) +
geom_point()+
stat_poly_line(formula = y~x, se=FALSE) +
stat_poly_eq(aes(label = paste(after_stat(eq.label), after_stat(rr.label), sep = "*\", \"*")), formula= y~x, label.y=0.25)
}
s1 <- lapply(targets_adj, target_plots_adj)
s2 <- lapply(seq_along(s1), function(i) {
s1[[i]] + ggtitle(names(targets_adj)[i])
})
#s2 #careful-this command makes a lot of plots!
#faceted plot with all adjusted standard curves for all targets
t<-ggplot(TAC_adj, aes(x = Quantity_log, y = CT)) +
geom_point() +
facet_wrap(~ Target.Name)+
stat_poly_line(formula = y~x, se=FALSE) +
stat_poly_eq(aes(label = paste(after_stat(eq.label), after_stat(rr.label), sep = "*\", \"*")), formula= y~x, label.y=20)
t
Summarize dilutions with detection and high and low CT values for each target FIX THIS
TAC_CT<-TAC %>% select(Target.Name, Sample.Name, CT, R.superscript.2.)
high_low<-TAC_CT %>%
group_by(Target.Name) %>%
mutate(
CT_min=min(CT, na.rm=T),
CT_max=max(CT, na.rm=T)) %>%
mutate(
Dilution_min=min(Sample.Name, na.rm=T),
Dilution_max=max(Sample.Name, na.rm=T)) %>%
select(-Sample.Name) %>%
select(-CT)Summarize adjusted (without lowest dilution) slope, efficiency, and R2 values for each target
R2<- TAC_adj %>%
group_by(Target.Name) %>%
dplyr::summarise(model = list(lm(CT ~ Quantity_log, data = cur_data())),
coef = list(coef(model[[1]])),
R2_adj = summary(model[[1]])$r.sq)%>%
unnest_wider(coef, names_repair = 'unique')
R2$Efficiency_adj<- 10^(-1/R2$Quantity_log)-1Summarize NTC results
TAC_NTC<-read.csv("sample_data.csv", na.strings=c("","NA"))
TAC_NTC <- TAC_NTC %>%
subset(Sample.Name=="NTC")
TAC_NTC$CT<-as.numeric(gsub(",","",TAC_NTC$CT))
TAC_NTC<- TAC_NTC%>%
select(Target.Name, CT)%>%
dplyr::rename(NTC_CT=CT)Summary table
summary<-merge(efficiency_summ, R2)
summary<-merge(high_low, summary)
summary<-merge(TAC_NTC, summary)
summary <- summary %>%
distinct(.keep_all = TRUE) %>%
select(-model) %>%
relocate(R.superscript.2., .before=Efficiency) %>%
relocate(NTC_CT, .before=Y.intercept)%>%
dplyr::rename(R2=R.superscript.2.) %>%
dplyr::rename(Y.intercept_adj='(Intercept)') %>%
dplyr::rename(Slope_adj=Quantity_log)%>%
mutate_if(is.numeric, round, digits=2)
summary## Target.Name CT_min CT_max Dilution_min Dilution_max NTC_CT
## 1 aaiC 14.07 32.61 PCP 10^1 PCP 10^6 NA
## 2 aatA 12.37 31.06 PCP 10^1 PCP 10^6 NA
## 3 afaB 15.97 35.63 PCP 10^1 PCP 10^6 NA
## 4 aggR 13.97 36.36 PCP 10^0 PCP 10^6 NA
## 5 AMR CTX-M1 14.99 33.21 PCP 10^1 PCP 10^6 NA
## 6 AMR CTX-M2-M74 17.77 37.26 PCP 10^1 PCP 10^6 NA
## 7 AMR CTX-M8-M25 15.03 33.90 PCP 10^1 PCP 10^6 NA
## 8 AMR CTX-M9 16.26 35.41 PCP 10^1 PCP 10^6 NA
## 9 AMR Intl1 15.62 34.34 PCP 10^1 PCP 10^6 NA
## 10 AMR KPC 14.57 32.94 PCP 10^1 PCP 10^6 NA
## 11 AMR NDM 15.10 34.74 PCP 10^0 PCP 10^6 NA
## 12 AMR qnrA 15.82 35.65 PCP 10^0 PCP 10^6 NA
## 13 AMR qnrB1 15.64 33.94 PCP 10^1 PCP 10^6 NA
## 14 AMR qnrB4 15.37 33.16 PCP 10^1 PCP 10^6 NA
## 15 AMR SHV 14.64 33.00 PCP 10^1 PCP 10^6 NA
## 16 AMR sul1 16.16 33.74 PCP 10^1 PCP 10^6 NA
## 17 AMR sul2 15.20 32.64 PCP 10^1 PCP 10^6 NA
## 18 AMR TEM 15.43 34.76 PCP 10^0 PCP 10^6 NA
## 19 AMR VIM 13.55 33.12 PCP 10^0 PCP 10^6 NA
## 20 Ascaris 12.90 31.42 PCP 10^1 PCP 10^6 NA
## 21 AstroV capsid 11.07 34.54 PCP 10^0 PCP 10^6 NA
## 22 bfpA 13.32 33.36 PCP 10^0 PCP 10^6 NA
## 23 C. hominus 13.84 32.51 PCP 10^1 PCP 10^6 NA
## 24 C. parvum 13.27 34.62 PCP 10^0 PCP 10^6 NA
## 25 cadF 13.05 32.07 PCP 10^1 PCP 10^6 NA
## 26 Crypto 18S 14.72 34.11 PCP 10^1 PCP 10^6 NA
## 27 Cyclospora 13.50 34.17 PCP 10^0 PCP 10^6 NA
## 28 eae 13.35 32.31 PCP 10^1 PCP 10^6 NA
## 29 Entamoeba 16.95 35.51 PCP 10^1 PCP 10^6 NA
## 30 EntV 5' UTR 10.21 32.61 PCP 10^0 PCP 10^6 NA
## 31 Fecal bacteria 16S 17.17 34.60 PCP 10^1 PCP 10^6 NA
## 32 Giardia 13.05 34.34 PCP 10^0 PCP 10^6 NA
## 33 GlyA 13.55 34.79 PCP 10^0 PCP 10^6 NA
## 34 HAdV fiber 15.37 34.80 PCP 10^1 PCP 10^6 NA
## 35 HAdV hexon 13.26 32.87 PCP 10^1 PCP 10^6 NA
## 36 hipO 15.69 37.43 PCP 10^0 PCP 10^6 NA
## 37 ipaH 12.98 33.33 PCP 10^0 PCP 10^6 NA
## 38 LT 13.26 31.89 PCP 10^1 PCP 10^6 NA
## 39 MS2g1 11.46 35.07 PCP 10^0 PCP 10^6 NA
## 40 Noro GII 11.43 35.23 PCP 10^0 PCP 10^6 NA
## 41 PhHV gB 13.12 32.97 PCP 10^1 PCP 10^6 NA
## 42 Rotavirus 13.71 35.65 PCP 10^0 PCP 10^6 NA
## 43 Sapovirus 17.88 38.02 PCP 10^1 PCP 10^6 NA
## 44 SARS-CoV-2 N1 12.83 34.63 PCP 10^0 PCP 10^6 10.66
## 45 SARS-CoV-2 N2 12.91 37.00 PCP 10^0 PCP 10^6 15.00
## 46 STh 13.98 35.67 PCP 10^0 PCP 10^6 NA
## 47 STp 14.46 35.35 PCP 10^0 PCP 10^6 NA
## 48 stx1 15.32 33.46 PCP 10^1 PCP 10^6 NA
## 49 stx2 14.74 36.41 PCP 10^0 PCP 10^6 NA
## 50 STY0201 14.15 33.11 PCP 10^1 PCP 10^6 NA
## 51 Trichuris 13.40 31.99 PCP 10^1 PCP 10^6 NA
## 52 ttr 13.37 32.00 PCP 10^1 PCP 10^6 NA
## 53 tviB 9.30 27.70 PCP 10^1 PCP 10^6 NA
## 54 uidA 13.79 31.71 PCP 10^1 PCP 10^6 NA
## 55 virF 15.04 34.82 PCP 10^1 PCP 10^6 NA
## Y.intercept Slope R2 Efficiency Y.intercept_adj Slope_adj R2_adj
## 1 36.67 -3.63 0.99 0.88 36.67 -3.63 0.99
## 2 34.86 -3.59 0.99 0.90 34.86 -3.59 0.99
## 3 39.42 -3.83 1.00 0.82 39.42 -3.83 1.00
## 4 36.50 -3.64 1.00 0.88 36.63 -3.67 1.00
## 5 36.98 -3.45 0.98 0.95 36.98 -3.45 0.98
## 6 41.08 -3.82 1.00 0.83 41.08 -3.82 1.00
## 7 37.53 -3.66 1.00 0.88 37.53 -3.66 1.00
## 8 39.17 -3.86 1.00 0.82 39.17 -3.86 1.00
## 9 38.04 -3.66 1.00 0.88 38.04 -3.66 1.00
## 10 37.26 -3.63 0.99 0.89 37.26 -3.63 0.99
## 11 35.88 -3.26 0.98 1.03 37.94 -3.71 0.99
## 12 36.75 -3.32 0.99 1.00 37.71 -3.54 0.99
## 13 37.77 -3.57 1.00 0.91 37.77 -3.57 1.00
## 14 37.06 -3.47 0.99 0.94 37.06 -3.47 0.99
## 15 37.06 -3.59 0.99 0.90 37.06 -3.59 0.99
## 16 37.91 -3.54 0.99 0.92 37.91 -3.54 0.99
## 17 36.90 -3.57 0.99 0.91 36.90 -3.57 0.99
## 18 36.05 -3.24 0.98 1.04 37.18 -3.50 0.99
## 19 34.11 -3.22 0.98 1.04 35.89 -3.62 0.99
## 20 35.44 -3.64 1.00 0.88 35.44 -3.64 1.00
## 21 34.83 -3.93 1.00 0.80 35.08 -3.99 1.00
## 22 34.52 -3.35 0.99 0.99 35.52 -3.58 0.99
## 23 36.45 -3.68 1.00 0.87 36.45 -3.68 1.00
## 24 34.89 -3.46 1.00 0.95 35.13 -3.52 0.99
## 25 36.14 -3.73 0.99 0.85 36.14 -3.73 0.99
## 26 37.83 -3.78 1.00 0.84 37.83 -3.78 1.00
## 27 34.99 -3.38 0.99 0.98 36.47 -3.71 0.99
## 28 36.47 -3.73 0.99 0.86 36.47 -3.73 0.99
## 29 38.99 -3.57 0.99 0.90 38.99 -3.57 0.99
## 30 33.06 -3.71 1.00 0.86 33.45 -3.80 1.00
## 31 38.88 -3.47 0.99 0.94 38.88 -3.47 0.99
## 32 34.76 -3.47 1.00 0.94 35.13 -3.56 0.99
## 33 36.28 -3.60 0.98 0.90 37.58 -3.90 1.00
## 34 38.58 -3.80 1.00 0.83 38.58 -3.80 1.00
## 35 36.50 -3.79 1.00 0.84 36.50 -3.79 1.00
## 36 38.88 -3.72 0.97 0.86 40.51 -4.09 0.99
## 37 34.52 -3.38 0.99 0.98 35.56 -3.61 0.99
## 38 36.05 -3.67 0.99 0.87 36.05 -3.67 0.99
## 39 35.01 -3.79 1.00 0.84 34.96 -3.77 0.99
## 40 35.11 -3.85 1.00 0.82 35.01 -3.83 1.00
## 41 36.73 -3.83 1.00 0.82 36.73 -3.83 1.00
## 42 36.22 -3.65 1.00 0.88 36.70 -3.76 1.00
## 43 42.25 -4.00 1.00 0.78 42.25 -4.00 1.00
## 44 35.17 -3.61 1.00 0.89 35.64 -3.72 1.00
## 45 37.26 -4.00 1.00 0.78 37.48 -4.05 0.99
## 46 36.28 -3.57 0.99 0.91 36.81 -3.69 0.99
## 47 36.38 -3.52 0.99 0.92 37.27 -3.73 1.00
## 48 37.57 -3.58 0.99 0.90 37.57 -3.58 0.99
## 49 36.96 -3.56 1.00 0.91 37.44 -3.66 1.00
## 50 36.94 -3.67 0.99 0.87 36.94 -3.67 0.99
## 51 36.02 -3.67 1.00 0.87 36.02 -3.67 1.00
## 52 36.00 -3.64 0.99 0.88 36.00 -3.64 0.99
## 53 31.89 -3.49 0.97 0.93 31.89 -3.49 0.97
## 54 35.75 -3.51 0.99 0.93 35.75 -3.51 0.99
## 55 38.45 -3.81 1.00 0.83 38.45 -3.81 1.00
## Efficiency_adj
## 1 0.88
## 2 0.90
## 3 0.82
## 4 0.87
## 5 0.95
## 6 0.83
## 7 0.88
## 8 0.82
## 9 0.88
## 10 0.89
## 11 0.86
## 12 0.92
## 13 0.91
## 14 0.94
## 15 0.90
## 16 0.92
## 17 0.91
## 18 0.93
## 19 0.89
## 20 0.88
## 21 0.78
## 22 0.90
## 23 0.87
## 24 0.93
## 25 0.85
## 26 0.84
## 27 0.86
## 28 0.86
## 29 0.90
## 30 0.83
## 31 0.94
## 32 0.91
## 33 0.81
## 34 0.83
## 35 0.84
## 36 0.76
## 37 0.89
## 38 0.87
## 39 0.84
## 40 0.82
## 41 0.82
## 42 0.84
## 43 0.78
## 44 0.86
## 45 0.77
## 46 0.87
## 47 0.85
## 48 0.90
## 49 0.87
## 50 0.87
## 51 0.87
## 52 0.88
## 53 0.93
## 54 0.93
## 55 0.83write.csv(summary,file="TAC_sample_linearity_summary.csv")