Files
plateflow/appPLAscraper.R

591 lines
25 KiB
R
Raw Normal View History

2026-01-11 18:15:51 +01:00
library(shiny)
library(shinyjs)
library(shinybusy)
library(gslnls)
library(tidyverse)
library(ggplot2)
library(reshape2)
library(pdftools)
library(stringr)
library(openxlsx)
library(evd)
# library(fitdistrplus)
RestrMF <- function(data, M) {
all_l <- melt(data, id.vars="uniDoses", variable.name="replname", value.name="readout")
all_l$sample <- c(rep(0,nrow(all_l)/2), rep(1,nrow(all_l)/2))
all_l$isRef <- c(rep(1,nrow(all_l)/2), rep(0,nrow(all_l)/2))
all_l$log_dose <- log(all_l$uniDoses)
CORro <- cor(all_l$uniDoses[1:nrow(all_l)], all_l$readout[1:nrow(all_l)])
if (CORro<0) B <- -1 else B <- 1
startlist <- list(a= min(all_l$readout), b = B, cs = mean(all_l$log_dose),
r = 0.0001, d = max(all_l$readout))
mr <- tryCatch({
gsl_nls(fn=readout ~ a + (d-a)/(1+exp(b*((cs-r*sample)-log_dose))),
data = all_l,
start= startlist, control = gsl_nls_control(xtol=1e-6, ftol=1e-6, gtol=1e-6))
}, error = function(msg) {
return(0)
})
if (length(mr)>1) {
s_mr_coef <- summary(mr)$coefficients
s_mr <- summary(mr)
} else {
s_mr_coef <- rep(0,4); s_mr <- 0
}
if (M == "R") {
return(s_mr_coef)
} else if (M == "CI") {
if (length(mr) == 1) return(rep(0,4))
RelPot_mr <- exp(s_mr_coef[4,1])
POTr_CI <- exp(confint(mr, "r", method="asymptotic"))
# Pure error
FitAnova <- anova(lm(readout ~ factor(uniDoses)*sample, all_l))
meanPureErr <- FitAnova[4,3]
SEsPure <- sqrt(diag(vcov(mr)/s_mr$sigma^2)*meanPureErr)
lCIpure <- exp(s_mr_coef[4,1]-qt(0.975,nrow(all_l)-16)*SEsPure[4])
uCIpure <- exp(s_mr_coef[4,1]+qt(0.975,nrow(all_l)-16)*SEsPure[4])
if (is.infinite(uCIpure)) uCIpure <- NA
if (is.infinite(POTr_CI[2])) POTr_CI[2] <- NA
relCIpure <- c(lCIpure/RelPot_mr, uCIpure/RelPot_mr)
return(c(RelPot_mr, POTr_CI, lCIpure, uCIpure))
} else if (M == "U") {
if (length(mr) == 1) return(rep(0,4))
#browser()
startlistmu <- list(as = min(all_l$readout), bs = B, cs = mean(all_l$log_dose), ds = max(all_l$readout),
at = min(all_l$readout), bt = B, ct = mean(all_l$log_dose), dt = min(all_l$readout))
mu <- tryCatch({
gsl_nls(readout ~ as*isRef + at*sample + (ds*isRef + dt*sample - as*isRef - at*sample)/
(1 + isRef*exp(bs*(cs-log_dose)) + sample*exp(bt*(ct-log_dose))),
data = all_l,
start = startlistmu,
control = gsl_nls_control(xtol=1e-10, ftol=1e-10, gtol=1e-10))
}, error = function(msg) {
return(0)
})
if (length(mu)>1) {
smu_coef <- tryCatch({
summary(mu)$coefficients
}, error = function(msg) {
return(rep(NA,8))
})
return(smu_coef)
} else return(rep(NA,8))
} else return()
}
ui <- fluidPage(
useShinyjs(),
titlePanel("PDF Scraper"),
tabsetPanel(id="toptab",
tabPanel("Data loading",
wellPanel(
fluidRow(
# use_busy_spinner(
# spin = "cube_grid",
# color = "#122440",
# position = "top-left",
# margins = c(200,330),
# height = "50px",
# width = "50px"
# ),
column(2,
tags$img(src = "Plateflow.png", class="adv_logo"),
h4("Highlight all pdfs (at least 3)"),
fileInput(inputId = "iFile", label = "", multiple = T, accept = "application/pdf"),
downloadButton('downloadPLA', "Download scraped data"),
# h5("Author:")
),
column(6, style = "background: lightgrey",
h4("Sigmoid plots (legend numbers refer to the order of uploaded pdfs - see list below)"),
plotOutput("sigPlotREF"),
plotOutput("sigPlotTEST"),
plotOutput("histCIs"),
plotOutput("histCIsPure"),
plotOutput("GumbelPlot"),
tableOutput("CalcPotDF")
),
column(4,
plotOutput("histEC50REF"),
plotOutput("histLasREF"),
plotOutput("histUasREF"),
tableOutput("ListofPDFsF")
)
)
))
)
)
#### Server ----
server <- function(input, output, session) {
Dat <- reactiveValues()
PLOTS <- reactiveValues()
observe({
if (!is.null(input$iFile)) {
inFile <- input$iFile
filenames <- inFile$name
ldf <- lapply(inFile$datapath, pdf_text)
Dat$PDFs <- ldf
Dat$FNames <- filenames
}
})
observeEvent(input$iFile, {
show_spinner()
})
observe({
if (!is.null(input$iFile)) {
PDFs <- Dat$PDFs
filenames <- Dat$FNames
POTlist <- FINALresponses <- SampleNamesL <- list()
DilSerVec <- vector()
hide_spinner()
for (Nfile in 1:length(PDFs)) {
tt = PDFs[[Nfile]]
PageI <- strg <- POT <- vector(); resIND <- matrix(NA, ncol=2,nrow=3)
T1flag <- T2flag <- StFlag <- 0
for (Seite in 1:length(tt)) {
Zeilen <- strsplit(tt[Seite], "\n")
REGEXloc <- grep(pattern='Dose .alue', Zeilen[[1]])
if (length(REGEXloc)>0) {
#### get Standard ----
for (n in 1:length(REGEXloc)) {
print(paste("Response page no.:", Seite, "location", REGEXloc[n]))
PageI <- c(PageI, Seite)
strg <- c(strg, REGEXloc[n])
if (StFlag == 0) {
grepStS <- grep("Standard .ample", Zeilen[[1]])
if (length(grepStS)>0) {
for (z in 1:length(grepStS)) {
if (abs(REGEXloc[n]-grepStS[z]) <5) {
resIND[1,] <- c(Seite, REGEXloc[n])
REF <- gsub("Standard .ample: ","",Zeilen[[1]][grepStS[z]])
StFlag <- 1
}
}
}
}
#### get Test ----
if (T1flag == 0) {
grepT3 <- grep("Test .ample", Zeilen[[1]])
if (length(grepT3)>0) {
for (z in 1:length(grepT3)) {
if (abs(REGEXloc[n]-grepT3[z]) <5) {
resIND[2,] <- c(Seite, REGEXloc[n])
TEST1 <- gsub("Test .ample: ","",Zeilen[[1]][grepT3[z]])
T1flag <- 1
}
}
} else TESTflag <- 1
}
if (T1flag == 1) {
grepT4 <- grep("Test .ample", Zeilen[[1]])
if (length(grepT4)>0) {
for (z in 1:length(grepT4)) {
if (abs(REGEXloc[n]-grepT4[z]) <5) {
resIND[3,] <- c(Seite, REGEXloc[n])
TEST2 <- gsub("Test .ample: ","",Zeilen[[1]][grepT4[z]])
T2flag <- 1
}
}
} else TESTflag <- 1
}
}
}
#### Potency from PDF ----
POTloc <- grep(pattern='Potency .atio', Zeilen[[1]])
if (length(POTloc)>0) {
for (n in 1:length(POTloc)) {
POT <- c(POT, Zeilen[[1]][POTloc:(POTloc+5)])
}
}
} # tt
if (sum(resIND[2,] == resIND[3,]) ==2) { resIND[3,] <- c(NA,NA) }
POTlist[[Nfile]] <- POT
#### geet responses ----
rm(ALLE3)
if (all(is.na(resIND[3,]))) ProdsinPDF <- 1:2 else if (all(is.na(resIND[2,]))) ProdsinPDF <- c(1,3) else ProdsinPDF <- 1:3
if (all(is.na(resIND[2:3,]))) next
for (x in ProdsinPDF) {
ZeilenR <- strsplit(tt[resIND[x,1]], "\n")[[1]]
ZeilenREF <- ZeilenR[resIND[x,2]:length(ZeilenR)]
grepDoses <- grep("Dose .alue", ZeilenREF)
Doses <- unlist(regmatches(c(ZeilenREF[grepDoses[1]], ZeilenREF[grepDoses[2]]),
gregexpr("([0-9]+)\\.?([0-9]+)",
c(ZeilenREF[grepDoses[1]], ZeilenREF[grepDoses[2]]))))
uniDoses <- as.numeric(unique(Doses))
Ndoses <- length(uniDoses)
grepResponses <- grep(".esponse", ZeilenREF)
grepMean <- grep(".ean", ZeilenREF)
NoResp <- grepMean[1]-grepResponses[1]
IND <- c(grepDoses)
StepIX <- IND+1
INDresp1 <- grepResponses
INDresp2 <- grepMean-1
Responses <- list()
for (n in 1:length(IND)) {
Matches <- unlist(regmatches(ZeilenREF[INDresp1[n]:INDresp2[n]],
gregexpr("([0-9]+)\\.?([0-9]+)",
ZeilenREF[INDresp1[n]:INDresp2[n]])))
Responses[[n]] <- as.numeric(Matches)
}
Step <- unlist(regmatches(ZeilenREF[StepIX[1]], gregexpr("([0-9]+)", ZeilenREF[StepIX[1]])))
maxStep <- max(as.numeric(Step))
resStep <- Ndoses-maxStep
if (maxStep== (Ndoses-1)) {
ref_data <- data.frame(c(Responses[[1]][1:maxStep], Responses[[2]][1]),
c(Responses[[1]][(maxStep+1):(maxStep*2)], Responses[[2]][2]),
c(Responses[[1]][(maxStep*2+1):(maxStep*3)], Responses[[2]][3]))
} else if (maxStep< (Ndoses-1)) {
ref_data <- data.frame(c(Responses[[1]][1:maxStep], Responses[[2]][1:resStep]),
c(Responses[[1]][(maxStep+1):(maxStep*2)], Responses[[2]][(resStep+1):(resStep*2)]),
c(Responses[[1]][(maxStep*2+1):(maxStep*3)], Responses[[2]][(resStep*2+1):(resStep*3)]))
} else {
ref_data <- data.frame(c(Responses[[1]][1:maxStep]),
c(Responses[[1]][(maxStep+1):(maxStep*2)]),
c(Responses[[1]][(maxStep*2+1):(maxStep*3)]))
}
if (x==1) colnames(ref_data) <- c(paste0(REF, "_1"), paste0(REF, "_2"), paste0(REF, "_3"))
if (x==2) colnames(ref_data) <- c(paste0(TEST1, "_1"), paste0(TEST1, "_2"), paste0(TEST1, "_3"))
if (x==3) colnames(ref_data) <- c(paste0(TEST2, "_1"), paste0(TEST2, "_2"), paste0(TEST2, "_3"))
if (!exists("ALLE3")) {
ALLE3 <- ref_data
} else { ALLE3 <- cbind(ALLE3, ref_data) }
} # loop ProdsinPDF
ALLE3_ <- Filter(function(x) !all(is.na(x)), ALLE3)
ALLE3 <- cbind(uniDoses, ALLE3_)
FINALresponses[[Nfile]] <- ALLE3
SampleNamesL[[Nfile]] <- c(REF, TEST1, TEST2)
DilSerVec <- c(DilSerVec, NoResp)
} # loop Nfile
Dat$FINALreponses <- FINALresponses
Dat$POTlist <- POTlist
Dat$DilSerVec <- DilSerVec
Dat$SampleNamesL <- SampleNamesL
} # if (iFile)
})
observe({
if (!is.null(Dat$FINALresponses)) {
filenames <- Dat$FNames
FINALresponses <- Dat$FINALresponses
POTlist <- Dat$POTlist
DilSerVec <- Dat$DilSerVec
rm(SIGrefDF, SIGtestDF)
CalcPot <- RestrM <- UnrestrM <- list()
for (pdfInd in 1:length(FINALresponses)) {
Potratio <- POTlist[[pdfInd]]
getPotratio <- unlist(regmatches(Potratio[1],
gregexpr("([0-9]+)\\.?([0-9]+)",
Potratio[1])))
Alle3 <- FINALresponses[[pdfInd]]
if (is.null(Alle3)) next
if (DilSerVec[pdfInd] == 3) {
StS <- Alle3[,1:4]
TestSs <- Alle3[,c(5:length(Alle3))]
StTest1 <- cbind(StS, TestSs[,1:3])
if (ncol(TestSs) >3) StTest2 <- cbind(StS, TestSs[,4:6])
} else if (DilSerVec[pdfInd] == 2) {
StS <- Alle3[,1:3]
TestSs <- Alle3[,c(4:length(Alle3))]
StTest1 <- cbind(StS, TestSs[,1:2])
if (ncol(TestSs) >2) StTest2 <- cbind(StS, TestSs[,3:4])
}
PotM1_1 <- RestrMF(StTest1, M="CI")
if (ncol(TestSs)/DilSerVec[pdfInd] ==2) {
PotM1_2 <- RestrMF(StTest2, M="CI")
PotM1 <- data.frame(rbind(unlist(PotM1_1), unlist(PotM1_2)))
} else { PotM1 <- data.frame(t(PotM1_1)) }
colnames(PotM1) <- c("rel. potency", "lowerCI", "upperCI","lowerCIpureEr", "upperCIpureEr")
PotM1_ <- as.matrix(PotM1)
colnames(PotM1_) <- names(PotM1)
PotM1 <- cbind(rep(filenames[pdfInd], nrow(PotM1_)), PotM1_)
CalcPot[[pdfInd]] <- PotM1_
if (!exists("CalcPotDF")) CalcPotDF <- PotM1_ else CalcPotDF <- rbind(CalcPotDF, PotM1_)
RMcoefs1 <- RestrMF(StTest1, M="R")
if (ncol(TestSs)/DilSerVec[pdfInd] ==2) {
RMcoefs2 <- RestrMF(StTest2, M="R")
RMcoefs <- rbind(c(basename(filenames[pdfInd]), NA, NA, NA), RMcoefs1, rep(NA,4), RMcoefs2)
} else { RMcoefs <- rbind(c(basename(filenames[pdfInd]), NA, NA, NA), RMcoefs1) }
RestrM[[pdfInd]] <- RMcoefs
if (!exists("RMcoefsDF")) RMcoefsDF <- RMcoefs[,1] else RMcoefsDF <- rbind(RMcoefsDF, RMcoefs[,1])
URMcoefs1 <- RestrMF(StTest1, M="U")
if (ncol(TestSs)/DilSerVec[pdfInd] ==2) {
URMcoefs2 <- RestrMF(StTest2, M="U")
URMcoefs <- rbind(c(basename(filenames[pdfInd]), NA, NA, NA), URMcoefs1, rep(NA,4), URMcoefs2)
} else { URMcoefs <- rbind(c(basename(filenames[pdfInd]), NA, NA, NA), URMcoefs1) }
UnrestrM[[pdfInd]] <- URMcoefs
if (!exists("URMcoefsDF")) URMcoefsDF <- URMcoefs[,1] else URMcoefsDF <- rbind(URMcoefsDF, URMcoefs[,1])
X <- seq(min(log(Alle3$uniDoses)), max(log(Alle3$uniDoses)), 0.1)
sigRef <- URMcoefs1[1,1] + (URMcoefs1[4,1]-URMcoefs1[1,1])/(1+exp(URMcoefs1[2,1]*(URMcoefs1[3,1]-X)))
sigTest1 <- URMcoefs1[5,1] + (URMcoefs1[8,1]-URMcoefs1[5,1])/(1+exp(URMcoefs1[6,1]*(URMcoefs1[7,1]-X)))
if (exists("URMcoefs2")) {
if (!all(is.na(URMcoefs2)) & sum(URMcoefs2) != 2) {
sigTest1 <- URMcoefs2[5,1] + (URMcoefs2[8,1]-URMcoefs2[5,1])/(1+exp(URMcoefs2[6,1]*(URMcoefs2[7,1]-X)))
}
}
dfPlotsigRef <- data.frame(X=X, sigRef = sigRef, Prod = pdfInd)
dfPlotsigTest1 <- data.frame(X=X, sigTest = sigTest1, Prod = paste0(pdfInd, "_1"))
if (exists("URMcoefs2")) {
if (!all(is.na(URMcoefs2))) {
dfPlotsigTest2 <- data.frame(X=X, sigTest = sigTest2, Prod = paste0(pdfInd, "_2"))
}
}
if (!exists("SIGrefDF")) SIGrefDF <- dfPlotsigRef else SIGrefDF <- rbind(SIGrefDF, dfPlotsigRef)
if (!exists("SIGtestDF")) {
if (exists("dfPlotsigTest2")) SIGtestDF <- rbind(dfPlotsigTest1, dfPlotsigTest2) else SIGtestDF <- dfPlotsigTest1
} else {
if (exists("dfPlotsigTest2")) SIGtestDF <- rbind(SIGtestDF,dfPlotsigTest1, dfPlotsigTest2) else SIGtestDF <- rbind(dfPlotsigTest1, dfPlotsigTest2)
}
} # pdfInd
if (ncol(URMcoefsDF)>11) EC50TEST <- as.numeric(c(URMcoefsDF[,8], URMcoefsDF[,17])) else EC50TEST <- as.numeric(c(URMcoefsDF[,8]))
EC50TEST <- EC50TEST[!EC50TEST %in% boxplot.stats(EC50TEST)$out]
EC50REF <- as.numeric(URMcoefsDF[,4])
EC50REF <- EC50REF[!EC50REF %in% boxplot.stats(EC50REF)$out]
UasREF <- as.numeric(URMcoefsDF[,5])
UasREF <- UasREF[!UasREF %in% boxplot.stats(UasREF)$out]
LasREF <- as.numeric(URMcoefsDF[,2])
LasREF <- LasREF[!LasREF %in% boxplot.stats(LasREF)$out]
Dat$URMcoefsDF <- URMcoefsDF
Dat$RestrM <- RestrM
Dat$CalcPot <- CalcPot
#### sigmoid plots ----
Slope <- as.numeric(URMcoefsDF[1,3])
if (Slope > 0) {
x_UA <- max(X); x_LA <- min(X)
} else { x_UA <- min(X); x_LA <- max(X) }
p1 <- ggplot(SIGrefDF, aes(x_X, y=sigRef, col=as.factor(Prod))) +
geom_line() +
annotate("text", label="x", x=x_UA, y=UasREF, alpha=0.2) +
annotate("text", label="o", x=x_LA, y=LasREF, alpha=0.2) +
geom_vline(xintercept = EC50REF, alpha = 0.2) +
xlab("dilutions") +
ggtitle("Plot of all calculated reference fits (unrestricted model, in gray vertical lines: EC50)") +
theme_bw() +
theme(axis.text = element_text(face = "bold", size = 15),
plot.title = element_text(size = 15, face = "bold"))
output$sigPlotREF <- renderPlot({ p1 })
PLOTS$sigPlotREF <- p1
p2 <- ggplot(SIGtestDF, aes(x_X, y=sigTest, col=as.factor(Prod))) +
geom_line() +
#annotate("text", label="x", x=x_UA, y=UasREF, alpha=0.2) +
#annotate("text", label="o", x=x_LA, y=LasREF, alpha=0.2) +
geom_vline(xintercept = EC50TEST, alpha = 0.2) +
xlab("dilutions") +
ggtitle("Plot of all calculated reference fits (unrestricted model, in gray vertical lines: EC50)") +
theme_bw() +
theme(axis.text = element_text(face = "bold", size = 15),
plot.title = element_text(size = 15, face = "bold"))
output$sigPlotREF <- renderPlot({ p2 })
PLOTS$sigPlotTEST <- p2
#### histograms ----
output$histEC50REF <- renderPlot({
hist(EC50REF, main = 'Histogram of EC50REF')
})
output$histLasREF <- renderPlot({
hist(Las0REF, main = 'Histogram of lower asymptotes REF')
})
output$histUasREF <- renderPlot({
hist(UasREF, main = 'Histogram of upper asymptotes REF')
})
PLOTS$histEO50REF <- hist(EC50REF, main = 'Histogram of EC50REF')
PLOTS$histLasREF <- hist(LasREF, main = 'Histogram of EC50REF')
PLOTS$histUasREF <- hist(UasREF, main = 'Histogram of EC50REF')
#### Calculated Potencies table ----
colnames(CalcPotDF) <- c("pdf","rel. potency","lowerCI","upperCI","lower CI puerErr","upper CI pureErr")
SampleNamesL <- Dat$SampleNamesL
SampleNames <- unlist(SampleNamesL)
if (length(SampleNames)/length(SampleNamesL) == 3) {
REFIndices <- seq(1,length(SampleNamesL))*3-2
} else { REFIndices <- seq(1,length(SampleNamesL))*2-1 }
Test_Sample_Names <- SampleNames[-c(REFIndices)]
CalcPotDF <- as.data.frame(CalcPotDF)
CalcPotDF[,2:6] <- apply(CalcPotDF[,2:6],2,as.numeric)
CalcPotDF$relLCI <- CalcPotDF$lowerCI/CalcPotDF$`rel. potency`
CalcPotDF$relUCI <- CalcPotDF$upperCI/CalcPotDF$`rel. potency`
CalcPotDF$`relLCI pureErr` <- CalcPotDF$`lowerCI pureErr`/CalcPotDF$`rel. potency`
CalcPotDF$`relUCI pureErr` <- CalcPotDF$`upperCI pureErr`/CalcPotDF$`rel. potency`
CalcPotDF2 <- cbind(CalcPotDF, Test_Sample_Names)
Dat$CalcPotDF2 <- CalcPotDF2
output$CalcPotDF <- renderTable({ CalcPotDF2 })
#### CI plots ----
CalcPotDF_ <- CalcPotDF[,c(1,7,8)]
all_lPot <- melt(as.data.frame(CalcPotDF_, id.vars="pdf",variable.name="var", value.name="Climit"))
all_lPot[,3] <- as.numeric(all_lPot[,3])
all_lPot[,3][all_lPot[,3] > 5] <- NA
all_lPot[,3][all_lPot[,3] < 0.1] <- NA
P_histCI <- ggplot(all_lPot, aes(x=Climit, fill=var)) +
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$Climit, na.rm=T)*10)/10, max(all_lPot$Climit, na.rm=T)*1.1, by=0.2),
) +
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$histCIsPure <- renderPlot({ P_histCI })
CalcPotDF_Pure <- CalcPotDF[,c(1,9,10)]
all_lPotPure <- melt(as.data.frame(CalcPotDF_Pure, id.vars="pdf",variable.name="var", value.name="Climit"))
all_lPotPure[,3][all_lPotPure[,3] == 0] <- NA
all_lPotPure[,3][all_lPotPure[,3] > 5] <- NA
all_lPotPure[,3][all_lPotPure[,3] < 0.1] <- NA
P_histCIPure <- ggplot(all_lPotPure, aes(x=Climit, fill=var)) +
geom_histogram(color="#e9ecef", alpha=0.6, position = "identity") +
scale-fill-manual(values=c("darkgreen","darkblue","salmon2","tomato3")) +
ggtitle("Histogram of relative potencies, pure error") +
scale_x_continuous(
breaks=seq(trunc(min(all_lPotPure$Climit, na.rm=T)*10)/10, max(all_lPotPure$Climit, na.rm=T)*1.1, by=0.2),
) +
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$histCIsPure <- renderPlot({ P_histCIPure })
ListofPDFs <- sapply(RestrM, function(x) rbind(x[1]))
output$ListofPDFs <- renderTable({ ListofPDFs}, rownames=T)
logCPDF <- log(CalcPotDF[,9:10])
logCPDF1 <- cbind(CalcPotDF[,1], logCPDF)
logCPDF1 <- logCPDF1[!duplicated(logCPDF1[,2:3]),]
logCPDF_ <- logCPDF1[!is.na(logCPDF1[,3]),]
colnames(logCPDF_) <- c("pdf", "relLCI","relUCI")
dg <- function(x,loc,scale,log) {
r <- try(dgumbel(x,loc,scale,log), silent=T)
if (inherits(r, "try-error")) return(NA)
return(r)
}
uCIGumbel <- fitdistr(logCPDF_[,3], dg,start=list(loc=1,scale=1))
lCIGumbel <- fitdistr(logCPDF_[,2], dg,start=list(loc=1,scale=1))
Gumbel <- c(rev(dgumbel(seq(0,0.5,0.01), loc=uCIGumbel$extimate[1], scale=uCIGumbel$extimate[2])), dgumbel(seq(0,0.5,0.01), loc=uCIGumbel$extimate[1], scale=uCIGumbel$extimate[2]))
GumbelDF <- as.data.frame(cbind(X=seq(-0.5,0.5,0.01), Gumbel))
all_lPotPure <- melt(as.data.frame(logCPDF_), id.vars="pdf",variable.name = "var",value.name="Climit")
all_lPotPure[,3][all_lPotPure[,3] == 0] <- NA
uEAC099 <- exp(qgumbel(0.99, uCIGumbel$estimate[1],uCIGumbel$estimate[2]))*100
lEAC001 <- exp(-qgumbel(0.99, uCIGumbel$estimate[1],uCIGumbel$estimate[2]))*100
uEAC095 <- exp(qgumbel(0.95, uCIGumbel$estimate[1],uCIGumbel$estimate[2]))*100
lEAC005 <- exp(-qgumbel(0.95, uCIGumbel$estimate[1],uCIGumbel$estimate[2]))*100
output$GumbelPlot <- renderPlot({
hist(all_lPotPure$Climit, breaks = 50, xlab="log(relative CLs)",
main=paste("calculated EACs 99%:",round(lEAC001,0)), "-", round(uEAC099,0), " 95%", round(lEAC005,0), "-",round(uEAC095,0))
lines(GumbelDF$X, GumbelDF$Gumbel)
})
}
})
output$downloadPLA <- downloadHandler(
filename="PLAcontent.zip",
content=function(fname) {
fs <- c()
tmpdir <- tempdir()
fileOutUNRM=paste(paste0(tmpdir,sep='/', "Coefficients_UnrestrModel"), sep='','.xlsx')
fs=c(fs, fileOutUNRM)
URMcoelsDF <- Dat$URMcoefsDF
write.xlsx(URMcoefsDF, fileOutUNRM, rowNames=T)
fileOutPOT=paste(paste0(tmpdir,sep='/', "relPot_CIs"), sep='','.xlsx')
fs=c(fs, fileOutPOT)
CalcPotDF <- Dat$CalcPotDF
write.xlsx(CalcPotDF, fileOutPOT, rowNames=T)
sigPlotREF <- PLOTS$sigPlotREF
sigPlotTEST <- PLOTS$sigPlotTEST
fileOutSigREF=paste(paste0(tmpdir,sep='/', "_SigmoidsREF"), sep='','.png')
fs=c(fs, fileOutSigREF)
png(fileOutSigREF, width=600, height=400)
plot(sigPlotREF)
dev.off()
fileOutSigTEST=paste(paste0(tmpdir,sep='/', "_SigmoidsTEST"), sep='','.png')
fs=c(fs, fileOutSigTEST)
png(fileOutSigTEST, width=600, height=400)
plot(sigPlotTEST)
dev.off()
fileOutRM=paste(paste0(tmpdir,sep='/', "_RestrictedModel"), sep='','.xlsx')
fileOutCalcPot=paste(paste0(tmpdir,sep='/', "_CalculatedPotencies"), sep='','.xlsx')
fileOutResponses=paste(paste0(tmpdir,sep='/', "_Responses"), sep='','.xlsx')
fs=c(fs, fileOutRM,fileOutCalcPot,fileOutResponses)
RestrM <- Dat$RestrM
CalcPot <- Dat$CalcPot
FINALresponses <- Dat$FINALresponses
write.xlsx(RestrM, fileOutRM, rowNames=T)
write.xlsx(CalcPot, fileOutCalcPot, rowNames=F)
write.xlsx(FINALresponses, fileOutResponses, rowNames=F)
zip::zipr(zipfile=fname, fileas=fs, include_directories = F)
}, contentType = "application/zip"
)
}
shinyApp(ui=ui, server=server)