Commit 851e0209 authored by GAREL Marc's avatar GAREL Marc
Browse files

For xmid et Asym, round git =3

parent 03e53da2
##################################################
## Project: Shiny app allowing plot and parameters for microbial growth dataset
## Script purpose: Perform logistic regression to estimate growth rate and maximum cells density.
## Date: October 2017
## Author: Marc Garel, Severine Martini, Christian Tamburini
##################################################
### libraries
library(shiny)
library(DT)
library(rgl)
library(quantreg)
library(SparseM)
library(rmarkdown)
# server script
shinyServer(
function(input, output, session){
####### Upload data #######
getData <- reactive({
inFile <- input$file1
if (is.null(inFile)) return(NULL)
df <- read.csv(inFile$datapath, header = input$header, sep = input$sep, dec = input$dec)
return(df)
})
####### Calcul des estimateurs #######
calculate <- reactive({
df <- getData()
req(df)
# calcul coef pour affichage dans en dessous du graph
x <- rep(df[,1], ncol(df)-1)
y <- unlist(df[, 2:ncol(df)])
crssce <- data.frame(x, y)
fit <- nls(y ~ SSlogis(x, Asym, xmid, scal), crssce)
param <- summary(fit)$parameters
# cat(param,'\n\n')
scal <- param[3,1]
sdscal <- param[3,2]
mu <- round(1 / scal, digits = 3)
sdmu <- round((1 / scal) - (1 / (scal + sdscal)), digits = 3)
asym <- param[1,1]
sdasym <- param[1,2]
xmid <- param[2,1]
sdxmid <- param[2,2]
params <- as.data.frame(cbind(scal, sdscal, mu, sdmu, asym, sdasym, xmid, sdxmid))
colnames(params) <- c("scal", "sdscal", "Growth rate (Gr)", "Gr Standard deviation", "Maximum population (Mp)", "Mp Standard deviation",
"Start expo. phase (Ste)", "Ste Standard deviation")
return(params)
})
####### Raw_plot #######
output$raw_plot <- renderPlot({
df <- getData()
req(df)
x <- rep(df[, 1], ncol(df) - 1)
y <- unlist(df[, 2:ncol(df)])
crssce <- data.frame(y, x)
# limit axes extended
maxx <- max(x, na.rm = T) + (max(x, na.rm = T) / 4)
maxy <- max(y, na.rm = T) + (max(y, na.rm = T) / 4)
w0 <- seq(0, maxx)
coef <- getInitial(y ~ SSlogis(x, asym, xmid, scal), data = cbind.data.frame(x, y))
mu <- 1 / coef[3]
# output plot
plot(x, y, xlab = input$xlabel, ylab = input$ylabel, xlim = c (0, maxx),
ylim = c(0, maxy), las = 1, pch = 22, bg = 1)
lines(w0, SSlogis(w0, coef[1], coef[2], coef[3]), lwd = 1.5)
grid()
q1b <- nlrq(y ~ SSlogis(x, K, xmid, r), tau = 0.25, data = crssce)
q3b <- nlrq(y ~ SSlogis(x, K, xmid, r), tau = 0.75, data = crssce)
q05b <- nlrq(y ~ SSlogis(x, K, xmid, r), tau = 0.025, data = crssce)
q95b <- nlrq(y ~ SSlogis(x, K, xmid, r), tau = 0.975, data = crssce)
lines(w0, predict(q1b, newdata = list(x = w0)), col = "blue", lty=2)
lines(w0, predict(q3b, newdata = list(x = w0)), col = "blue", lty=2)
lines(w0, predict(q05b, newdata = list(x = w0)), col = "blue", lty=3)
lines(w0, predict(q95b, newdata = list(x = w0)), col = "blue", lty=3)
#
})
####### Affichage du tableau des parametres #######
output$parameters <- DT::renderDataTable({
params <- calculate()
paste("The table below gives the parameters estimation.")
DT::datatable(params[,3:8])
#
})
####### Affichage du tableau initial #######
output$raw_data <- DT::renderDataTable({
df <- getData()
DT::datatable(df)
#
})
####### Affichage de la citation #######
output$citation <- renderText({
paste("Martini, S., B. Al Ali, M. Garel, and others. 2013. Effects of Hydrostatic Pressure on Growth and Luminescence of a Moderately-Piezophilic Luminous Bacteria\nPhotobacterium phosphoreum ANT-2200 A. Driks [ed.]. PLoS One 8: e66580. doi:10.1371/journal.pone.0066580")
#
})
####### Fonction de cout #######
cost <- function(Asym, scal, xmid, xdata, ydata){
ymod <- SSlogis(xdata, Asym, xmid, scal)
SSE <- sum((ymod - ydata)^2, na.rm = T)
return(SSE)
}
####### Verify_plot #######
output$verify_plot_2D <- renderPlot({
# on recupereles donnes de l'utilisateur
df <- getData()
req(df)
# on calcule les parametres
params <- calculate()
req(params)
# on assigne les valeurs de params a nos differentes variables
scal <- params[1,1]
mu <- params[1,3]
Asym <- params[1,5]
xmid <- params[1,7]
# cat(Asym, '\n\n')
# cat(mu, '\n\n')
# cat(scal, '\n\n')
# cat(xmid, '\n\n')
x <- rep(df[, 1], ncol(df) - 1)
y <- unlist(df[, 2:ncol(df)])
Asymmin <- Asym - (5*Asym) / 100
Asymmax <- Asym + (5*Asym) / 100
scal_min <- scal - (5*scal) / 100
scal_max <- scal + (5*scal) / 100
mu_min <- 1/scal_max
mu_max <- 1/scal_min
# cat(Asymmin, '\n\n')
# cat(Asymmax, '\n\n')
#
# cat(scal_min, '\n\n')
# cat(scal_max, '\n\n')
#
# cat(mu_min, '\n\n')
# cat(mu_max, '\n\n')
SSE <- cost(Asym, scal, xmid, x, y)
# cat(SSE, '\n\n')
scal_Vector <- seq(scal_min, scal_max, length=100)
mu_Vector <- seq(mu_min, mu_max, length=100)
Asym_vector <- seq(Asymmin, Asymmax, length=100)
costV <- Vectorize(cost, c("scal", "Asym"))
z <- outer(Asym_vector, scal_Vector, FUN = "costV", xmid = xmid, xdata = x, ydata = y)
# print(costV, '\n\n')
# cat(Asym, '\n\n')
# cat(mu, '\n\n')
# cat(head(z), '\n\n')
contour(Asym_vector, mu_Vector, z, nlevels = 30, xlab = "K (en cellules/ml)", ylab = "Growth rate (h-1)",
xlim = range(Asym_vector, finite = TRUE),
ylim = range(mu_Vector, finite = TRUE),
zlim = range(z, finite = TRUE))
points(Asym, mu, pch = 19, col = "lightseagreen")
})
####### Verify plot 3D #######
output$verify_plot_3D <- renderRglwidget({
df <- getData()
req(df)
params <- calculate()
req(params)
scal <- params[1,1]
mu <- params[1,3]
Asym <- params[1,5]
xmid <- params[1,7]
x <- rep(df[, 1], ncol(df) - 1)
y <- unlist(df[, 2:ncol(df)])
Asymmin <- Asym - (5*Asym) / 100
Asymmax <- Asym + (5*Asym) / 100
scal_min <- scal - (5*scal) / 100
scal_max <- scal + (5*scal) / 100
mu_min <- 1/scal_max
mu_max <- 1/scal_min
SSE <- cost(Asym, scal, xmid, x, y)
scal_Vector <- seq(scal_min, scal_max, length=100)
mu_Vector <- seq(mu_min, mu_max, length=100)
Asym_vector <- seq(Asymmin, Asymmax, length=100)
costV <- Vectorize(cost, c("scal", "Asym"))
z <- outer(Asym_vector, scal_Vector, FUN = "costV", xmid = xmid, xdata = x, ydata = y)
rgl.open(useNULL=T)
bg3d("white")
view3d( theta = 0, phi = -45)
plot3d(rep(Asym, length(Asym)),
rep(mu, each = length(mu)), z,
xlim = range(Asym_vector, finite = TRUE),
ylim = range(mu_Vector, finite = TRUE),
zlim = range(z, finite = TRUE),
type = "n", xlab = "Maximum population", ylab = "Growth rate", zlab = "SSE")
surface3d(Asym_vector, mu_Vector, z, col = "lightseagreen")
points3d(Asym, mu, SSE, size = 5, col = "black")
rglwidget()
})
output$downloadData <- downloadHandler(
filename = function() {
paste('data-', Sys.Date(), '.txt', sep='')
},
content = function(con) {
param <- calculate()
write.table((param[,3:8]), con, row.names = F, col.names = T, sep = '\t')
}
)
##################################################
## Project: Shiny app allowing plot and parameters for microbial growth dataset
## Script purpose: Perform logistic regression to estimate growth rate and maximum cells density.
## Date: October 2017
## Author: Marc Garel, Severine Martini, Christian Tamburini
##################################################
### libraries
library(shiny)
library(DT)
library(rgl)
library(quantreg)
library(SparseM)
library(rmarkdown)
# server script
shinyServer(
function(input, output, session){
####### Upload data #######
getData <- reactive({
inFile <- input$file1
if (is.null(inFile)) return(NULL)
df <- read.csv(inFile$datapath, header = input$header, sep = input$sep, dec = input$dec)
return(df)
})
####### Calcul des estimateurs #######
calculate <- reactive({
df <- getData()
req(df)
# calcul coef pour affichage dans en dessous du graph
x <- rep(df[,1], ncol(df)-1)
y <- unlist(df[, 2:ncol(df)])
crssce <- data.frame(x, y)
fit <- nls(y ~ SSlogis(x, Asym, xmid, scal), crssce)
param <- summary(fit)$parameters
# cat(param,'\n\n')
scal <- param[3,1]
sdscal <- param[3,2]
mu <- round(1 / scal, digits = 3)
sdmu <- round((1 / scal) - (1 / (scal + sdscal)), digits = 3)
asym <- round(param[1,1], digit=3)
sdasym <- round(param[1,2], digit=3)
xmid <- round(param[2,1], digit=3)
sdxmid <- round(param[2,2], digit=3)
params <- as.data.frame(cbind(scal, sdscal, mu, sdmu, asym, sdasym, xmid, sdxmid))
colnames(params) <- c("scal", "sdscal", "Growth rate (Gr)", "Gr Standard deviation", "Maximum population (Mp)", "Mp Standard deviation",
"Start expo. phase (Ste)", "Ste Standard deviation")
return(params)
})
####### Raw_plot #######
output$raw_plot <- renderPlot({
df <- getData()
req(df)
x <- rep(df[, 1], ncol(df) - 1)
y <- unlist(df[, 2:ncol(df)])
crssce <- data.frame(y, x)
# limit axes extended
maxx <- max(x, na.rm = T) + (max(x, na.rm = T) / 4)
maxy <- max(y, na.rm = T) + (max(y, na.rm = T) / 4)
w0 <- seq(0, maxx)
coef <- getInitial(y ~ SSlogis(x, asym, xmid, scal), data = cbind.data.frame(x, y))
mu <- 1 / coef[3]
# output plot
plot(x, y, xlab = input$xlabel, ylab = input$ylabel, xlim = c (0, maxx),
ylim = c(0, maxy), las = 1, pch = 22, bg = 1)
lines(w0, SSlogis(w0, coef[1], coef[2], coef[3]), lwd = 1.5)
grid()
q1b <- nlrq(y ~ SSlogis(x, K, xmid, r), tau = 0.25, data = crssce)
q3b <- nlrq(y ~ SSlogis(x, K, xmid, r), tau = 0.75, data = crssce)
q05b <- nlrq(y ~ SSlogis(x, K, xmid, r), tau = 0.025, data = crssce)
q95b <- nlrq(y ~ SSlogis(x, K, xmid, r), tau = 0.975, data = crssce)
lines(w0, predict(q1b, newdata = list(x = w0)), col = "blue", lty=2)
lines(w0, predict(q3b, newdata = list(x = w0)), col = "blue", lty=2)
lines(w0, predict(q05b, newdata = list(x = w0)), col = "blue", lty=3)
lines(w0, predict(q95b, newdata = list(x = w0)), col = "blue", lty=3)
#
})
####### Affichage du tableau des parametres #######
output$parameters <- DT::renderDataTable({
params <- calculate()
paste("The table below gives the parameters estimation.")
DT::datatable(params[,3:8])
#
})
####### Affichage du tableau initial #######
output$raw_data <- DT::renderDataTable({
df <- getData()
DT::datatable(df)
#
})
####### Affichage de la citation #######
output$citation <- renderText({
paste("Martini, S., B. Al Ali, M. Garel, and others. 2013. Effects of Hydrostatic Pressure on Growth and Luminescence of a Moderately-Piezophilic Luminous Bacteria\nPhotobacterium phosphoreum ANT-2200 A. Driks [ed.]. PLoS One 8: e66580. doi:10.1371/journal.pone.0066580")
#
})
####### Fonction de cout #######
cost <- function(Asym, scal, xmid, xdata, ydata){
ymod <- SSlogis(xdata, Asym, xmid, scal)
SSE <- sum((ymod - ydata)^2, na.rm = T)
return(SSE)
}
####### Verify_plot #######
output$verify_plot_2D <- renderPlot({
# on recupereles donnes de l'utilisateur
df <- getData()
req(df)
# on calcule les parametres
params <- calculate()
req(params)
# on assigne les valeurs de params a nos differentes variables
scal <- params[1,1]
mu <- params[1,3]
Asym <- params[1,5]
xmid <- params[1,7]
# cat(Asym, '\n\n')
# cat(mu, '\n\n')
# cat(scal, '\n\n')
# cat(xmid, '\n\n')
x <- rep(df[, 1], ncol(df) - 1)
y <- unlist(df[, 2:ncol(df)])
Asymmin <- Asym - (5*Asym) / 100
Asymmax <- Asym + (5*Asym) / 100
scal_min <- scal - (5*scal) / 100
scal_max <- scal + (5*scal) / 100
mu_min <- 1/scal_max
mu_max <- 1/scal_min
# cat(Asymmin, '\n\n')
# cat(Asymmax, '\n\n')
#
# cat(scal_min, '\n\n')
# cat(scal_max, '\n\n')
#
# cat(mu_min, '\n\n')
# cat(mu_max, '\n\n')
SSE <- cost(Asym, scal, xmid, x, y)
# cat(SSE, '\n\n')
scal_Vector <- seq(scal_min, scal_max, length=100)
mu_Vector <- seq(mu_min, mu_max, length=100)
Asym_vector <- seq(Asymmin, Asymmax, length=100)
costV <- Vectorize(cost, c("scal", "Asym"))
z <- outer(Asym_vector, scal_Vector, FUN = "costV", xmid = xmid, xdata = x, ydata = y)
# print(costV, '\n\n')
# cat(Asym, '\n\n')
# cat(mu, '\n\n')
# cat(head(z), '\n\n')
contour(Asym_vector, mu_Vector, z, nlevels = 30, xlab = "K (en cellules/ml)", ylab = "Growth rate (h-1)",
xlim = range(Asym_vector, finite = TRUE),
ylim = range(mu_Vector, finite = TRUE),
zlim = range(z, finite = TRUE))
points(Asym, mu, pch = 19, col = "lightseagreen")
})
####### Verify plot 3D #######
output$verify_plot_3D <- renderRglwidget({
df <- getData()
req(df)
params <- calculate()
req(params)
scal <- params[1,1]
mu <- params[1,3]
Asym <- params[1,5]
xmid <- params[1,7]
x <- rep(df[, 1], ncol(df) - 1)
y <- unlist(df[, 2:ncol(df)])
Asymmin <- Asym - (5*Asym) / 100
Asymmax <- Asym + (5*Asym) / 100
scal_min <- scal - (5*scal) / 100
scal_max <- scal + (5*scal) / 100
mu_min <- 1/scal_max
mu_max <- 1/scal_min
SSE <- cost(Asym, scal, xmid, x, y)
scal_Vector <- seq(scal_min, scal_max, length=100)
mu_Vector <- seq(mu_min, mu_max, length=100)
Asym_vector <- seq(Asymmin, Asymmax, length=100)
costV <- Vectorize(cost, c("scal", "Asym"))
z <- outer(Asym_vector, scal_Vector, FUN = "costV", xmid = xmid, xdata = x, ydata = y)
rgl.open(useNULL=T)
bg3d("white")
view3d( theta = 0, phi = -45)
plot3d(rep(Asym, length(Asym)),
rep(mu, each = length(mu)), z,
xlim = range(Asym_vector, finite = TRUE),
ylim = range(mu_Vector, finite = TRUE),
zlim = range(z, finite = TRUE),
type = "n", xlab = "Maximum population", ylab = "Growth rate", zlab = "SSE")
surface3d(Asym_vector, mu_Vector, z, col = "lightseagreen")
points3d(Asym, mu, SSE, size = 5, col = "black")
rglwidget()
})
output$downloadData <- downloadHandler(
filename = function() {
paste('data-', Sys.Date(), '.txt', sep='')
},
content = function(con) {
param <- calculate()
write.table((param[,3:8]), con, row.names = F, col.names = T, sep = '\t')
}
)
})
\ No newline at end of file
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment