route of administration

absorption

delay

distribution

elimination

output

pkmodelCode.R

ui.R

Download

library(shinydashboard)
library(mlxR)
# source("C:/Users/Marc/mlxtoolbox/mlxR220/initMlxR220.R")

sidebar <- dashboardSidebar(
  hr(),
  sidebarMenu(
    menuItem("Plot", tabName = "plot", icon = icon("line-chart")),
    menuItem("Codes",  icon = icon("file-text-o"),
             menuSubItem("pkmodelCode.R", tabName = "pkmodel", icon = icon("angle-right")),
             menuSubItem("ui.R", tabName = "ui", icon = icon("angle-right")),
             menuSubItem("server.R", tabName = "server", icon = icon("angle-right"))
    ),
    menuItem("ReadMe", tabName = "readme", icon = icon("mortar-board"))
    #     menuItem("a", tabName = "perso", icon = icon("user"))
  ),
  hr(),
  sidebarMenu(
    menuItem("Administration", tabName = "admin", icon = icon("chevron-circle-right")),
    menuItem("Parameterization", icon = icon("chevron-circle-right"),
             fluidRow(
               column(1),
               column(10,
             radioButtons("parameterization","",c("rate constant" = "1","clearance" = "2"))
   )) ),
    menuItem("Model", icon = icon("chevron-circle-right"),
             menuSubItem("Absorption", tabName = "absorption", icon = icon("angle-right")),
             menuSubItem("Distribution", tabName = "distribution", icon = icon("angle-right")),
             menuSubItem("Elimination", tabName = "elimination", icon = icon("angle-right"))
    ),
    menuItem("Output", tabName = "output", icon = icon("chevron-circle-right"))
  ),
  hr()
#   conditionalPanel(condition = "input.distribution==2 | input.distribution==3" ,
#                    checkboxInput("cmt1", "Compartment 1 (central)", TRUE),
#                    checkboxInput("cmt2", "Compartment 2 (peripheral)", FALSE)
#   ),
#   conditionalPanel(condition = "input.distribution==3" ,
#                    checkboxInput("cmt3", "Compartment 3 (peripheral)", FALSE)
#   )
)

body <- dashboardBody(
  tabItems(
    tabItem(tabName = "plot",
            fluidRow( 
              box(width = 12, status = "primary",
#                   verbatimTextOutput("dynamic_value")
                   plotOutput("plot",  height="400px")
              ),
              tabItems(
                tabItem(tabName = "admin",
                        box(width = 4, status = "primary", solidHeader = TRUE, title="route of administration",          
                            radioButtons("administration", "",
                                         c("iv bolus" = "bolus","iv infusion" = "infusion","oral" = "oral"),selected="oral")
                        ),
                        box(width = 4, status = "primary",          
                            sliderInput("tfd", "Time of first dose:", value=0, min=0, max = 20, step=1),
                            sliderInput("nd", "Number of doses:", value=4, min=0, max = 10, step=1)
                        ),
                        box(width = 4, status = "primary",          
                            sliderInput("ii", "Interdose interval:", value = 9, min = 0.5, max = 15, step=0.5),
                            sliderInput("amt", "Amount:", value = 5, min = 0, max = 20, step=1),
                            conditionalPanel(condition = "input.administration == 'infusion'",
                                             sliderInput("tinf", "Infusion time:", value = 1, min = 0, max = 5, step=0.2)
                            )
                        )
                ),
                tabItem(tabName = "absorption",
                        conditionalPanel(condition = "input.administration == 'oral'",
                                         box(width = 4, status = "primary", solidHeader = TRUE, title="absorption",       
                                             radioButtons("absorption","",c("first order" = "1","zero order" = "0")),
                                             
                                             conditionalPanel(condition = "input.absorption == '1'",
                                                              sliderInput("ka", "rate constant ka:", value = 0.8, min = 0, max = 4, step=0.1)
                                             ),
                                             
                                             conditionalPanel(condition = "input.absorption == '0'",
                                                              sliderInput("tk0", "duration Tk0:", value = 3, min = 0, max = 10, step=0.5)
                                             )
                                         ),
                                         box(width = 4, status = "primary", solidHeader = TRUE, title="delay",          
                                             conditionalPanel(condition = "input.absorption == '0'",
                                                              selectInput("delay0","",c("none" = "0","lag time" = "1"))
                                                              
                                             ),                            
                                             conditionalPanel(condition = "input.absorption == '1'",
                                                              selectInput("delay1","",c("none" = "0","lag time" = "1","transit compartment" = "2"))
                                             ),
                                             #                                              conditionalPanel(condition = "input.delay1 == '1' | input.delay0 == '1' ",
                                             conditionalPanel(condition = "(input.absorption == '1' & input.delay1 == '1') | (input.absorption == '0' & input.delay0 == '1')",
                                                              sliderInput("tlag", "lag time Tlag:", value = 1, min = 0, max = 10, step=0.5)
                                             ),
                                             
                                             conditionalPanel(condition = "input.absorption == '1' & input.delay1 == '2'",
                                                              sliderInput("mtt", "mean transit time Mtt:", value = 1, min = 0, max = 10, step=0.5),
                                                              sliderInput("ntr", "number of compartments Ntr:", value = 2, min = 0.5, max = 10, step=0.5)
                                             )
                                         )
                        )
                ),
                tabItem(tabName = "distribution",
                        #                         fluidRow( 
                        box(width = 4, status = "primary", solidHeader = TRUE, title="distribution",        
                            radioButtons("distribution","",c("1 compartment" = "1","2 compartments" = "2","3 compartments" = "3")),
                            sliderInput("v", "volume V:", value = 10, min = 1, max = 20, step=1)
                        ),
                        box(width = 4,   status = "primary",         
                            conditionalPanel(condition = "input.distribution != '1' & input.parameterization=='1'",
                                             sliderInput("k12", "transition rate constant k12:", value = 0.4, min = 0, max = 2, step=0.05),
                                             sliderInput("k21", "transition rate constant k21:", value = 0.2, min = 0, max = 2, step=0.05)
                            ),
                            conditionalPanel(condition = "input.distribution != '1' & input.parameterization=='2'",
                                             sliderInput("v2", "periph. comp. volume V2:", value = 5, min = 1, max = 10, step=0.5),
                                             sliderInput("q2", "inter comp. clearance Q2:", value = 0.5, min = 0.1, max = 2, step=0.1)
                            )
                        ),
                        box(width = 4,   status = "primary",         
                            conditionalPanel(condition = "input.distribution == '3' & input.parameterization=='1'",
                                             sliderInput("k13", "transition rate constant k13:", value = 0.4, min = 0, max = 2, step=0.05),
                                             sliderInput("k31", "transition rate constant k31:", value = 0.2, min = 0, max = 2, step=0.05)
                            ),
                            conditionalPanel(condition = "input.distribution == '3' & input.parameterization=='2'",
                                             sliderInput("v3", "periph. comp. volume V3:", value = 5, min = 1, max = 10, step=0.5),
                                             sliderInput("q3", "inter comp. clearance Q3:", value = 0.5, min = 0.1, max = 2, step=0.1)
                            )
                        )
                ),
                tabItem(tabName = "elimination",
                        box(width = 4, status = "primary", solidHeader = TRUE, title="elimination",          
                            radioButtons("elimination","",c("linear" = "1","Michaelis Menten" = "2"))
                        ),
                        box(width = 4, status = "primary",          
                            conditionalPanel(condition = "input.elimination=='1' & input.parameterization=='1'",
                                             sliderInput("k", "elimination rate constant k:", value = 0.2, min = 0, max = 2, step=0.05)
                                             #uiOutput('myk')
                            ),
                            conditionalPanel(condition = "input.elimination == '1' & input.parameterization=='2'",
                                             sliderInput("Cl", "Clearance Cl:", value = 1, min = 0, max = 2, step=0.05)
                                             #uiOutput('myk')
                            ),
                            conditionalPanel(condition = "input.elimination == '2'",
                                             sliderInput("vm", "Vm:", value = 1, min = 0, max = 5, step=0.05),
                                             sliderInput("km", "Km:", value = 1, min = 0, max = 5, step=0.05)
                            )
                        )
                ),
                tabItem(tabName = "output",
                        box(width = 4, status = "primary", solidHeader = TRUE, title="output",          
                            radioButtons("log", "", c("linear scale" = FALSE,"log scale" = TRUE))  
                        ),
                        box(width = 4, status = "primary",
                            sliderInput("range", "time range", min = -10, max = 200, value = c(-5,100), step=5),
                            sliderInput("ngp", "grid size", min = 10, max = 1000, value = 200, step=10),
                            sliderInput("lsize", "line width", min = 0, max = 2, value = 0.75, step=0.25)
                        )
                )
              )
            )
    ),
    tabItem(tabName = "pkmodel",
            box( width = NULL, status = "primary", solidHeader = TRUE, title="pkmodelCode.R",                
                 downloadButton('downloadData1', 'Download'),
                 br(),br(),
                 verbatimTextOutput("RFile")
            )
    ),
    tabItem(tabName = "ui",
            box( width = NULL, status = "primary", solidHeader = TRUE, title="ui.R",
                 downloadButton('downloadData2', 'Download'),
                 br(),br(),
                 pre(includeText("ui.R"))
            )
    ),
    tabItem(tabName = "server",
            box( width = NULL, status = "primary", solidHeader = TRUE, title="server.R",
                 downloadButton('downloadData3', 'Download'),
                 br(),br(),
                 pre(includeText("server.R"))
            )
    ),
    tabItem(tabName = "readme",
            withMathJax(), 
            includeMarkdown("ReadMe.Rmd")
    ),
    tabItem(tabName = "perso",
            h4("Marc Lavielle"), 
            "Inria Saclay, Popix team",br(),
            "April 24th, 2015"
    ) 
  )
)

dashboardPage(
  dashboardHeader(title = "Build your PK model"),
  sidebar,
  body
)

server.R

Download

library(shinydashboard)
library(mlxR)
# my.dir=getwd()
# source("C:/Users/Marc/mlxtoolbox/mlxR220/initMlxR220.R")
# setwd(my.dir)

shinyServer(function(input, output){
  l.text <- ("library(mlxR)\n")
  pkmodel.text <- ("
res <- pkmodel(t,adm,p)

print(ggplot(data=res, aes(x=time, y=cc)) + geom_line(size=1) +
  xlab('time (h)') + ylab('concentration (mg/L)'))
")
  
  
  #   output$myk <- renderUI({
  #     if (input$parameterization=="1"){
  #       sliderInput("k", "elimination rate constant k:", value = 0.2, min = 0, max = 2, step=0.05)        
  #     }else{
  #       sliderInput("Cl", "Clearance Cl:", value = 1.5, min = 0, max = 3, step=0.05)        
  #     }
  #   })
  
  r <- reactive({ 
    param.name=NULL
    param.value=NULL
    if (input$administration=="oral"){
      if (input$absorption=="1"){
        if (input$delay1=="1"){
          param.name=c(param.name,'Tlag')
          param.value=c(param.value,input$tlag)
        }
        if (input$delay1=="2"){
          param.name=c(param.name,'Mtt','Ktr')
          param.value=c(param.value,input$mtt,(input$ntr+1)/input$mtt)
        }
        param.name=c(param.name,'ka')
        param.value=c(param.value,input$ka)
      }else{
        param.name=c(param.name,'Tk0')
        param.value=c(param.value,input$tk0)        
        if (input$delay0=="1"){
          param.name=c(param.name,'Tlag')
          param.value=c(param.value,input$tlag)
        }
      }
    }
    param.name=c(param.name,'V')
    param.value=c(param.value,input$v)

    if (input$distribution!="1"){
      if (input$parameterization=="1"){
        param.name=c(param.name,'k12','k21')
        param.value=c(param.value,input$k12,input$k21)
      }else{
        param.name=c(param.name,'Q2','V2')
        param.value=c(param.value,input$q2,input$v2)
      }
    }
    if (input$distribution=="3"){
      if (input$parameterization=="1"){
        param.name=c(param.name,'k13','k31')
        param.value=c(param.value,input$k13,input$k31)
      }else{
        param.name=c(param.name,'Q3','V3')
        param.value=c(param.value,input$q3,input$v3)
      }
    }
    
    if (input$elimination=="1"){
      if (input$parameterization=="1"){
        param.name=c(param.name,'k')
        param.value=c(param.value,input$k)
      }else{
        param.name=c(param.name,'Cl')
        param.value=c(param.value,input$Cl)
      }
    }else{
      param.name=c(param.name,'Vm','Km')
      param.value=c(param.value,input$vm,input$km)      
    }
    p=list(name=param.name,value=param.value)
    t.value=seq(input$range[1],input$range[2],length.out=input$ngp)
    t1=input$tfd
    t2=input$ii*(input$nd-1)+t1
    if (t2>=t1){
      t.dose=seq(t1,t2,by=input$ii)
      adm <- list(time=t.dose, amount=input$amt)
    }else{
      adm <- list(time=t1, amount=0)
    }
    if (input$administration == 'infusion'){
      adm$rate <- input$amt/input$tinf
    }
    
    #----------------------------------------------------    
    res   <- pkmodel(time=t.value,treatment=adm,parameter=p)
    #----------------------------------------------------    
    
    a.text <- 'adm <- list('
    if (input$nd==0){
      a.text <- 'adm <- list(amount=0, time=0'
    } else if(input$nd==1){
      a.text <- paste0('adm <- list(amount=',input$amt,', time=',input$tfd)
    } else{
      t1=input$tfd
      t2=input$ii*(input$nd-1)+t1
      a.text <- paste0('adm <- list(amount=',input$amt,
                       ', time=seq(',t1,', ',t2,', by=',input$ii,')')
    }
    if (input$administration == 'infusion'){
      a.text <- paste0(a.text,', tinf=',input$tinf,')')
    }else{
      a.text <- paste0(a.text,')')      
    }
    K <- length(param.value)
    p.text <- 'p <- c('
    for (k in (1:(K-1))){
      p.text <- paste0(p.text,param.name[k],'=',param.value[k],', ')
    }
    p.text <- paste0(p.text,param.name[K],'=',param.value[K],')')
    t.text <- paste0('t <- seq(',input$range[1],', ',input$range[2],
                     ', length=',length.out=input$ngp,')')
    pkmodel.text <- paste0(l.text,a.text,'\n',p.text,'\n',t.text,'\n',pkmodel.text)
    
    out <- list(res,pkmodel.text)
    return(out)
  })
  
  output$plot <- renderPlot({
    pl=ggplotmlx(data=r()[[1]], aes(x=time, y=cc)) + geom_line(size=input$lsize) 
    if (input$log==TRUE)
      pl=pl + scale_y_log10()
    print(pl)
  })
  
  output$table <- renderTable({ r()[[1]] })
  output$RFile  <- renderText(r()[[2]])
  
  output$downloadData1 <- downloadHandler(
    filename = "pkmodelCode.R",
    content = function(file) {
      write(r()[[2]], file)
    }
  )
  
  output$downloadData2 <- downloadHandler(
    filename = "ui.R",
    content = function(file) {
      con        = file("ui.R", open = "r")
      lines      = readLines(con)
      close(con)
      write(lines, file)
    }
  )
  
  output$downloadData3 <- downloadHandler(
    filename = "server.R",
    content = function(file) {
      con        = file("server.R", open = "r")
      lines      = readLines(con)
      close(con)
      write(lines, file)
    }
  )
  
  #   output$dynamic_value <- renderPrint({
  #     str(input$k)
  #   })
  
})

Build your PK model

Marc Lavielle

April 24th, 2015


Select the administration route and the dosage regimen in the menu Administration:

  • iv bolus, iv infusion, or oral administration,
  • time of first dose,
  • number of doses,
  • interdose interval,
  • amount.
  • infusion time (only for infusion),

Define the PK model in the menu Model:

  • Define the absorption process for oral administration in the submenu Absorption:
    • zero-order or first-order absorption process with the associated paramameters (duration of absorption \(Tk0\) or absorption rate constant \(ka\))
    • lag time \(Tlag\)
    • transit compartment model with mean transit time \(Mtt\) and number of transit compartments \(Ntr\) (only for first-order absorption)


  • Define the distribution process in the submenu Distribution:
    • 1, 2 or 3 compartments
    • volume of distribution \(V\)
    • transfer rates constants \(k_{12}\) and \(k_{21}\) between compartments 1 and 2
    • transfer rates constants \(k_{13}\) and \(k_{31}\) between compartments 1 and 3


  • Define the elimination process in the submenu Elimination:
    • linear or Michaelis Menten elimination process
    • elimination rate constant \(k\), or
    • Michaelis Menten parameters \(V_m\) and \(K_m\).


Define the output in the menu Output:

  • select the time range where the predicted concentration is computed,
  • select the number of time points of the grid where the predicted concentration is computed,
  • choose between linear and semi-log scale.


Marc Lavielle

Inria Saclay, Popix team
April 24th, 2015