diff --git a/R/Global.R b/R/Global.R index 4bf9614..8591e82 100644 --- a/R/Global.R +++ b/R/Global.R @@ -40,6 +40,7 @@ library(scales) COR_FUNC <- function(vec1, vec2) { df <- data.frame(v1 = vec1, v2 = vec2) df2 <- df[complete.cases(df),] + #browser() return(cor(df2[,1],df2[,2])) } @@ -942,12 +943,13 @@ PlotLinPLA_FUNC <- function(circle, sigmoid, all_l2, pl_df, indS, indT) { truePL_df <- NULL } + p <- ggplot(all_l2, aes(x = log_dose, y = readout, color = factor(isRef))) + geom_point(size = 2) + # labs(title=paste("linear regression model", indS,indT), color="product") + scale_colour_manual(labels = c("test", "reference"), values = c("#C2173F", "#4545BA")) + - ylim(min(all_l2$readout), max(all_l2$readout)) + + scale_x_continuous(breaks = scales::pretty_breaks(n = 10)) + scale_y_continuous(breaks = scales::pretty_breaks(n = 10)) + theme_bw() @@ -981,6 +983,7 @@ PlotLinPLA_FUNC <- function(circle, sigmoid, all_l2, pl_df, indS, indT) { x = log_dose, y = readout, shape = factor(isRef), size = 5, alpha = 0.2 ), col = c("black"), inherit.aes = FALSE) + + ylim(min(all_l2$readout), max(all_l2$readout)) + scale_shape_manual(labels = c("test", "reference"), values = c(21, 21)) # fit intercept for test and ref and common slope @@ -1017,6 +1020,7 @@ PlotLinPLA_FUNC <- function(circle, sigmoid, all_l2, pl_df, indS, indT) { title = paste("restricted linear regression model"), subtitle = paste("Regression on highlighted points") ) + + ylim(min(all_l2$readout), max(all_l2$readout)) + theme(legend.position = "none", axis.text = element_text(size = 14)) pr3 <- pr2 + geom_point(circle, mapping = aes( x = log_dose, y = readout, shape = factor(isRef), @@ -1219,7 +1223,7 @@ tests_FUNC <- function(ro_new, Lim, PureErrFlag) { VCOVpure <- V_V * meanPureErr DFsPure <- FitAnova[4, 1] - +#browser() testPOTr <- logical() if (POTr_CI[1] * 100 > Lim[[9]] & POTr_CI[2] * 100 < Lim[[10]]) testPOTr <- 0 else testPOTr <- 1 diff --git a/Tests3Plates.numbers b/Tests3Plates.numbers index adf4ae2..72d3d18 100644 Binary files a/Tests3Plates.numbers and b/Tests3Plates.numbers differ diff --git a/Tests4Plates.xlsx b/Tests4Plates.xlsx new file mode 100644 index 0000000..5c00054 Binary files /dev/null and b/Tests4Plates.xlsx differ diff --git a/dev/.DS_Store b/dev/.DS_Store index 2c3cac7..532ee15 100644 Binary files a/dev/.DS_Store and b/dev/.DS_Store differ diff --git a/dev/app.R b/dev/app.R index 2ad4cb1..b3ed2c9 100644 --- a/dev/app.R +++ b/dev/app.R @@ -118,8 +118,32 @@ server <- function(input, output, session) { "The reference readout columns have to be before the test sample readout columns. The column names for reference and test are free to set, but different for all columns.", "The column name of the dilution concentrations can be as follows: concentration, dose, log_concentration, log_dose (first letter can be capital)", "If the concentrations are in logarithmized, any log base can be used.", + br(), br(), + "EXPLORE the 4pl function: visualize the meta data in the context of a 4 PL fit or a linear regression fit. ", + "Enter the 4 parameters of test and reference sample and see, what this means.", br(), + br(), + "OPTIMIZE the concentrations: plot all results you have and adjust the concentrations accordingly. ", + "Get help if you want to read pdfs in contacting us.", br(), ), - column(6, ) + column(2, + style = "background: #7FAEFF88", + "HERE: Enter the equivalence limits for 4PL suitability tests. If you need help to set them, contact us.", + numericInput("lEACratiola", "lower EAC ratio of LAs", 0.005, step = 0.001), + numericInput("uEACratiola", "upper EAC for ratio of LAs", 100, step = 1), + numericInput("lEACratioSlope", "lower EAC for ratio of slopes", 0.55, step = 0.01), + numericInput("uEACratioSlope", "upper EAC for ratio of slopes", 1.84, step = 0.1), + numericInput("lEACratioua", "lower EAC for ratio of UAs", 0.75, step = 0.1), + numericInput("uEACratioua", "upper EAC for ratio of UAs", 1.33, step = 0.1) + ), + column(2, + style = "background: #7FAEFF88", + numericInput("lowerPot", "lower EAC for potency", 75, step = 1), + numericInput("upperPot", "upper EAC for potency", 133, step = 1), + numericInput("lEACratioAdiff", "lower EAC of ratio of asymptote differences", 0.75, step = 0.01), + numericInput("uEACratioAdiff", "upper EAC of ratio of asymptote differences", 1.33, step = 0.01), + numericInput("lEACdiffla", "lower EAC for diff. of LA", -0.175, step = 0.001), + numericInput("uEACdiffla", "upper EAC for diff. of LA", 0.189, step = 0.001) + ) ), tabPanel( "Documentation", @@ -188,24 +212,6 @@ server <- function(input, output, session) { ), selected = c("1", "2", "3", "4", "5", "6", "7", "8") ) - ), - column(2, - style = "background: #7FAEFF88", - numericInput("lEACratiola", "lower EAC ratio of LAs", 0.005, step = 0.001), - numericInput("uEACratiola", "upper EAC for ratio of LAs", 100, step = 1), - numericInput("lEACratioSlope", "lower EAC for ratio of slopes", 0.55, step = 0.01), - numericInput("uEACratioSlope", "upper EAC for ratio of slopes", 1.84, step = 0.1), - numericInput("lEACratioua", "lower EAC for ratio of UAs", 0.75, step = 0.1), - numericInput("uEACratioua", "upper EAC for ratio of UAs", 1.33, step = 0.1) - ), - column(2, - style = "background: #7FAEFF88", - numericInput("lowerPot", "lower EAC for potency", 75, step = 1), - numericInput("upperPot", "upper EAC for potency", 133, step = 1), - numericInput("lEACratioAdiff", "lower EAC of ratio of asymptote differences", 0.75, step = 0.01), - numericInput("uEACratioAdiff", "upper EAC of ratio of asymptote differences", 1.33, step = 0.01), - numericInput("lEACdiffla", "lower EAC for diff. of LA", -0.175, step = 0.001), - numericInput("uEACdiffla", "upper EAC for diff. of LA", 0.189, step = 0.001) ) ), @@ -559,8 +565,8 @@ server <- function(input, output, session) { tabPanel( "Report", h4("Settings for report"), - downloadButton("downloadXLReport", label = "Download PDF report", class = "butt"), - tags$style(type = "text/css", "#downloadXLReport {background-color: orange; color: black;font-family: COurier New}"), + downloadButton("downloadXLReportMeta", label = "Download PDF report", class = "butt"), + tags$style(type = "text/css", "#downloadXLReportMeta {background-color: orange; color: black;font-family: COurier New}"), ) ) ) @@ -632,7 +638,11 @@ server <- function(input, output, session) { ), tabPanel("Histograms", h4("Histograms of parameters"), - plotOutput("histCIs"), + plotOutput("linerangeCIs"), + plotOutput("ratioSlope"), + plotOutput("ratioLas"), + plotOutput("ratioUas"), + plotOutput("widthCIs"), column(6, plotOutput("histEC50REF"), plotOutput("histLasREF"), @@ -647,7 +657,16 @@ server <- function(input, output, session) { ), tabPanel( "Report", - h4("Settings for report")) + h4("Settings for report"), + useShinyjs(), + actionButton("btn2", "Download PDF report", icon = icon("download")), + downloadButton("downloadWizardReport", label = "Download report", class = "butt"), + tags$style(type = "text/css", "#downloadWizardReport {background-color: #4FCBD9; color: black;font-family: Courier New}"), + # textInput("Author", "Author", value = ""), + # textInput("RepIdentifier", "Report name", value = ""), + # textInput("NoP", "Product name", value = ""), + # textInput("Assay", "Assay name", value = "") + ) ) ) # main panel @@ -2003,7 +2022,7 @@ server <- function(input, output, session) { pottab4_$`upper95%CI` <- round(as.numeric(pottab4[, 4]) * 100, 2) pottab4_$relative_lowerCL <- round(pottab4_[, 6] / pottab4_[, 5] * 100, 2) pottab4_$relative_upperCL <- round(pottab4_[, 7] / pottab4_[, 5] * 100, 2) - +#browser() if (as.numeric(pottab4_$relative_lowerCL[1]) > Lim[[9]] & as.numeric(pottab4_$relative_upperCL[1]) < Lim[[10]]) { test_potCI <- 0 } else { @@ -2129,7 +2148,7 @@ server <- function(input, output, session) { }) - #### Dilutions Simulator ---- + #### Meta plots all XL ---- observe({ if (!is.null(Dat$Mws)) { @@ -2194,8 +2213,9 @@ server <- function(input, output, session) { } } #for N_WS - #browser() +#browser() URMcoefsDF <- t(matrix(unlist(URMcoefsL),nrow=9)) + colnames(URMcoefsDF) <- c("WS name", "lowerAs REF","slope REF","upperAs REF","EC50 REF", "lowerAs TEST","slope TEST","upperAs TEST","EC50 Difference") EC50TEST <- as.numeric(URMcoefsDF[,5]) - as.numeric(URMcoefsDF[,9]) # EC50TEST <- EC50TEST[!EC50TEST %in% boxplot.stats(EC50TEST)$out] EC50REF <- as.numeric(URMcoefsDF[,5]) @@ -2204,8 +2224,15 @@ server <- function(input, output, session) { # UasREF <- UasREF[!UasREF %in% boxplot.stats(UasREF)$out] LasREF <- as.numeric(URMcoefsDF[,2]) # LasREF <- LasREF[!LasREF %in% boxplot.stats(LasREF)$out] - UasTEST <- as.numeric(URMcoefsDF[,4]) - LasTEST <- as.numeric(URMcoefsDF[,2]) + UasTEST <- as.numeric(URMcoefsDF[,8]) + LasTEST <- as.numeric(URMcoefsDF[,6]) + slopeREF <- as.numeric(URMcoefsDF[,3]) + slopeTEST <- as.numeric(URMcoefsDF[,7]) + + slopeRatio <- slopeTEST/slopeREF + LasRatio <- LasTEST/LasREF + UasRatio <- UasTEST/UasREF + ratioDF <- data.frame(WS_name = URMcoefsDF[,1], slopeRatio = slopeRatio, LasRatio = LasRatio, UasRatio = UasRatio) RMcoefsDF <- t(matrix(unlist(RMcoefsL),nrow=6)) @@ -2213,6 +2240,7 @@ server <- function(input, output, session) { Dat$RestrM <- RMcoefsDF CalcPotDF <- t(matrix(unlist(potEstL),nrow=3)) + colnames(CalcPotDF) <- c("rel_potency","lower_CI","upper_CI") Dat$CalcPot <- CalcPotDF # #### sigmoid plots ---- @@ -2241,7 +2269,7 @@ server <- function(input, output, session) { expand_limits(x = c(min(SIGrefDF$X) - 0.1 * diff(range(SIGrefDF$X)), max(SIGrefDF$X) + 0.1 * diff(range(SIGrefDF$X)))) + xlab("dilutions") + - #ggtitle("Plot of all calculated reference fits (unrestricted model, in gray vertical lines: EC50)") + + ggtitle("REF sample 4PL-fits (unrestricted model, gray vertical lines: EC50)") + theme_bw() + theme(axis.text = element_text(face = "bold", size = 15), plot.title = element_text(size = 15, face = "bold"), @@ -2261,7 +2289,7 @@ server <- function(input, output, session) { geom_hline(yintercept = c(LasParTolTEST[[4]], LasParTolTEST[[5]]), linetype=2, col="grey") + geom_vline(xintercept = EC50TEST, alpha = 0.2) + xlab("dilutions") + - ggtitle("Calculated test sample fits (unrestricted model, in gray vertical lines: EC50)") + + ggtitle("TEST sample 4PL-fits (unrestricted model, gray vertical lines: EC50)") + theme_bw() + theme(axis.text = element_text(face = "bold", size = 15), plot.title = element_text(size = 15, face = "bold")) @@ -2269,29 +2297,90 @@ server <- function(input, output, session) { output$sigPlotTEST <- renderPlot({ p2 }) Dat$sigPlotTEST <- p2 - + #### histograms right panel ---- - #browser() + - all_lPot <- data.frame(Cat_potency= c(rep("rel poteny",nrow(CalcPotDF)), rep("lower CI",nrow(CalcPotDF)),rep("upper CI",nrow(CalcPotDF))), + all_lPot <- data.frame(Cat_potency= c(rep("rel_poteny",nrow(CalcPotDF)), rep("lower_CI",nrow(CalcPotDF)),rep("upper_CI",nrow(CalcPotDF))), Potency_and_CI = c(CalcPotDF[,1], CalcPotDF[,2],CalcPotDF[,3])) all_lPot[,2][all_lPot[,2] > 5] <- NA all_lPot[,2][all_lPot[,2] < 0.1] <- NA - P_histCI <- ggplot(all_lPot, aes(x=Potency_and_CI, fill=Cat_potency)) + + CalcPotDF <- as.data.frame(CalcPotDF) + CalcPotDF$width_CI <- CalcPotDF$upper_CI - CalcPotDF$lower_CI + widthCLTol <- normtol.int(x = CalcPotDF$width_CI, alpha = 1-input$TolConf, P = input$TolPop, side = 1) + + P_linerangeCI <- ggplot(CalcPotDF, aes(x=seq(1,nrow(CalcPotDF)))) + #, aes(x=Potency_and_CI, fill=Cat_potency) + geom_linerange(aes(ymin=lower_CI, ymax=upper_CI), color="black") + + geom_point(aes(y=rel_potency), alpha = 0.1) + + #scale_fill_manual(values=c("darkgreen","darkblue","salmon2","tomato3")) + + ggtitle("CLs of relative potencies, standard RMSEs") + + # scale_x_continuous( + # breaks=seq(trunc(min(all_lPot$Potency_and_CI, na.rm=T)*10)/10, max(all_lPot$Potency_and_CI, na.rm=T)*1.1, by=0.4), + # ) + + theme_bw() + + theme(axis.text = element_text(face="bold", size=15), + axis.text.x = element_text(angle=90), + plot.title= element_text(size=15, face="bold")) + #P_linerangeCI + output$linerangeCIs <- renderPlot({ P_linerangeCI }) + + P_widthCIs <- ggplot(CalcPotDF, aes(x=width_CI, fill="blue")) + geom_histogram(color="#e9ecef", alpha=0.6, position = "identity") + - scale_fill_manual(values=c("darkgreen","darkblue","salmon2","tomato3")) + - ggtitle("Histogram of relative potencies, standard RMSEs") + - scale_x_continuous( - breaks=seq(trunc(min(all_lPot$Potency_and_CI, na.rm=T)*10)/10, max(all_lPot$Potency_and_CI, na.rm=T)*1.1, by=0.4), - ) + + geom_density(alpha = 0.1) + + #scale_fill_manual(values=c("darkgreen","darkblue","salmon2","tomato3")) + + labs(title = "Histogram of width of CLs", subtitle = paste("with upper",input$TolConf,input$TolPop, "tolerance interval")) + + geom_vline(xintercept = widthCLTol[[5]]) + + theme_bw() + + theme(axis.text = element_text(face="bold", size=15), + axis.text.x = element_text(angle=90), + plot.title= element_text(size=15, face="bold")) + + output$widthCIs <- renderPlot({ P_widthCIs }) + + SlopeTol <- normtol.int(x = slopeRatio, alpha = 1-input$TolConf, P = input$TolPop, side = 2) + LasTol <- normtol.int(x = LasRatio, alpha = 1-input$TolConf, P = input$TolPop, side = 2) + UasTol <- normtol.int(x = UasRatio, alpha = 1-input$TolConf, P = input$TolPop, side = 2) + + P_ratioSlope <- ggplot(ratioDF, aes(x=slopeRatio, fill="turquoise")) + + geom_histogram(color="#e9ecef", alpha=0.6, position = "identity") + + geom_density(alpha = 0.1) + + #scale_fill_manual(values=c("darkgreen","darkblue","salmon2","tomato3")) + + labs(title = "Histogram of Hill's slope ratios ", subtitle = paste("with",input$TolConf,input$TolPop, "tolerance interval")) + + geom_vline(xintercept = c(SlopeTol[[4]], SlopeTol[[5]]), linetype=2, col="grey") + theme_bw() + theme(axis.text = element_text(face="bold", size=15), axis.text.x = element_text(angle=90), plot.title= element_text(size=15, face="bold")) - output$histCIs <- renderPlot({ P_histCI }) + output$ratioSlope <- renderPlot({ P_ratioSlope }) + + P_ratioLas <- ggplot(ratioDF, aes(x=LasRatio, fill="turquoise")) + + geom_histogram(color="#e9ecef", alpha=0.6, position = "identity") + + geom_density(alpha = 0.1) + + #scale_fill_manual(values=c("darkgreen","darkblue","salmon2","tomato3")) + + labs(title = "Histogram of lower asymptote ratios ", subtitle = paste("with",input$TolConf,input$TolPop, "tolerance interval")) + + geom_vline(xintercept = c(LasTol[[4]], LasTol[[5]]), linetype=2, col="grey") + + theme_bw() + + theme(axis.text = element_text(face="bold", size=15), + axis.text.x = element_text(angle=90), + plot.title= element_text(size=15, face="bold")) + + output$ratioLas <- renderPlot({ P_ratioLas }) + + P_ratioUas <- ggplot(ratioDF, aes(x=UasRatio, fill="turquoise")) + + geom_histogram(color="#e9ecef", alpha=0.6, position = "identity") + + geom_density(alpha = 0.1) + + #scale_fill_manual(values=c("darkgreen","darkblue","salmon2","tomato3")) + + labs(title = "Histogram of upper asymptote ratios ", subtitle = paste("with",input$TolConf,input$TolPop, "tolerance interval")) + + geom_vline(xintercept = c(UasTol[[4]], UasTol[[5]]), linetype=2, col="grey") + + theme_bw() + + theme(axis.text = element_text(face="bold", size=15), + axis.text.x = element_text(angle=90), + plot.title= element_text(size=15, face="bold")) + + output$ratioUas <- renderPlot({ P_ratioUas }) output$histEC50REF <- renderPlot({ hist(EC50REF, col="steelblue", border="white", main = 'Histogram of EC50REF') @@ -2315,11 +2404,12 @@ server <- function(input, output, session) { Dat$histEC50REF <- hist(EC50REF, col="steelblue", border="white", main = 'Histogram of EC50REF') Dat$histLasREF <- hist(LasREF, col="violet", border="white", main = 'Histogram of EC50REF') Dat$histUasREF <- hist(UasREF, col="darkturquoise", border="white", main = 'Histogram of EC50REF') - + + ##### Dilutions Simulator ---- tab <- AllXL[[1]] dils <- tab$log_dose - min_y <- min(tab[, 1:2]) - max_y <- max(tab[, 1:2]) + min_y <- min(tab[, 1:2], na.rm = T) + max_y <- max(tab[, 1:2], na.rm = T) #browser() # if (input$fixupper) { # dils_av <- dils - max(dils) @@ -2389,16 +2479,17 @@ server <- function(input, output, session) { ) DilsTable }) - +##### Plot for dilution slider ---- + if (!is.null(p2)) { #p2 <- Dat$p2 p_dil <- p2 + - annotate("pointrange", x = dils2, y = rep(min_y, length(dils2)), xmin = min(dils2), xmax = max(dils2)) + - annotate("text", x = dils2, y = rep(min_y + (max_y - min_y) * 0.05, length(dils2)), label = as.character(round(dils2, 3))) + + annotate("pointrange", x = dils2, y = rep(min_y, length(dils2)), xmin = min(dils2), xmax = max(dils2),colour = "red" ,linetype = 3) + + annotate("text", x = dils2, y = rep(min_y + (max_y - min_y) * 0.05, length(dils2)), label = as.character(round(dils2, 3)),colour = "red") + annotate("text", x = dils2[-1] + (max(dils2) - min(dils2)) * 0.05, y = rep(min_y + (max_y - min_y) * 0.1, length(dils2[-1])), - label = as.character(round(dilfactors[-1], 3))) + label = as.character(round(dilfactors[-1], 3)),colour = "red") # geom_line( # data = as.data.frame(pl_df), aes(x = dils2, y = SAMPLE50), color = "grey15", linetype = 2, # inherit.aes = F @@ -2512,7 +2603,7 @@ server <- function(input, output, session) { }) }) - #### simulations ---- + #### NOT SHOWN: simulations ---- observe({ observeEvent(input$goSim, { sd_fac_ <- as.numeric(input$sdfac) @@ -2625,7 +2716,7 @@ server <- function(input, output, session) { }) - #### simulation Histograms output ---- + #### NOT SHOWN: simulation Histograms output ---- output$plotHistuAs <- renderPlot({ if (!is.null(Dat$resHist)) { @@ -2756,6 +2847,76 @@ server <- function(input, output, session) { ) } ) + + + #### download Meta 4PL report---- + + observeEvent(input$btn2, { + if(!Dat$FITsFlag) { + runjs("$('#downloadXLReportMeta')[0].click();") + } else { + showModal(modalDialog( + title = "No 4PL model to Download", + "Please select other data before download.", + easyClose = TRUE, + footer = NULL + )) + } + }) + + output$downloadXLReportMeta <- downloadHandler( + filename = paste0("Report_4PLEvaluation", Dat$RepIdentifier, ".pdf"), + content = function(file) { + tpdr <- tempdir() + tempReport <- file.path(tpdr, "Doc_BioassayReport.Rmd") + file.copy("Doc_BioassayReport.Rmd", tempReport, overwrite = T) + + tempReportc <- file.path(tpdr, "logov2.png") + file.copy("logov2.png", tempReportc, overwrite = T) + + rmarkdown::render(tempReport, + output_file = file, + params = list( + FileName = Dat$FileName, + author = Dat$Author, + NoP = Dat$NoP, + Assay = Dat$Assay, + REP = REP, + coeffs = Dat$coeffs_UN + ), + envir = new.env(parent = globalenv()) + ) + } + ) + + + #### download Wizard report ---- + + output$downloadWizardReport <- downloadHandler( + filename = paste0("Report_Meta", Dat$nameRep, ".pdf"), + content = function(file) { + tpdr <- tempdir() + tempReport <- file.path(tpdr, "Doc_BioassayLinReport.Rmd") + file.copy("Doc_BioassayLinReport.Rmd", tempReport, overwrite = TRUE) + + tempReportc <- file.path(tpdr, "logov2.png") + file.copy("logov2.png", tempReportc, overwrite = TRUE) + + rmarkdown::render(tempReport, + output_file = file, + params = list( + FileName = Dat$FileName, + author = Dat$Author, + NoP = Dat$NoP, + Assay = Dat$Assay, + REP = REP, + REPlin = REPlin, + coeffsLin = Dat$coeffs_UN + ), + envir = new.env(parent = globalenv()) + ) + } + ) } shinyApp(ui, server)