server.R 31.4 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
Maciej Dobrzynski's avatar
Maciej Dobrzynski committed
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
Maciej Dobrzynski's avatar
Maciej Dobrzynski committed
23 24

library(sparcl) # sparse hierarchical and k-means
dmattek's avatar
Added:  
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
Added:  
dmattek committed
27
library(imputeTS) # for interpolating NAs
Maciej Dobrzynski's avatar
Maciej Dobrzynski committed
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

Maciej Dobrzynski's avatar
Maciej Dobrzynski committed
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
Added:  
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 ----
dmattek's avatar
dmattek committed
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
Added:  
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
Added:  
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
        accept = c("text/csv", 
                   "text/comma-separated-values,text/plain", 
                   "application/gzip", 
                   "application/bz2"), 
dmattek's avatar
Added:  
dmattek committed
147 148 149 150
      )
  })
  
  output$uiButLoadTrajRem = renderUI({
151
    if (DEB)
152
      cat(file = stdout(), 'server:uiButLoadTrajRem\n')
dmattek's avatar
Added:  
dmattek committed
153 154 155 156 157
    
    if(input$chBtrajRem)
      actionButton("inButLoadTrajRem", "Load Data")
  })

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

dmattek's avatar
dmattek committed
183
  
dmattek's avatar
dmattek committed
184
  # UI-side-panel-column-selection ----
185 186
  
  # Select a column with time series ID
dmattek's avatar
dmattek committed
187
  output$varSelTrackLabel = renderUI({
188
    if (DEB)
189
      cat(file = stdout(), 'server:varSelTrackLabel\n')
190
    
dmattek's avatar
dmattek committed
191
    locCols = getDataNucCols()
192
    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
193 194 195
    
    selectInput(
      'inSelTrackLabel',
dmattek's avatar
dmattek committed
196
      'Track ID column',
dmattek's avatar
dmattek committed
197 198 199 200 201 202
      locCols,
      width = '100%',
      selected = locColSel
    )
  })
  
203
  # Select a column with time
dmattek's avatar
dmattek committed
204
  output$varSelTime = renderUI({
205
    if (DEB)
206
      cat(file = stdout(), 'server:varSelTime\n')
207
    
dmattek's avatar
dmattek committed
208
    locCols = getDataNucCols()
209
    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
210 211 212
    
    selectInput(
      'inSelTime',
dmattek's avatar
dmattek committed
213
      'Time column',
dmattek's avatar
dmattek committed
214 215 216 217 218
      locCols,
      width = '100%',
      selected = locColSel
    )
  })
219

dmattek's avatar
dmattek committed
220
  # This is the main field to select plot facet grouping
dmattek's avatar
dmattek committed
221
  # It's typically a column with the entire experimental description,
dmattek's avatar
dmattek committed
222 223
  # 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
224
  output$varSelGroup = renderUI({
225
    if (DEB)
226
      cat(file = stdout(), 'server:varSelGroup\n')
dmattek's avatar
dmattek committed
227
    
dmattek's avatar
dmattek committed
228 229 230 231 232
    if (input$chBgroup) {
      
      locCols = getDataNucCols()
      
      if (!is.null(locCols)) {
233
        locColSel = locCols[grep('(G|g)roup|(S|s)tim|(S|s)timulation|(S|s)ite', locCols)[1]]
234 235

        #cat('UI varSelGroup::locColSel ', locColSel, '\n')
dmattek's avatar
dmattek committed
236 237
        selectInput(
          'inSelGroup',
dmattek's avatar
dmattek committed
238
          'Grouping columns',
dmattek's avatar
dmattek committed
239 240 241 242 243
          locCols,
          width = '100%',
          selected = locColSel,
          multiple = TRUE
        )
dmattek's avatar
dmattek committed
244 245 246 247
      }
    }
  })
  
248 249
  # UI for selecting grouping to add to track ID to make 
  # the track ID unique across entire dataset
dmattek's avatar
dmattek committed
250
  output$varSelSite = renderUI({
251
    if (DEB)
252
      cat(file = stdout(), 'server:varSelSite\n')
dmattek's avatar
dmattek committed
253
    
254
    if (input$chBtrackUni) {
dmattek's avatar
Added:  
dmattek committed
255
      locCols = getDataNucCols()
Maciej Dobrzynski's avatar
Maciej Dobrzynski committed
256
      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
Added:  
dmattek committed
257 258 259
      
      selectInput(
        'inSelSite',
dmattek's avatar
dmattek committed
260
        'Prepend track ID with',
dmattek's avatar
Added:  
dmattek committed
261 262
        locCols,
        width = '100%',
263 264
        selected = locColSel,
        multiple = T
dmattek's avatar
Added:  
dmattek committed
265 266
      )
    }
dmattek's avatar
dmattek committed
267 268
  })
  
269
  # Select a column with the measurement
dmattek's avatar
dmattek committed
270
  output$varSelMeas1 = renderUI({
271
    if (DEB)
272
      cat(file = stdout(), 'server:varSelMeas1\n')
dmattek's avatar
dmattek committed
273 274 275
    locCols = getDataNucCols()
    
    if (!is.null(locCols)) {
Maciej Dobrzynski's avatar
Maciej Dobrzynski committed
276
      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
277

dmattek's avatar
dmattek committed
278 279
      selectInput(
        'inSelMeas1',
dmattek's avatar
dmattek committed
280
        '1st measurement column',
dmattek's avatar
dmattek committed
281 282 283 284 285 286 287
        locCols,
        width = '100%',
        selected = locColSel
      )
    }
  })
  
288 289
  # Select a column with the 2nd measurement.
  # Some simple operations can be performed betwee the two meaurements
dmattek's avatar
dmattek committed
290
  output$varSelMeas2 = renderUI({
291
    if (DEB)
292
      cat(file = stdout(), 'server:varSelMeas2\n')
293
    
dmattek's avatar
dmattek committed
294 295 296 297
    locCols = getDataNucCols()
    
    if (!is.null(locCols) &&
        !(input$inSelMath %in% c('', '1 / '))) {
Maciej Dobrzynski's avatar
Maciej Dobrzynski committed
298
      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
299

dmattek's avatar
dmattek committed
300 301
      selectInput(
        'inSelMeas2',
dmattek's avatar
dmattek committed
302
        '2nd measurement column',
dmattek's avatar
dmattek committed
303 304 305 306 307 308 309
        locCols,
        width = '100%',
        selected = locColSel
      )
    }
  })
  
dmattek's avatar
dmattek committed
310
  # UI-side-panel-trim x-axis (time) ----
311
  
dmattek's avatar
dmattek committed
312
  output$uiSlTimeTrim = renderUI({
313
    if (DEB)
314
      cat(file = stdout(), 'server:uiSlTimeTrim\n')
dmattek's avatar
dmattek committed
315 316 317 318 319 320 321 322 323 324 325 326
    
    if (input$chBtimeTrim) {
      locTpts  = getDataTpts()
      
      if(is.null(locTpts))
        return(NULL)
      
      locRTmin = min(locTpts)
      locRTmax = max(locTpts)
      
      sliderInput(
        'slTimeTrim',
dmattek's avatar
dmattek committed
327
        label = 'Use time range',
dmattek's avatar
dmattek committed
328 329 330 331 332 333 334
        min = locRTmin,
        max = locRTmax,
        value = c(locRTmin, locRTmax),
        step = 1
      )
      
    }
335 336 337 338 339 340 341
  }) 
  
  # Return the value of slider for trimming time; 
  # output delayed by MILLIS
  returnValSlTimeTrim = reactive({
    return(input$slTimeTrim)
  }) %>% debounce(MILLIS)
dmattek's avatar
dmattek committed
342
  
343 344 345 346 347 348 349 350 351 352 353 354 355 356 357 358 359 360
  # UI-side-panel-interpolation ----
  # Provide interval between 2 time points (for interpolation of NAs and missing time points)
  output$varSelTimeFreq = renderUI({
    if (DEB)
      cat(file = stdout(), 'server:varSelTimeFreq\n')
    
    if (input$chBtrajInter) {
      numericInput(
        'inSelTimeFreq',
        'Interval between 2 time points',
        min = 0,
        step = 1,
        width = '100%',
        value = 1
      )
    }
  })
  
dmattek's avatar
dmattek committed
361
  # UI-side-panel-normalization ----
362 363 364 365
  
  # 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
366
  output$uiChBnorm = renderUI({
367
    if (DEB)
368
      cat(file = stdout(), 'server:uiChBnorm\n')
dmattek's avatar
dmattek committed
369 370
    
    if (input$chBnorm) {
371
      tagList(
dmattek's avatar
dmattek committed
372 373 374
      radioButtons(
        'rBnormMeth',
        label = 'Select method',
375 376 377
        choices = list('fold-change' = 'mean', 'z-score' = 'z.score'),
        width = "40%"
      ),
dmattek's avatar
dmattek committed
378
      bsTooltip('rBnormMeth', helpText.server[["rBnormMeth"]], placement = "top", trigger = "hover", options = NULL)
dmattek's avatar
dmattek committed
379 380 381 382
      )
    }
  })
  
383
  # select the region of the time series for normalisation
dmattek's avatar
dmattek committed
384
  output$uiSlNorm = renderUI({
385
    if (DEB)
386
      cat(file = stdout(), 'server:uiSlNorm\n')
dmattek's avatar
dmattek committed
387 388 389 390 391 392 393 394 395 396
    
    if (input$chBnorm) {
      locTpts  = getDataTpts()
      
      if(is.null(locTpts))
        return(NULL)
      
      locRTmin = min(locTpts)
      locRTmax = max(locTpts)
      
397
      tagList(
dmattek's avatar
dmattek committed
398 399
      sliderInput(
        'slNormRtMinMax',
400
        label = 'Time span',
dmattek's avatar
dmattek committed
401 402
        min = locRTmin,
        max = locRTmax,
dmattek's avatar
dmattek committed
403 404
        value = c(locRTmin, 0.1 * locRTmax), 
        step = 1
405
      ),
dmattek's avatar
dmattek committed
406
      bsTooltip('slNormRtMinMax', helpText.server[["slNormRtMinMax"]], placement = "top", trigger = "hover", options = NULL)
dmattek's avatar
dmattek committed
407 408 409 410
      )
    }
  })
  
411 412 413 414 415 416 417
  # Return the value of slider for normalisation time; 
  # output delayed by MILLIS
  returnValSlNormRtMinMax = reactive({
    return(input$slNormRtMinMax)
  }) %>% debounce(MILLIS)
  
  
418
  # use robust stats (median instead of mean, mad instead of sd)
dmattek's avatar
dmattek committed
419
  output$uiChBnormRobust = renderUI({
420
    if (DEB)
421
      cat(file = stdout(), 'server:uiChBnormRobust\n')
dmattek's avatar
dmattek committed
422 423
    
    if (input$chBnorm) {
424
      tagList(
dmattek's avatar
dmattek committed
425 426
      checkboxInput('chBnormRobust',
                    label = 'Robust stats',
427 428
                    FALSE, 
                    width = "40%"),
dmattek's avatar
dmattek committed
429
      bsTooltip('chBnormRobust', helpText.server[["chBnormRobust"]], placement = "top", trigger = "hover", options = NULL)
430
      )
dmattek's avatar
dmattek committed
431 432 433
    }
  })
  
434
  # choose whether normalisation should be calculated for the entire dataset, group, or trajectory
dmattek's avatar
dmattek committed
435
  output$uiChBnormGroup = renderUI({
436
    if (DEB)
437
      cat(file = stdout(), 'server:uiChBnormGroup\n')
dmattek's avatar
dmattek committed
438 439
    
    if (input$chBnorm) {
440
      tagList(
dmattek's avatar
dmattek committed
441
      radioButtons('chBnormGroup',
dmattek's avatar
Mod:  
dmattek committed
442
                   label = 'Normalisation grouping',
443 444
                   choices = list('Entire dataset' = 'none', 'Per group' = 'group', 'Per trajectory' = 'id'), 
                   width = "40%"),
dmattek's avatar
dmattek committed
445
      bsTooltip('chBnormGroup', helpText.server[["chBnormGroup"]], placement = "top", trigger = "hover", options = NULL)
446
      )
dmattek's avatar
dmattek committed
447 448 449 450
    }
  })
  
  
451 452 453 454 455
  # Pop-overs ----
  addPopover(session, 
             "alDataFormat",
             title = "Data format",
             content = helpText.server[["alDataFormat"]],
dmattek's avatar
dmattek committed
456
             trigger = "click")
dmattek's avatar
dmattek committed
457
  
dmattek's avatar
dmattek committed
458

dmattek's avatar
dmattek committed
459
  # Processing-data ----
dmattek's avatar
dmattek committed
460
  
461
  # Obtain data either from an upload or by generating a synthetic dataset
462 463 464 465 466 467 468 469 470 471 472 473
  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
    
474
    # Don't wrap around if(DEB) !!!
475
    cat(
476
      "server:dataInBoth\n   inGen1: ",
477
      locInGen1,
478
      "      prev=",
479
      isolate(counter$dataGen1),
480
      "\n   inDataNuc: ",
481 482 483 484 485 486 487 488 489 490
      locInLoadNuc,
      "   prev=",
      isolate(counter$dataLoadNuc),
      # "\ninDataStim: ",
      # locInLoadStim,
      # "   prev=",
      # isolate(counter$dataLoadStim),
      "\n"
    )
    
491
    # isolate the checks of the counter reactiveValues
492 493
    # as we set the values in this same reactive
    if (locInGen1 != isolate(counter$dataGen1)) {
494
      cat("server:dataInBoth if inDataGen1\n")
495 496 497 498
      dm = dataGen1()
      # no need to isolate updating the counter reactive values!
      counter$dataGen1 <- locInGen1
    } else if (locInLoadNuc != isolate(counter$dataLoadNuc)) {
499
      cat("server:dataInBoth if inDataLoadNuc\n")
500
      dm = dataLoadNuc()
501 502 503
      
      # convert to long format if radio box set to "wide"
      # the input data in long format should contain:
504
      # - the first row with a header: group, track id, time points as columns with numeric header
505 506
      # - consecutive rows with time series, where columns are time points
      if (input$inRbutLongWide == 1) {
507 508 509 510 511 512 513 514 515 516 517 518
        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")
          
519

520 521 522 523 524 525 526 527 528 529 530 531 532 533 534 535 536 537 538 539 540 541 542 543 544 545 546 547 548 549 550 551 552
        } 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")
          }
        }
553 554
      }
      
555 556 557
      # no need to isolate updating the counter reactive values!
      counter$dataLoadNuc <- locInLoadNuc
    } else {
558
      cat("server:dataInBoth else\n")
559 560
      dm = NULL
    }
561
    
562 563 564
    return(dm)
  })
  
565 566 567
  # Return a dt with mods depending on UI options::
  # - an added column with unique track object label created from the existing track id and prepended with columns chosen in the UI
  # - removed track IDs based on a separate file uploaded; the file should contain a single column with a header and unique track IDs
568
  dataMod <- reactive({
569
    if (DEB)
570
      cat(file = stdout(), 'server:dataMod\n')
571
    
572 573
    loc.dt = dataInBoth()
    
dmattek's avatar
dmattek committed
574
    if (is.null(loc.dt))
575 576
      return(NULL)
    
577
    if (input$chBtrackUni) {
578
      # create unique track ID based on columns specified in input$inSelSite field and combine with input$inSelTrackLabel
579
      loc.dt[, (COLIDUNI) := do.call(paste, c(.SD, sep = "_")), .SDcols = c(input$inSelSite, input$inSelTrackLabel) ]
dmattek's avatar
Added:  
dmattek committed
580
    } else {
581
      # Leave track ID provided in the loaded dataset; has to be unique
582
      loc.dt[, (COLIDUNI) := get(input$inSelTrackLabel)]
dmattek's avatar
Added:  
dmattek committed
583 584
    }
    
dmattek's avatar
Added:  
dmattek committed
585 586
    # remove trajectories based on uploaded csv
    if (input$chBtrajRem) {
587
      if (DEB)
588
        cat(file = stdout(), 'server:dataMod: trajRem not NULL\n')
dmattek's avatar
Added:  
dmattek committed
589 590
      
      loc.dt.rem = dataLoadTrajRem()
dmattek's avatar
dmattek committed
591
      loc.dt = loc.dt[!(trackObjectsLabelUni %in% loc.dt.rem[[1]])]
dmattek's avatar
Added:  
dmattek committed
592 593
    }
    
594 595 596 597 598
    return(loc.dt)
  })
  
  # prepare data for plotting time courses
  # returns dt with these columns:
dmattek's avatar
dmattek committed
599
  #    realtime - selected from input
dmattek's avatar
dmattek committed
600
  #    y        - measurement selected from input
dmattek's avatar
dmattek committed
601
  #               (can be a single column or result of an operation on two cols)
602 603
  #    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
604 605
  #    group    - grouping variable for facetting from input
  #    mid.in   - column with trajectory selection status from the input file or
606 607 608 609
  #               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
610
  dataLong <- reactive({
611
    if (DEB)
612
      cat(file = stdout(), 'server:dataLong\n')
613 614
    
    loc.dt = dataMod()
dmattek's avatar
dmattek committed
615
    if (is.null(loc.dt))
616 617
      return(NULL)
    
618
    # create expression for 'y' column based on measurements and math operations selected in UI
dmattek's avatar
dmattek committed
619
    if (input$inSelMath == '')
620 621 622 623 624 625
      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)
    
626
    # create expression for 'group' column
627 628
    # creates a merged column based on other columns from input
    # used for grouping of plot facets
dmattek's avatar
dmattek committed
629 630 631 632 633 634 635 636 637 638 639
    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')"
    }
640
    
dmattek's avatar
dmattek committed
641 642

    # column name with time
643 644
    loc.s.rt = input$inSelTime
    
dmattek's avatar
dmattek committed
645 646
    # Assign tracks selected for highlighting in UI
    loc.tracks.highlight = input$inSelHighlight
647
    locButHighlight = input$chBhighlightTraj
dmattek's avatar
dmattek committed
648
    
dmattek's avatar
Added:  
dmattek committed
649 650
    
    # Find column names with position
651
    loc.s.pos.x = names(loc.dt)[grep('(L|l)ocation.*X|(P|p)os.x|(P|p)osx', names(loc.dt))[1]]
652
    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
Added:  
dmattek committed
653
    
654
    if (DEB)
655
      cat('server:dataLong:\n   Position columns: ', loc.s.pos.x, loc.s.pos.y, '\n')
656 657
    
    if (!is.na(loc.s.pos.x) & !is.na(loc.s.pos.y))
dmattek's avatar
Added:  
dmattek committed
658 659 660 661
      locPos = TRUE
    else
      locPos = FALSE
    
662 663 664 665
    
    # Find column names with ObjectNumber
    # This is different from TrackObject_Label and is handy to keep
    # because labels on segmented images are typically ObjectNumber
666
    loc.s.objnum = names(loc.dt)[grep('(O|o)bject(N|n)umber', names(loc.dt))[1]]
667
    #cat('dataLong::loc.s.objnum ', loc.s.objnum, '\n')
668 669 670 671
    if (is.na(loc.s.objnum)) {
      locObjNum = FALSE
    }
    else {
dmattek's avatar
dmattek committed
672
      loc.s.objnum = loc.s.objnum[1]
673
      locObjNum = TRUE
dmattek's avatar
dmattek committed
674
    }
675 676
    
    
677 678
    # if dataset contains column mid.in with trajectory filtering status,
    # then, include it in plotting
dmattek's avatar
dmattek committed
679
    if (sum(names(loc.dt) %in% COLIN) > 0)
680 681 682 683 684 685
      locMidIn = TRUE
    else
      locMidIn = FALSE
    
    ## Build expression for selecting columns from loc.dt
    # Core columns
dmattek's avatar
dmattek committed
686 687 688 689
    s.colexpr = paste0('.(',  COLY, ' = ', loc.s.y,
                       ', ', COLID, ' = ', COLIDUNI, 
                       ', ', COLGR, ' = ', loc.s.gr,
                       ', ', COLRT, ' = ', loc.s.rt)
690 691
    
    # account for the presence of 'mid.in' column in uploaded data
dmattek's avatar
dmattek committed
692
    # future: choose this column in UI
693 694
    if(locMidIn)
      s.colexpr = paste0(s.colexpr, 
dmattek's avatar
dmattek committed
695
                         ',', COLIN, ' = ', COLIN)
696 697 698 699
    
    # include position x,y columns in uploaded data
    if(locPos)
      s.colexpr = paste0(s.colexpr, 
dmattek's avatar
dmattek committed
700 701
                         ', ', COLPOSX, '= ', loc.s.pos.x,
                         ', ', COLPOSY, '= ', loc.s.pos.y)
702 703 704 705
    
    # include ObjectNumber column
    if(locObjNum)
      s.colexpr = paste0(s.colexpr, 
dmattek's avatar
dmattek committed
706
                         ', ', COLOBJN, ' = ', loc.s.objnum)
707 708 709 710 711 712 713
    
    # 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))]
    
714 715 716 717 718
    # 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))]
    
719 720 721 722 723 724
    
    # 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
725
        loc.out[, mid.in := ifelse(get(COLID) %in% loc.tracks.highlight, 'SELECTED', get(COLIN))]
726
      else
dmattek's avatar
Mod:  
dmattek committed
727
        # add a column with status of track selection
dmattek's avatar
dmattek committed
728
        loc.out[, mid.in := ifelse(get(COLID) %in% loc.tracks.highlight, 'SELECTED', 'NOT SEL')]
729
    }
730
      
dmattek's avatar
dmattek committed
731

732
    ## Interpolate missing data and NA data points
733
    # From: https://stackoverflow.com/questions/28073752/r-how-to-add-rows-for-missing-values-for-unique-group-sequences
734
    # Tracks are interpolated only within first and last time points of every track id
735 736
    # 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.
737
    # Therefore, we cannot rely on that info to get time frequency; user must provide this number!
738
    
739 740
    # required for subsetting downstream
    setkeyv(loc.out, c(COLGR, COLID, COLRT))
741

dmattek's avatar
dmattek committed
742 743
    if (input$chBtrajInter) {
      # check if time between 2 time points provided and greater than 0
744
      if (input$inSelTimeFreq > 0) {
dmattek's avatar
dmattek committed
745 746
        closeAlert(session, "alertTimeFreq0")
        
747 748 749 750 751 752 753 754 755 756 757 758 759 760 761 762 763 764 765
        # NA's may be already present in the dataset'.
        # Interpolate (linear) them with na.interpolate as well
        if(locPos)
          s.cols = c(COLY, COLPOSX, COLPOSY)
        else
          s.cols = c(COLY)
        
        loc.out = LOCinterpolate(loc.out, COLGR, COLID, COLRT, s.cols, input$inSelTimeFreq, T)
        
        # !!! 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
dmattek committed
766 767 768 769 770 771 772 773 774
      } else {
        closeAlert(session, "alertNAsPresent")
        createAlert(session = session, 
                    anchorId = "alertAnchorSidePanelNAsPresent", 
                    alertId = "alertTimeFreq0", 
                    title = "Error",
                    content = helpText.server[["alertTimeFreq0"]], 
                    append = T,
                    style = "danger")
775
      }
dmattek's avatar
dmattek committed
776 777 778 779 780 781 782 783 784 785 786 787 788 789 790 791
    } else
      closeAlert(session, "alertTimeFreq0")
    
    # Check for explicit NA's in the measurement columns
    # Has to be here (and not in dataMod()) because we need to know the name of the measurement column (COLY)
    if (sum(is.na(loc.out[[COLY]])))
      createAlert(session, "alertAnchorSidePanelNAsPresent", 
                  "alertNAsPresent", 
                  title = "Warning",
                  content = helpText.server[["alertNAsPresent"]], 
                  append = T,
                  style = "warning")
    else
      closeAlert(session, "alertNAsPresent")
    
    
dmattek's avatar
Mod:  
dmattek committed
792
    
793
    ## Trim x-axis (time)
dmattek's avatar
dmattek committed
794
    if(input$chBtimeTrim) {
795
      loc.out = loc.out[get(COLRT) >= returnValSlTimeTrim()[[1]] & get(COLRT) <= returnValSlTimeTrim()[[2]] ]
dmattek's avatar
dmattek committed
796
    }
dmattek's avatar
dmattek committed
797
    
798
    ## Normalization
dmattek's avatar
dmattek committed
799
    # F-n normTraj adds additional column with .norm suffix
dmattek's avatar
dmattek committed
800
    if (input$chBnorm) {
dmattek's avatar
dmattek committed
801
      loc.out = LOCnormTraj(
dmattek's avatar
dmattek committed
802
        in.dt = loc.out,
dmattek's avatar
dmattek committed
803 804
        in.meas.col = COLY,
        in.rt.col = COLRT,
805 806
        in.rt.min = returnValSlNormRtMinMax()[1],
        in.rt.max = returnValSlNormRtMinMax()[2],
dmattek's avatar
dmattek committed
807 808 809 810 811
        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
812
      # Column with normalized data is renamed to the original name
813
      # Further code assumes column name y produced by dataLong
Maciej Dobrzynski's avatar
Maciej Dobrzynski committed
814 815
      
      loc.out[, c(COLY) := NULL]
dmattek's avatar
dmattek committed
816
      setnames(loc.out, 'y.norm', COLY)
dmattek's avatar
dmattek committed
817 818 819
    }
    
    return(loc.out)
dmattek's avatar
dmattek committed
820 821
  })
  
dmattek's avatar
dmattek committed
822
  
823 824 825 826 827
  # Prepare data in wide format, ready for distance calculation in clustering
  # Return a matrix with:
  # - time series as rows
  # - time points as columns
  dataWide <- reactive({
828
    if (DEB)  
829
      cat(file = stdout(), 'server:dataWide\n')
dmattek's avatar
dmattek committed
830
    
831
    loc.dt = dataLongNoOut()
dmattek's avatar
dmattek committed
832 833 834
    if (is.null(loc.dt))
      return(NULL)
    
835 836 837 838
    # 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
839
    
840 841
    # store row names for later
    loc.rownames = loc.dt.wide[[COLID]]
dmattek's avatar
Mod:  
dmattek committed
842
    
843 844
    # omit first column that contains row names
    loc.m.out = as.matrix(loc.dt.wide[, -1])
dmattek's avatar
Added:  
dmattek committed
845
    
846 847
    # assign row names to the matrix
    rownames(loc.m.out) = loc.rownames
dmattek's avatar
Added:  
dmattek committed
848
    
849 850 851
    # Check for missing time points
    # Missing rows in the long format give rise to NAs during dcast
    # Here, we are not checking for explicit NAs in COLY column
dmattek's avatar
dmattek committed
852 853 854
    if ((sum(is.na(loc.dt[[COLY]])) == 0) & (sum(is.na(loc.dt.wide)) > 0))
      cat(stderr(), helpText.server[["alertNAsPresentLong2WideConv"]])

855
    return(loc.m.out)
dmattek's avatar
Mod:  
dmattek committed
856
  }) 
dmattek's avatar
dmattek committed
857
  
dmattek's avatar
dmattek committed
858
  
859 860 861
  # Prepare data with stimulation pattern
  # This dataset is displayed underneath of trajectory plot (modules/trajPlot.R) as geom_segment
  dataStim <- reactive({
862
    if (DEB)  
863
      cat(file = stdout(), 'server:dataStim\n')
864 865
    
    if (input$chBstim) {
866
      if (DEB)  
867
        cat(file = stdout(), 'server:dataStim: stim not NULL\n')
868 869 870 871
      
      loc.dt.stim = dataLoadStim()
      return(loc.dt.stim)
    } else {
872
      if (DEB)  
873
        cat(file = stdout(), 'server:dataStim: stim is NULL\n')
874
      
875 876 877 878
      return(NULL)
    }
  })
  
879 880 881 882 883 884 885 886 887 888 889 890 891
  # Return all unique track object labels (created in dataMod)
  # Used to display track IDs in UI for trajectory highlighting
  getDataTrackObjLabUni <- reactive({
    if (DEB)
      cat(file = stdout(), 'server:getDataTrackObjLabUni\n')
    
    loc.dt = dataMod()
    
    if (is.null(loc.dt))
      return(NULL)
    else
      return(unique(loc.dt[[COLIDUNI]]))
  })
dmattek's avatar
dmattek committed
892
  
893 894 895 896 897 898 899
  
  # Return all unique time points (real time)
  # Used to set limits of sliders for trimming time and for normalisation
  # These timepoints are from the original dt and aren't affected by trimming of x-axis
  getDataTpts <- reactive({
    if (DEB)
      cat(file = stdout(), 'server:getDataTpts\n')
dmattek's avatar
dmattek committed
900
    
901
    loc.dt = dataMod()
dmattek's avatar
dmattek committed
902
    
903 904 905 906 907 908 909 910 911 912 913 914
    if (is.null(loc.dt))
      return(NULL)
    else
      return(unique(loc.dt[[input$inSelTime]]))
  })
  
  
  # Return column names of the main dt
  # Used to fill UI input fields with a choice of column names
  getDataNucCols <- reactive({
    if (DEB)
      cat(file = stdout(), 'server:getDataNucCols: in\n')
dmattek's avatar
dmattek committed
915
    
916 917 918 919 920 921
    loc.dt = dataInBoth()
    
    if (is.null(loc.dt))
      return(NULL)
    else
      return(colnames(loc.dt))
dmattek's avatar
dmattek committed
922 923
  })
  
924 925 926 927 928 929 930
  # Unfinished f-n!
  # prepare y-axis label in time series plots, depending on UI setting
  createYaxisLabel = reactive({
    locLabel = input$inSelMeas1
    
    return(locLabel)
  })
dmattek's avatar
Added:  
dmattek committed
931
  
dmattek's avatar
dmattek committed
932 933 934
  # Plotting-trajectories ----

  # UI for selecting trajectories
935
  # The output data table of dataLong is modified based on inSelHighlight field
dmattek's avatar
dmattek committed
936
  output$varSelHighlight = renderUI({
937
    if (DEB)  
938
      cat(file = stdout(), 'server:varSelHighlight\n')
dmattek's avatar
dmattek committed
939
    
dmattek's avatar
dmattek committed
940 941 942
    locBut = input$chBhighlightTraj
    if (!locBut)
      return(NULL)
dmattek's avatar
dmattek committed
943
    
dmattek's avatar
dmattek committed
944
    loc.v = getDataTrackObjLabUni()
dmattek's avatar
dmattek committed
945
    if (!is.null(loc.v)) {
946
      selectInput(
dmattek's avatar
dmattek committed
947
        'inSelHighlight',
948
        'Select one or more trajectories:',
dmattek's avatar
dmattek committed
949
        loc.v,
950
        width = '100%',
dmattek's avatar
dmattek committed
951
        multiple = TRUE
952
      )
dmattek's avatar
dmattek committed
953 954 955
    }
  })
  
956 957 958 959 960 961 962 963 964 965 966
  # Modules within main window ----
  
  # download data as prepared for plotting
  # after all modification
  output$downloadDataClean <- downloadHandler(
    filename = FCSVTCCLEAN,
    content = function(file) {
      write.csv(dataLongNoOut(), file, row.names = FALSE)
    }
  )

dmattek's avatar
dmattek committed
967
  # Taking out outliers 
968
  dataLongNoOut = callModule(modSelOutliers, 'returnOutlierIDs', dataLong)
dmattek's avatar
dmattek committed
969
  
dmattek's avatar
dmattek committed
970 971
  # Trajectory plotting - ribbon
  callModule(modTrajRibbonPlot, 'modTrajRibbon', 
972 973
             in.data = dataLongNoOut,
             in.data.stim = dataStim,
974 975
             in.group = COLGR,
             in.group.color = NULL,
dmattek's avatar
dmattek committed
976
             in.fname = function() return(FPDFTCMEAN))
dmattek's avatar
dmattek committed
977
  
dmattek's avatar
dmattek committed
978
  # Trajectory plotting - individual
dmattek's avatar
dmattek committed
979
  callModule(modTrajPlot, 'modTrajPlot', 
980 981
             in.data = dataLongNoOut, 
             in.data.stim = dataStim,
dmattek's avatar
dmattek committed
982 983
             in.fname = function() {return(FPDFTCSINGLE)},
             in.ylab = createYaxisLabel)
dmattek's avatar
dmattek committed
984
  
985 986
  # Trajectory plotting - PSD
  callModule(modPSDPlot, 'modPSDPlot',
987
             in.data = dataLongNoOut,
988 989
             in.fname = function() {return(FPDFTCPSD)})
  
dmattek's avatar
dmattek committed
990 991
  
  # Tabs ----
992
  ###### AUC calculation and plotting
993 994 995
  callModule(tabAUCplot, 'tabAUC', 
             dataLongNoOut, 
             in.fname = function() return(FPDFBOXAUC))
dmattek's avatar
Added:  
dmattek committed
996
  
dmattek's avatar
Added:  
dmattek committed
997
  ###### Box-plot
998