server.R 29.9 KB
Newer Older
dmattek's avatar
dmattek committed
1
#
dmattek's avatar
dmattek committed
2 3 4 5
# Time Course Inspector: Shiny app for plotting time series data
# Author: Maciej Dobrzynski
#
# This is the server logic for a Shiny web application.
dmattek's avatar
dmattek committed
6 7 8 9
#

library(shiny)
library(shinyjs) #http://deanattali.com/shinyjs/
dmattek's avatar
dmattek committed
10 11
library(shinyBS) # for tooltips
library(shinycssloaders) # for loader animations
dmattek's avatar
dmattek committed
12 13
library(data.table)
library(ggplot2)
dmattek's avatar
dmattek committed
14
library(gplots) # for heatmap.2
15 16 17
library(plotly) # interactive plot
library(DT) # interactive tables

dmattek's avatar
dmattek committed
18
library(dendextend) # for color_branches
19
library(colorspace) # for palettes (used to colour dendrogram)
dmattek's avatar
dmattek committed
20 21
library(RColorBrewer)
library(scales) # for percentages on y scale
dmattek's avatar
dmattek committed
22
library(ggthemes) # nice colour palettes
23 24

library(sparcl) # sparse hierarchical and k-means
dmattek's avatar
dmattek committed
25
library(dtw) # for dynamic time warping
dmattek's avatar
dmattek committed
26
library(factoextra) # extract and visualize the output of multivariate data analyses 
dmattek's avatar
dmattek committed
27
library(imputeTS) # for interpolating NAs
28 29
library(robust) # for robust linear regression
library(MASS)
dmattek's avatar
dmattek committed
30 31
library(pracma) # for trapz used in AUC calculation

32

dmattek's avatar
dmattek committed
33

34
# Global parameters ----
dmattek's avatar
dmattek committed
35
# change to increase the limit of the upload file size
dmattek's avatar
dmattek committed
36
options(shiny.maxRequestSize = 200 * 1024 ^ 2)
dmattek's avatar
dmattek committed
37

dmattek's avatar
dmattek committed
38 39 40 41 42
# Important when joining, grouping or ordering numeric (i.e. double, POSIXct) columns.
# https://stackoverflow.com/questions/58230619/xy-join-of-keyed-data-table-fails-when-key-on-numeric-column-and-data-fread-fr
setNumericRounding(2)


dmattek's avatar
dmattek committed
43 44 45
# colour of loader spinner (shinycssloaders)
options(spinner.color="#00A8AA")

dmattek's avatar
dmattek committed
46
# Server logic ----
47
shinyServer(function(input, output, session) {
48
  useShinyjs()
dmattek's avatar
dmattek committed
49
  
50
  # This is only set at session start
dmattek's avatar
dmattek committed
51
  # We use this as a way to determine which input was
52 53
  # clicked in the dataInBoth reactive
  counter <- reactiveValues(
dmattek's avatar
dmattek committed
54 55 56
    # The value of actionButton is the number of times the button is pressed
    dataGen1        = isolate(input$inDataGen1),
    dataLoadNuc     = isolate(input$inButLoadNuc),
57 58
    dataLoadTrajRem = isolate(input$inButLoadTrajRem),
    dataLoadStim    = isolate(input$inButLoadStim)
dmattek's avatar
dmattek committed
59
  )
dmattek's avatar
dmattek committed
60 61 62 63 64 65 66 67 68

  nCellsCounter <- reactiveValues(
    nCellsOrig = 0,
    nCellsAfterOutlierTrim = 0
  )
    
  myReactVals = reactiveValues(
    outlierIDs = NULL
  )
dmattek's avatar
dmattek committed
69
  
dmattek's avatar
dmattek committed
70
  # UI-side-panel-data-load ----
dmattek's avatar
dmattek committed
71
  
dmattek's avatar
dmattek committed
72
  # Generate random dataset
73
  dataGen1 <- eventReactive(input$inDataGen1, {
74
    if (DEB)
75
      cat("server:dataGen1\n")
76
    
dmattek's avatar
dmattek committed
77
    return(LOCgenTraj2(n_perGroup = 20, sd_noise = 0.01, sampleFreq = 0.4, endTime = 40))
78 79
  })
  
dmattek's avatar
dmattek committed
80
  # Load main data file
81
  dataLoadNuc <- eventReactive(input$inButLoadNuc, {
82
    if (DEB)
83
      cat("server:dataLoadNuc\n")
84

85 86 87 88 89 90 91
    locFilePath = input$inFileLoadNuc$datapath
    
    counter$dataLoadNuc <- input$inButLoadNuc - 1
    
    if (is.null(locFilePath) || locFilePath == '')
      return(NULL)
    else {
92
      return(fread(locFilePath, strip.white = T))
93 94 95
    }
  })
  
dmattek's avatar
dmattek committed
96 97 98 99
  # This button will reset the inFileLoad
  observeEvent(input$butReset, {
    reset("inFileLoadNuc")  # reset is a shinyjs function
  })
100

dmattek's avatar
dmattek committed
101
  # Load data with trajectories to remove
102
  dataLoadTrajRem <- eventReactive(input$inButLoadTrajRem, {
103
    if (DEB)
104
      cat(file = stdout(), "server:dataLoadTrajRem\n")
105
    
106 107 108 109 110 111 112 113 114 115
    locFilePath = input$inFileLoadTrajRem$datapath
    
    counter$dataLoadTrajRem <- input$inButLoadTrajRem - 1
    
    if (is.null(locFilePath) || locFilePath == '')
      return(NULL)
    else {
      return(fread(locFilePath))
    }
  })
dmattek's avatar
dmattek committed
116
  
dmattek's avatar
dmattek committed
117
  # Load data with stimulation pattern
118
  dataLoadStim <- eventReactive(input$inButLoadStim, {
119
    if (DEB)
120
      cat(file = stdout(), "server:dataLoadStim\n")
121
    
122 123 124 125 126 127 128 129 130 131 132 133
    locFilePath = input$inFileLoadStim$datapath
    
    counter$dataLoadStim <- input$inButLoadStim - 1
    
    if (is.null(locFilePath) || locFilePath == '')
      return(NULL)
    else {
      return(fread(locFilePath))
    }
  })
  
    
dmattek's avatar
dmattek committed
134 135
  # UI for loading csv with cell IDs for trajectory removal
  output$uiFileLoadTrajRem = renderUI({
136
    if (DEB)
137
      cat(file = stdout(), 'server:uiFileLoadTrajRem\n')
dmattek's avatar
dmattek committed
138 139 140 141
    
    if(input$chBtrajRem) 
      fileInput(
        'inFileLoadTrajRem',
dmattek's avatar
dmattek committed
142
        "Select file and click Load Data",
dmattek's avatar
dmattek committed
143 144 145 146 147
        accept = c('text/csv', 'text/comma-separated-values,text/plain')
      )
  })
  
  output$uiButLoadTrajRem = renderUI({
148
    if (DEB)
149
      cat(file = stdout(), 'server:uiButLoadTrajRem\n')
dmattek's avatar
dmattek committed
150 151 152 153 154
    
    if(input$chBtrajRem)
      actionButton("inButLoadTrajRem", "Load Data")
  })

155 156
  # UI for loading csv with stimulation pattern
  output$uiFileLoadStim = renderUI({
157
    if (DEB)
158
      cat(file = stdout(), 'server:uiFileLoadStim\n')
dmattek's avatar
dmattek committed
159
    
160 161 162
    if(input$chBstim) 
      fileInput(
        'inFileLoadStim',
dmattek's avatar
dmattek committed
163
        "Select file and click Load Data",
164 165 166 167 168
        accept = c('text/csv', 'text/comma-separated-values,text/plain')
      )
  })
  
  output$uiButLoadStim = renderUI({
169
    if (DEB)
170
      cat(file = stdout(), 'server:uiButLoadStim\n')
dmattek's avatar
dmattek committed
171
    
172 173
    if(input$chBstim)
      actionButton("inButLoadStim", "Load Data")
dmattek's avatar
dmattek committed
174 175
  })
  
176

dmattek's avatar
dmattek committed
177
  
dmattek's avatar
dmattek committed
178
  # UI-side-panel-column-selection ----
dmattek's avatar
dmattek committed
179
  output$varSelTrackLabel = renderUI({
180
    if (DEB)
181
      cat(file = stdout(), 'server:varSelTrackLabel\n')
182
    
dmattek's avatar
dmattek committed
183
    locCols = getDataNucCols()
184
    locColSel = locCols[grep('(T|t)rack|ID|id', locCols)[1]] # index 1 at the end in case more matches; select 1st; matches TrackLabel, tracklabel, Track Label etc
dmattek's avatar
dmattek committed
185 186 187
    
    selectInput(
      'inSelTrackLabel',
dmattek's avatar
dmattek committed
188
      'Track ID column',
dmattek's avatar
dmattek committed
189 190 191 192 193 194 195
      locCols,
      width = '100%',
      selected = locColSel
    )
  })
  
  output$varSelTime = renderUI({
196
    if (DEB)
197
      cat(file = stdout(), 'server:varSelTime\n')
198
    
dmattek's avatar
dmattek committed
199
    locCols = getDataNucCols()
200
    locColSel = locCols[grep('(T|t)ime|Metadata_T', locCols)[1]] # index 1 at the end in case more matches; select 1st; matches RealTime, realtime, real time, etc.
dmattek's avatar
dmattek committed
201 202 203
    
    selectInput(
      'inSelTime',
dmattek's avatar
dmattek committed
204
      'Time column',
dmattek's avatar
dmattek committed
205 206 207 208 209
      locCols,
      width = '100%',
      selected = locColSel
    )
  })
210 211

  output$varSelTimeFreq = renderUI({
212
    if (DEB)
213
      cat(file = stdout(), 'server:varSelTimeFreq\n')
214
    
215 216 217
    if (input$chBtrajInter) {
      numericInput(
        'inSelTimeFreq',
dmattek's avatar
dmattek committed
218
        'Interval between 2 time points',
219 220 221 222 223 224
        min = 1,
        step = 1,
        width = '100%',
        value = 1
      )
    }
225
  })
dmattek's avatar
dmattek committed
226
  
dmattek's avatar
dmattek committed
227
  # This is the main field to select plot facet grouping
dmattek's avatar
dmattek committed
228
  # It's typically a column with the entire experimental description,
dmattek's avatar
dmattek committed
229 230
  # e.g.1 Stim_All_Ch or Stim_All_S.
  # e.g.2 a combination of 3 columns called Stimulation_...
dmattek's avatar
dmattek committed
231
  output$varSelGroup = renderUI({
232
    if (DEB)
233
      cat(file = stdout(), 'server:varSelGroup\n')
dmattek's avatar
dmattek committed
234
    
dmattek's avatar
dmattek committed
235 236 237 238 239
    if (input$chBgroup) {
      
      locCols = getDataNucCols()
      
      if (!is.null(locCols)) {
240
        locColSel = locCols[grep('(G|g)roup|(S|s)tim|(S|s)timulation|(S|s)ite', locCols)[1]]
241 242

        #cat('UI varSelGroup::locColSel ', locColSel, '\n')
dmattek's avatar
dmattek committed
243 244
        selectInput(
          'inSelGroup',
dmattek's avatar
dmattek committed
245
          'Grouping columns',
dmattek's avatar
dmattek committed
246 247 248 249 250
          locCols,
          width = '100%',
          selected = locColSel,
          multiple = TRUE
        )
dmattek's avatar
dmattek committed
251 252 253 254
      }
    }
  })
  
255 256
  # UI for selecting grouping to add to track ID to make 
  # the track ID unique across entire dataset
dmattek's avatar
dmattek committed
257
  output$varSelSite = renderUI({
258
    if (DEB)
259
      cat(file = stdout(), 'server:varSelSite\n')
dmattek's avatar
dmattek committed
260
    
261
    if (input$chBtrackUni) {
dmattek's avatar
dmattek committed
262
      locCols = getDataNucCols()
Maciej Dobrzynski's avatar
Maciej Dobrzynski committed
263
      locColSel = locCols[grep('(S|s)ite|(S|s)eries|(F|f)ov|(G|g)roup', locCols)[1]] # index 1 at the end in case more matches; select 1st
dmattek's avatar
dmattek committed
264 265 266
      
      selectInput(
        'inSelSite',
dmattek's avatar
dmattek committed
267
        'Prepend track ID with',
dmattek's avatar
dmattek committed
268 269
        locCols,
        width = '100%',
270 271
        selected = locColSel,
        multiple = T
dmattek's avatar
dmattek committed
272 273
      )
    }
dmattek's avatar
dmattek committed
274 275 276 277
  })
  
  
  output$varSelMeas1 = renderUI({
278
    if (DEB)
279
      cat(file = stdout(), 'server:varSelMeas1\n')
dmattek's avatar
dmattek committed
280 281 282
    locCols = getDataNucCols()
    
    if (!is.null(locCols)) {
Maciej Dobrzynski's avatar
Maciej Dobrzynski committed
283
      locColSel = locCols[grep('(R|r)atio|(I|i)ntensity|(Y|y)|(M|m)eas', locCols)[1]] # index 1 at the end in case more matches; select 1st
dmattek's avatar
dmattek committed
284

dmattek's avatar
dmattek committed
285 286
      selectInput(
        'inSelMeas1',
dmattek's avatar
dmattek committed
287
        '1st measurement column',
dmattek's avatar
dmattek committed
288 289 290 291 292 293 294 295 296
        locCols,
        width = '100%',
        selected = locColSel
      )
    }
  })
  
  
  output$varSelMeas2 = renderUI({
297
    if (DEB)
298
      cat(file = stdout(), 'server:varSelMeas2\n')
299
    
dmattek's avatar
dmattek committed
300 301 302 303
    locCols = getDataNucCols()
    
    if (!is.null(locCols) &&
        !(input$inSelMath %in% c('', '1 / '))) {
Maciej Dobrzynski's avatar
Maciej Dobrzynski committed
304
      locColSel = locCols[grep('(R|r)atio|(I|i)ntensity|(Y|y)|(M|m)eas', locCols)[1]] # index 1 at the end in case more matches; select 1st
dmattek's avatar
dmattek committed
305

dmattek's avatar
dmattek committed
306 307
      selectInput(
        'inSelMeas2',
dmattek's avatar
dmattek committed
308
        '2nd measurement column',
dmattek's avatar
dmattek committed
309 310 311 312 313 314 315
        locCols,
        width = '100%',
        selected = locColSel
      )
    }
  })
  
dmattek's avatar
dmattek committed
316
  # UI-side-panel-trim x-axis (time) ----
dmattek's avatar
dmattek committed
317
  output$uiSlTimeTrim = renderUI({
318
    if (DEB)
319
      cat(file = stdout(), 'server:uiSlTimeTrim\n')
dmattek's avatar
dmattek committed
320 321 322 323 324 325 326 327 328 329 330 331
    
    if (input$chBtimeTrim) {
      locTpts  = getDataTpts()
      
      if(is.null(locTpts))
        return(NULL)
      
      locRTmin = min(locTpts)
      locRTmax = max(locTpts)
      
      sliderInput(
        'slTimeTrim',
dmattek's avatar
dmattek committed
332
        label = 'Use time range',
dmattek's avatar
dmattek committed
333 334 335 336 337 338 339 340
        min = locRTmin,
        max = locRTmax,
        value = c(locRTmin, locRTmax),
        step = 1
      )
      
    }
  })
dmattek's avatar
dmattek committed
341
  
dmattek's avatar
dmattek committed
342
  # UI-side-panel-normalization ----
343 344 345 346
  
  # select normalisation method
  # - fold-change calculates fold change with respect to the mean
  # - z-score calculates z-score of the selected regione of the time series
dmattek's avatar
dmattek committed
347
  output$uiChBnorm = renderUI({
348
    if (DEB)
349
      cat(file = stdout(), 'server:uiChBnorm\n')
dmattek's avatar
dmattek committed
350 351
    
    if (input$chBnorm) {
352
      tagList(
dmattek's avatar
dmattek committed
353 354 355
      radioButtons(
        'rBnormMeth',
        label = 'Select method',
356 357 358
        choices = list('fold-change' = 'mean', 'z-score' = 'z.score'),
        width = "40%"
      ),
dmattek's avatar
dmattek committed
359
      bsTooltip('rBnormMeth', helpText.server[["rBnormMeth"]], placement = "top", trigger = "hover", options = NULL)
dmattek's avatar
dmattek committed
360 361 362 363
      )
    }
  })
  
364
  # select the region of the time series for normalisation
dmattek's avatar
dmattek committed
365
  output$uiSlNorm = renderUI({
366
    if (DEB)
367
      cat(file = stdout(), 'server:uiSlNorm\n')
dmattek's avatar
dmattek committed
368 369 370 371 372 373 374 375 376 377
    
    if (input$chBnorm) {
      locTpts  = getDataTpts()
      
      if(is.null(locTpts))
        return(NULL)
      
      locRTmin = min(locTpts)
      locRTmax = max(locTpts)
      
378
      tagList(
dmattek's avatar
dmattek committed
379 380
      sliderInput(
        'slNormRtMinMax',
381
        label = 'Time span',
dmattek's avatar
dmattek committed
382 383
        min = locRTmin,
        max = locRTmax,
dmattek's avatar
dmattek committed
384 385
        value = c(locRTmin, 0.1 * locRTmax), 
        step = 1
386
      ),
dmattek's avatar
dmattek committed
387
      bsTooltip('slNormRtMinMax', helpText.server[["slNormRtMinMax"]], placement = "top", trigger = "hover", options = NULL)
dmattek's avatar
dmattek committed
388 389 390 391
      )
    }
  })
  
392
  # use robust stats (median instead of mean, mad instead of sd)
dmattek's avatar
dmattek committed
393
  output$uiChBnormRobust = renderUI({
394
    if (DEB)
395
      cat(file = stdout(), 'server:uiChBnormRobust\n')
dmattek's avatar
dmattek committed
396 397
    
    if (input$chBnorm) {
398
      tagList(
dmattek's avatar
dmattek committed
399 400
      checkboxInput('chBnormRobust',
                    label = 'Robust stats',
401 402
                    FALSE, 
                    width = "40%"),
dmattek's avatar
dmattek committed
403
      bsTooltip('chBnormRobust', helpText.server[["chBnormRobust"]], placement = "top", trigger = "hover", options = NULL)
404
      )
dmattek's avatar
dmattek committed
405 406 407
    }
  })
  
408
  # choose whether normalisation should be calculated for the entire dataset, group, or trajectory
dmattek's avatar
dmattek committed
409
  output$uiChBnormGroup = renderUI({
410
    if (DEB)
411
      cat(file = stdout(), 'server:uiChBnormGroup\n')
dmattek's avatar
dmattek committed
412 413
    
    if (input$chBnorm) {
414
      tagList(
dmattek's avatar
dmattek committed
415
      radioButtons('chBnormGroup',
dmattek's avatar
Mod:  
dmattek committed
416
                   label = 'Normalisation grouping',
417 418
                   choices = list('Entire dataset' = 'none', 'Per group' = 'group', 'Per trajectory' = 'id'), 
                   width = "40%"),
dmattek's avatar
dmattek committed
419
      bsTooltip('chBnormGroup', helpText.server[["chBnormGroup"]], placement = "top", trigger = "hover", options = NULL)
420
      )
dmattek's avatar
dmattek committed
421 422 423 424
    }
  })
  
  
425 426 427 428 429
  # Pop-overs ----
  addPopover(session, 
             "alDataFormat",
             title = "Data format",
             content = helpText.server[["alDataFormat"]],
dmattek's avatar
dmattek committed
430
             trigger = "click")
dmattek's avatar
dmattek committed
431
  
dmattek's avatar
dmattek committed
432

dmattek's avatar
dmattek committed
433
  # Processing-data ----
dmattek's avatar
dmattek committed
434
  
435 436 437 438 439 440 441 442 443 444 445 446
  dataInBoth <- reactive({
    # Without direct references to inDataGen1,2 and inFileLoad, inDataGen2
    #    does not trigger running this reactive once inDataGen1 is used.
    # This is one of the more nuanced areas of reactive programming in shiny
    #    due to the if else logic, it isn't fetched once inDataGen1 is available
    # The morale is use direct retrieval of inputs to guarantee they are available
    #    for if else logic checks!
    
    locInGen1 = input$inDataGen1
    locInLoadNuc = input$inButLoadNuc
    #locInLoadStim = input$inButLoadStim
    
447
    # Don't wrap around if(DEB) !!!
448
    cat(
449
      "server:dataInBoth\n   inGen1: ",
450
      locInGen1,
451
      "      prev=",
452
      isolate(counter$dataGen1),
453
      "\n   inDataNuc: ",
454 455 456 457 458 459 460 461 462 463
      locInLoadNuc,
      "   prev=",
      isolate(counter$dataLoadNuc),
      # "\ninDataStim: ",
      # locInLoadStim,
      # "   prev=",
      # isolate(counter$dataLoadStim),
      "\n"
    )
    
464
    # isolate the checks of the counter reactiveValues
465 466
    # as we set the values in this same reactive
    if (locInGen1 != isolate(counter$dataGen1)) {
467
      cat("server:dataInBoth if inDataGen1\n")
468 469 470 471
      dm = dataGen1()
      # no need to isolate updating the counter reactive values!
      counter$dataGen1 <- locInGen1
    } else if (locInLoadNuc != isolate(counter$dataLoadNuc)) {
472
      cat("server:dataInBoth if inDataLoadNuc\n")
473
      dm = dataLoadNuc()
474 475 476
      
      # convert to long format if radio box set to "wide"
      # the input data in long format should contain:
477
      # - the first row with a header: group, track id, time points as columns with numeric header
478 479
      # - consecutive rows with time series, where columns are time points
      if (input$inRbutLongWide == 1) {
480 481 482 483 484 485 486 487 488 489 490 491
        print(length(names(dm)))
        
        # data in wide format requires at least 3 columns: grouping, track id, 1 time point
        if (length(names(dm)) < 3) {
          dm = NULL
          
          createAlert(session, "alertAnchorSidePanelDataFormat", "alertWideTooFewColumns", 
                      title = "Error",
                      content = helpText.server[["alertWideTooFewColumns"]], 
                      append = FALSE,
                      style = "danger")
          
492

493 494 495 496 497 498 499 500 501 502 503 504 505 506 507 508 509 510 511 512 513 514 515 516 517 518 519 520 521 522 523 524 525
        } else {
          closeAlert(session, "alertWideTooFewColumns")

          # obtain column headers from the wide format data
          # headers for grouping and track id columns
          loc.cols.idvars = names(dm)[1:2]
          
          # headers for time columns
          loc.cols.time = names(dm)[c(-1, -2)]
          
          # check if time columns are numeric
          # from https://stackoverflow.com/a/21154566/1898713
          loc.cols.time.numres = grepl("[-]?[0-9]+[.]?[0-9]*|[-]?[0-9]+[L]?|[-]?[0-9]+[.]?[0-9]*[eE][0-9]+", loc.cols.time)
          
          # melt the table only if time columns are numeric
          if (sum(!loc.cols.time.numres) == 0) {
            closeAlert(session, "alertWideMissesNumericTime")
            
            # long to wide
            dm = melt(dm, id.vars = loc.cols.idvars, variable.name = COLRT, value.name = COLY)
            
            # convert column names with time points to a number
            dm[, (COLRT) := as.numeric(levels(get(COLRT)))[get(COLRT)]]
            
          } else {
            dm = NULL

            createAlert(session, "alertAnchorSidePanelDataFormat", "alertWideMissesNumericTime", title = "Error",
                        content = helpText.server[["alertWideMissesNumericTime"]], 
                        append = FALSE,
                        style = "danger")
          }
        }
526 527
      }
      
528 529 530
      # no need to isolate updating the counter reactive values!
      counter$dataLoadNuc <- locInLoadNuc
    } else {
531
      cat("server:dataInBoth else\n")
532 533
      dm = NULL
    }
534
    
535 536 537 538
    return(dm)
  })
  
  # return column names of the main dt
dmattek's avatar
dmattek committed
539
  getDataNucCols <- reactive({
540
    if (DEB)
541
      cat(file = stdout(), 'server:getDataNucCols: in\n')
542
    
543 544 545 546 547 548 549 550 551 552
    loc.dt = dataInBoth()
    
    if (is.null(loc.dt))
      return(NULL)
    else
      return(colnames(loc.dt))
  })
  
  # return dt with an added column with unique track object label
  dataMod <- reactive({
553
    if (DEB)
554
      cat(file = stdout(), 'server:dataMod\n')
555
    
556 557
    loc.dt = dataInBoth()
    
dmattek's avatar
dmattek committed
558
    if (is.null(loc.dt))
559 560
      return(NULL)
    
561
    if (input$chBtrackUni) {
562
      # create unique track ID based on columns specified in input$inSelSite field and combine with input$inSelTrackLabel
563
      loc.dt[, (COLIDUNI) := do.call(paste, c(.SD, sep = "_")), .SDcols = c(input$inSelSite, input$inSelTrackLabel) ]
dmattek's avatar
dmattek committed
564
    } else {
565
      # stay with track ID provided in the loaded dataset; has to be unique
566
      loc.dt[, (COLIDUNI) := get(input$inSelTrackLabel)]
dmattek's avatar
dmattek committed
567 568
    }
    
dmattek's avatar
dmattek committed
569
    
dmattek's avatar
dmattek committed
570 571
    # remove trajectories based on uploaded csv
    if (input$chBtrajRem) {
572
      if (DEB)
573
        cat(file = stdout(), 'server:dataMod: trajRem not NULL\n')
dmattek's avatar
dmattek committed
574 575
      
      loc.dt.rem = dataLoadTrajRem()
dmattek's avatar
dmattek committed
576
      loc.dt = loc.dt[!(trackObjectsLabelUni %in% loc.dt.rem[[1]])]
dmattek's avatar
dmattek committed
577 578
    }
    
dmattek's avatar
dmattek committed
579 580 581
    # check if NAs present
    
    
582 583 584
    return(loc.dt)
  })
  
dmattek's avatar
dmattek committed
585 586 587
  # return all unique track object labels (created in dataMod)
  # This will be used to display in UI for trajectory highlighting
  getDataTrackObjLabUni <- reactive({
588
    if (DEB)
589
      cat(file = stdout(), 'server:getDataTrackObjLabUni\n')
590
    
dmattek's avatar
dmattek committed
591
    loc.dt = dataMod()
592
    
dmattek's avatar
dmattek committed
593 594 595 596
    if (is.null(loc.dt))
      return(NULL)
    else
      return(unique(loc.dt$trackObjectsLabelUni))
597 598
  })
  
dmattek's avatar
Mod:  
dmattek committed
599
  
dmattek's avatar
dmattek committed
600 601 602
  # return all unique time points (real time)
  # This will be used to display in UI for box-plot
  # These timepoints are from the original dt and aren't affected by trimming of x-axis
dmattek's avatar
dmattek committed
603
  getDataTpts <- reactive({
604
    if (DEB)
605
      cat(file = stdout(), 'server:getDataTpts\n')
606
    
dmattek's avatar
dmattek committed
607
    loc.dt = dataMod()
608
    
dmattek's avatar
dmattek committed
609 610 611 612
    if (is.null(loc.dt))
      return(NULL)
    else
      return(unique(loc.dt[[input$inSelTime]]))
613 614
  })
  
dmattek's avatar
dmattek committed
615
  
616 617 618
  
  # prepare data for plotting time courses
  # returns dt with these columns:
dmattek's avatar
dmattek committed
619
  #    realtime - selected from input
dmattek's avatar
dmattek committed
620
  #    y        - measurement selected from input
dmattek's avatar
dmattek committed
621
  #               (can be a single column or result of an operation on two cols)
622 623
  #    id       - trackObjectsLabelUni; created in dataMod based on TrackObjects_Label
  #               and FOV column such as Series or Site (if TrackObjects_Label not unique across entire dataset)
dmattek's avatar
dmattek committed
624 625
  #    group    - grouping variable for facetting from input
  #    mid.in   - column with trajectory selection status from the input file or
626 627 628 629
  #               highlight status from UI 
  #               (column created if mid.in present in uploaded data or tracks are selected in the UI)
  #    obj.num  - created if ObjectNumber column present in the input data 
  #    pos.x,y  - created if columns with x and y positions present in the input data
630
  data4trajPlot <- reactive({
631
    if (DEB)
632
      cat(file = stdout(), 'server:data4trajPlot\n')
633 634
    
    loc.dt = dataMod()
dmattek's avatar
dmattek committed
635
    if (is.null(loc.dt))
636 637
      return(NULL)
    
638
    # create expression for 'y' column based on measurements and math operations selected in UI
dmattek's avatar
dmattek committed
639
    if (input$inSelMath == '')
640 641 642 643 644 645
      loc.s.y = input$inSelMeas1
    else if (input$inSelMath == '1 / ')
      loc.s.y = paste0(input$inSelMath, input$inSelMeas1)
    else
      loc.s.y = paste0(input$inSelMeas1, input$inSelMath, input$inSelMeas2)
    
646
    # create expression for 'group' column
647 648
    # creates a merged column based on other columns from input
    # used for grouping of plot facets
dmattek's avatar
dmattek committed
649 650 651 652 653 654 655 656 657 658 659
    if (input$chBgroup) {
      if(length(input$inSelGroup) == 0)
        return(NULL)
      
      loc.s.gr = sprintf("paste(%s, sep=';')",
                         paste(input$inSelGroup, sep = '', collapse = ','))
    } else {
      # if no grouping required, fill 'group' column with 0
      # because all the plotting relies on the presence of the group column
      loc.s.gr = "paste('0')"
    }
660
    
dmattek's avatar
dmattek committed
661 662

    # column name with time
663 664
    loc.s.rt = input$inSelTime
    
dmattek's avatar
dmattek committed
665 666
    # Assign tracks selected for highlighting in UI
    loc.tracks.highlight = input$inSelHighlight
667
    locButHighlight = input$chBhighlightTraj
dmattek's avatar
dmattek committed
668
    
dmattek's avatar
dmattek committed
669 670
    
    # Find column names with position
671
    loc.s.pos.x = names(loc.dt)[grep('(L|l)ocation.*X|(P|p)os.x|(P|p)osx', names(loc.dt))[1]]
672
    loc.s.pos.y = names(loc.dt)[grep('(L|l)ocation.*Y|(P|p)os.y|(P|p)osy', names(loc.dt))[1]]
dmattek's avatar
dmattek committed
673
    
674
    if (DEB)
675
      cat('server:data4trajPlot:\n   Position columns: ', loc.s.pos.x, loc.s.pos.y, '\n')
676 677
    
    if (!is.na(loc.s.pos.x) & !is.na(loc.s.pos.y))
dmattek's avatar
dmattek committed
678 679 680 681
      locPos = TRUE
    else
      locPos = FALSE
    
682 683 684 685
    
    # Find column names with ObjectNumber
    # This is different from TrackObject_Label and is handy to keep
    # because labels on segmented images are typically ObjectNumber
686 687 688 689 690 691
    loc.s.objnum = names(loc.dt)[grep('(O|o)bject(N|n)umber', names(loc.dt))[1]]
    #cat('data4trajPlot::loc.s.objnum ', loc.s.objnum, '\n')
    if (is.na(loc.s.objnum)) {
      locObjNum = FALSE
    }
    else {
dmattek's avatar
dmattek committed
692
      loc.s.objnum = loc.s.objnum[1]
693
      locObjNum = TRUE
dmattek's avatar
dmattek committed
694
    }
695 696
    
    
697 698
    # if dataset contains column mid.in with trajectory filtering status,
    # then, include it in plotting
dmattek's avatar
dmattek committed
699
    if (sum(names(loc.dt) %in% COLIN) > 0)
700 701 702 703 704 705
      locMidIn = TRUE
    else
      locMidIn = FALSE
    
    ## Build expression for selecting columns from loc.dt
    # Core columns
dmattek's avatar
dmattek committed
706 707 708 709
    s.colexpr = paste0('.(',  COLY, ' = ', loc.s.y,
                       ', ', COLID, ' = ', COLIDUNI, 
                       ', ', COLGR, ' = ', loc.s.gr,
                       ', ', COLRT, ' = ', loc.s.rt)
710 711
    
    # account for the presence of 'mid.in' column in uploaded data
dmattek's avatar
dmattek committed
712
    # future: choose this column in UI
713 714
    if(locMidIn)
      s.colexpr = paste0(s.colexpr, 
dmattek's avatar
dmattek committed
715
                         ',', COLIN, ' = ', COLIN)
716 717 718 719
    
    # include position x,y columns in uploaded data
    if(locPos)
      s.colexpr = paste0(s.colexpr, 
dmattek's avatar
dmattek committed
720 721
                         ', ', COLPOSX, '= ', loc.s.pos.x,
                         ', ', COLPOSY, '= ', loc.s.pos.y)
722 723 724 725
    
    # include ObjectNumber column
    if(locObjNum)
      s.colexpr = paste0(s.colexpr, 
dmattek's avatar
dmattek committed
726
                         ', ', COLOBJN, ' = ', loc.s.objnum)
727 728 729 730 731 732 733
    
    # close bracket, finish the expression
    s.colexpr = paste0(s.colexpr, ')')
    
    # create final dt for output based on columns selected above
    loc.out = loc.dt[, eval(parse(text = s.colexpr))]
    
734 735 736 737 738
    # Convert track ID to a factor.
    # This is necessary for, e.g. merging data with cluster assignments.
    # If input dataset has track ID as a number, such a merge would fail.
    loc.out[, (COLID) := as.factor(get(COLID))]
    
739 740 741 742 743 744
    
    # if track selection ON
    if (locButHighlight){
      # add a 3rd level with status of track selection
      # to a column with trajectory filtering status in the uploaded file
      if(locMidIn)
dmattek's avatar
dmattek committed
745
        loc.out[, mid.in := ifelse(get(COLID) %in% loc.tracks.highlight, 'SELECTED', get(COLIN))]
746
      else
dmattek's avatar
Mod:  
dmattek committed
747
        # add a column with status of track selection
dmattek's avatar
dmattek committed
748
        loc.out[, mid.in := ifelse(get(COLID) %in% loc.tracks.highlight, 'SELECTED', 'NOT SEL')]
749
    }
750
      
dmattek's avatar
dmattek committed
751

752
    ## Interpolate missing data and NA data points
753
    # From: https://stackoverflow.com/questions/28073752/r-how-to-add-rows-for-missing-values-for-unique-group-sequences
754
    # Tracks are interpolated only within first and last time points of every track id
755 756
    # Datasets can have different realtime frequency (e.g. every 1', 2', etc),
    # or the frame number metadata can be missing, as is the case for tCourseSelected files that already have realtime column.
757
    # Therefore, we cannot rely on that info to get time frequency; user must provide this number!
758
    
dmattek's avatar
dmattek committed
759 760 761 762 763 764 765 766 767
    # check if NA's present
    if (sum(is.na(loc.out[[COLY]])))
      createAlert(session, "alertAnchorSidePanelNAsPresent", "alertNAsPresent", title = "Warning",
                  content = helpText.server[["alertNAsPresent"]], 
                  append = FALSE,
                  style = "warning")
    else
      closeAlert(session, "alertNAsPresent")
    
768
    setkeyv(loc.out, c(COLGR, COLID, COLRT))
769

770
    if (input$chBtrajInter) {
dmattek's avatar
dmattek committed
771 772 773 774 775 776
      # here we fill missing rows with NA's
      loc.out = loc.out[setkeyv(loc.out[, 
                                        .(seq(min(get(COLRT), na.rm = T), 
                                              max(get(COLRT), na.rm = T), 
                                              input$inSelTimeFreq)), 
                                        by = c(COLGR, COLID)], c(COLGR, COLID, 'V1'))]
777 778
      
      # x-check: print all rows with NA's
779
      if (DEB) {
780
        cat(file = stdout(), 'server:data4trajPlot: Rows with NAs:\n')
781 782
        print(loc.out[rowSums(is.na(loc.out)) > 0, ])
      }
783 784 785 786
      
      # NA's may be already present in the dataset'.
      # Interpolate (linear) them with na.interpolate as well
      if(locPos)
dmattek's avatar
dmattek committed
787
        s.cols = c(COLY, COLPOSX, COLPOSY)
788
      else
dmattek's avatar
dmattek committed
789
        s.cols = c(COLY)
790
      
791 792 793 794 795 796 797 798
      # Interpolated columns should be of type numeric (float)
      # This is to ensure that interpolated columns are of porper type.
      
      # Apparently the loop is faster than lapply+SDcols
      for(col in s.cols) {
        #loc.out[, (col) := as.numeric(get(col))]
        data.table::set(loc.out, j = col, value = as.numeric(loc.out[[col]]))

dmattek's avatar
dmattek committed
799
        loc.out[, (col) := na_interpolation(get(col)), by = c(COLID)]        
800 801 802
      }
      
      # loc.out[, (s.cols) := lapply(.SD, na.interpolation), by = c(COLID), .SDcols = s.cols]
803 804 805 806 807 808 809 810 811 812 813 814 815 816
      
      
      # !!! Current issue with interpolation:
      # The column mid.in is not taken into account.
      # If a trajectory is selected in the UI,
      # the mid.in column is added (if it doesn't already exist in the dataset),
      # and for the interpolated point, it will still be NA. Not really an issue.
      #
      # Also, think about the current option of having mid.in column in the uploaded dataset.
      # Keep it? Expand it?
      # Create a UI filed for selecting the column with mid.in data.
      # What to do with that column during interpolation (see above)
      
    }    
dmattek's avatar
Mod:  
dmattek committed
817
    
818
    ## Trim x-axis (time)
dmattek's avatar
dmattek committed
819
    if(input$chBtimeTrim) {
dmattek's avatar
dmattek committed
820
      loc.out = loc.out[get(COLRT) >= input$slTimeTrim[[1]] & get(COLRT) <= input$slTimeTrim[[2]] ]
dmattek's avatar
dmattek committed
821
    }
dmattek's avatar
dmattek committed
822
    
823
    ## Normalization
dmattek's avatar
dmattek committed
824
    # F-n normTraj adds additional column with .norm suffix
dmattek's avatar
dmattek committed
825
    if (input$chBnorm) {
dmattek's avatar
dmattek committed
826
      loc.out = LOCnormTraj(
dmattek's avatar
dmattek committed
827
        in.dt = loc.out,
dmattek's avatar
dmattek committed
828 829
        in.meas.col = COLY,
        in.rt.col = COLRT,
dmattek's avatar
dmattek committed
830 831 832 833 834 835 836
        in.rt.min = input$slNormRtMinMax[1],
        in.rt.max = input$slNormRtMinMax[2],
        in.type = input$rBnormMeth,
        in.robust = input$chBnormRobust,
        in.by.cols = if(input$chBnormGroup %in% 'none') NULL else input$chBnormGroup
      )
      
dmattek's avatar
dmattek committed
837 838
      # Column with normalized data is renamed to the original name
      # Further code assumes column name y produced by data4trajPlot
Maciej Dobrzynski's avatar
Maciej Dobrzynski committed
839 840
      
      loc.out[, c(COLY) := NULL]
dmattek's avatar
dmattek committed
841
      setnames(loc.out, 'y.norm', COLY)
dmattek's avatar
dmattek committed
842 843 844
    }
    
    return(loc.out)
dmattek's avatar
dmattek committed
845 846
  })
  
dmattek's avatar
dmattek committed
847 848
  
  # prepare data for clustering
dmattek's avatar
dmattek committed
849
  # convert from long to wide; return a matrix with:
dmattek's avatar
dmattek committed
850 851 852
  # cells as columns
  # time points as rows
  data4clust <- reactive({
853
    if (DEB)  
854
      cat(file = stdout(), 'server:data4clust\n')
dmattek's avatar
dmattek committed
855
    
dmattek's avatar
dmattek committed
856
    loc.dt = data4trajPlotNoOut()
dmattek's avatar
dmattek committed
857 858 859
    if (is.null(loc.dt))
      return(NULL)
    
860 861 862 863
    # convert from long to wide format
    loc.dt.wide = dcast(loc.dt, 
                    reformulate(response = COLID, termlabels = COLRT), 
                    value.var = COLY)
dmattek's avatar
dmattek committed
864
    
865 866
    # store row names for later
    loc.rownames = loc.dt.wide[[COLID]]
dmattek's avatar
Mod:  
dmattek committed
867
    
868 869
    # omit first column that contains row names
    loc.m.out = as.matrix(loc.dt.wide[, -1])
dmattek's avatar
dmattek committed
870
    
871 872
    # assign row names to the matrix
    rownames(loc.m.out) = loc.rownames
dmattek's avatar
dmattek committed
873
    
874
    return(loc.m.out)
dmattek's avatar
Mod:  
dmattek committed
875
  }) 
876
  
dmattek's avatar
dmattek committed
877
  
878 879 880
  # prepare data with stimulation pattern
  # this dataset is displayed underneath of trajectory plot (modules/trajPlot.R) as geom_segment
  data4stimPlot <- reactive({
881
    if (DEB)  
882
      cat(file = stdout(), 'server:data4stimPlot\n')
883 884
    
    if (input$chBstim) {
885
      if (DEB)  
886
        cat(file = stdout(), 'server:data4stimPlot: stim not NULL\n')
887 888 889 890
      
      loc.dt.stim = dataLoadStim()
      return(loc.dt.stim)
    } else {
891
      if (DEB)  
892
        cat(file = stdout(), 'server:data4stimPlot: stim is NULL\n')
893
      
894 895 896 897
      return(NULL)
    }
  })
  
dmattek's avatar
dmattek committed
898 899 900 901 902 903 904 905 906 907
  # prepare y-axis label in time series plots, depending on UI setting
  
  createYaxisLabel = reactive({
    locLabel = input$inSelMeas1
    
    
    
    return(locLabel)
  })
  
dmattek's avatar
dmattek committed
908 909 910
  # download data as prepared for plotting
  # after all modification
  output$downloadDataClean <- downloadHandler(
dmattek's avatar
dmattek committed
911
    filename = FCSVTCCLEAN,
dmattek's avatar
dmattek committed
912
    content = function(file) {
dmattek's avatar
dmattek committed
913
      write.csv(data4trajPlotNoOut(), file, row.names = FALSE)
dmattek's avatar
dmattek committed
914 915 916
    }
  )
  
dmattek's avatar
dmattek committed
917 918 919
  # Plotting-trajectories ----

  # UI for selecting trajectories
920
  # The output data table of data4trajPlot is modified based on inSelHighlight field
dmattek's avatar
dmattek committed
921
  output$varSelHighlight = renderUI({
922
    if (DEB)  
923
      cat(file = stdout(), 'server:varSelHighlight\n')
dmattek's avatar
dmattek committed
924
    
dmattek's avatar
dmattek committed
925 926 927
    locBut = input$chBhighlightTraj
    if (!locBut)
      return(NULL)
dmattek's avatar
dmattek committed
928
    
dmattek's avatar
dmattek committed
929
    loc.v = getDataTrackObjLabUni()
dmattek's avatar
dmattek committed
930
    if (!is.null(loc.v)) {
931
      selectInput(
dmattek's avatar
dmattek committed
932
        'inSelHighlight',
933
        'Select one or more trajectories:',
dmattek's avatar
dmattek committed
934
        loc.v,
935
        width = '100%',
dmattek's avatar
dmattek committed
936
        multiple = TRUE
937
      )
dmattek's avatar
dmattek committed
938 939 940
    }
  })
  
dmattek's avatar
dmattek committed
941 942 943
  # Taking out outliers 
  data4trajPlotNoOut = callModule(modSelOutliers, 'returnOutlierIDs', data4trajPlot)
  
dmattek's avatar
dmattek committed
944 945
  # Trajectory plotting - ribbon
  callModule(modTrajRibbonPlot, 'modTrajRibbon', 
dmattek's avatar
dmattek committed
946
             in.data = data4trajPlotNoOut,
dmattek's avatar
dmattek committed
947
             in.data.stim = data4stimPlot,
dmattek's avatar
dmattek committed
948
             in.fname = function() return(FPDFTCMEAN))
dmattek's avatar
dmattek committed
949
  
dmattek's avatar
dmattek committed
950
  # Trajectory plotting - individual
dmattek's avatar
dmattek committed
951
  callModule(modTrajPlot, 'modTrajPlot', 
dmattek's avatar
dmattek committed
952
             in.data = data4trajPlotNoOut, 
dmattek's avatar
dmattek committed
953
             in.data.stim = data4stimPlot,
dmattek's avatar
dmattek committed
954 955
             in.fname = function() {return(FPDFTCSINGLE)},
             in.ylab = createYaxisLabel)
dmattek's avatar
dmattek committed
956
  
957 958 959 960 961
  # Trajectory plotting - PSD
  callModule(modPSDPlot, 'modPSDPlot',
             in.data = data4trajPlotNoOut,
             in.fname = function() {return(FPDFTCPSD)})
  
dmattek's avatar
dmattek committed
962 963
  
  # Tabs ----
964
  ###### AUC calculation and plotting
dmattek's avatar
dmattek committed
965
  callModule(tabAUCplot, 'tabAUC', data4trajPlotNoOut, in.fname = function()