server.R 24.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 10 11
#

library(shiny)
library(shinyjs) #http://deanattali.com/shinyjs/
library(data.table)
library(ggplot2)
dmattek's avatar
dmattek committed
12
library(gplots) # for heatmap.2
Maciej Dobrzynski's avatar
Maciej Dobrzynski committed
13 14 15
library(plotly) # interactive plot
library(DT) # interactive tables

dmattek's avatar
dmattek committed
16
library(dendextend) # for color_branches
17
library(colorspace) # for palettes (used to colour dendrogram)
dmattek's avatar
dmattek committed
18 19
library(RColorBrewer)
library(scales) # for percentages on y scale
Maciej Dobrzynski's avatar
Maciej Dobrzynski committed
20 21

library(sparcl) # sparse hierarchical and k-means
dmattek's avatar
Added:  
dmattek committed
22 23
library(dtw) # for dynamic time warping
library(imputeTS) # for interpolating NAs
Maciej Dobrzynski's avatar
Maciej Dobrzynski committed
24 25 26 27
library(robust) # for robust linear regression
library(MASS)
library(pracma) # for trapz

dmattek's avatar
dmattek committed
28

29
# Global parameters ----
dmattek's avatar
dmattek committed
30
# change to increase the limit of the upload file size
dmattek's avatar
Added:  
dmattek committed
31
options(shiny.maxRequestSize = 200 * 1024 ^ 2)
dmattek's avatar
dmattek committed
32

dmattek's avatar
dmattek committed
33
# Server logic ----
dmattek's avatar
dmattek committed
34
shinyServer(function(input, output, session) {
35
  useShinyjs()
dmattek's avatar
dmattek committed
36
  
37
  # This is only set at session start
dmattek's avatar
dmattek committed
38
  # We use this as a way to determine which input was
39 40
  # clicked in the dataInBoth reactive
  counter <- reactiveValues(
dmattek's avatar
dmattek committed
41 42 43
    # The value of actionButton is the number of times the button is pressed
    dataGen1        = isolate(input$inDataGen1),
    dataLoadNuc     = isolate(input$inButLoadNuc),
44 45
    dataLoadTrajRem = isolate(input$inButLoadTrajRem),
    dataLoadStim    = isolate(input$inButLoadStim)
dmattek's avatar
dmattek committed
46
  )
dmattek's avatar
dmattek committed
47 48 49 50 51 52 53 54 55

  nCellsCounter <- reactiveValues(
    nCellsOrig = 0,
    nCellsAfterOutlierTrim = 0
  )
    
  myReactVals = reactiveValues(
    outlierIDs = NULL
  )
dmattek's avatar
dmattek committed
56
  
dmattek's avatar
dmattek committed
57
  # UI-side-panel-data-load ----
dmattek's avatar
dmattek committed
58
  
dmattek's avatar
dmattek committed
59
  # Generate random dataset
60
  dataGen1 <- eventReactive(input$inDataGen1, {
61 62
    if (DEB)
      cat("dataGen1\n")
63
    
dmattek's avatar
dmattek committed
64
    return(LOCgenTraj(in.nwells = 3, in.addout = 3))
65 66
  })
  
dmattek's avatar
dmattek committed
67
  # Load main data file
68
  dataLoadNuc <- eventReactive(input$inButLoadNuc, {
69 70 71
    if (DEB)
      cat("dataLoadNuc\n")

72 73 74 75 76 77 78 79 80 81 82
    locFilePath = input$inFileLoadNuc$datapath
    
    counter$dataLoadNuc <- input$inButLoadNuc - 1
    
    if (is.null(locFilePath) || locFilePath == '')
      return(NULL)
    else {
      return(fread(locFilePath))
    }
  })
  
dmattek's avatar
dmattek committed
83 84 85 86
  # This button will reset the inFileLoad
  observeEvent(input$butReset, {
    reset("inFileLoadNuc")  # reset is a shinyjs function
  })
87

dmattek's avatar
dmattek committed
88
  # Load data with trajectories to remove
89
  dataLoadTrajRem <- eventReactive(input$inButLoadTrajRem, {
90 91 92
    if (DEB)
      cat(file = stdout(), "dataLoadTrajRem\n")
    
93 94 95 96 97 98 99 100 101 102
    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
103
  
dmattek's avatar
dmattek committed
104
  # Load data with stimulation pattern
105
  dataLoadStim <- eventReactive(input$inButLoadStim, {
106 107 108
    if (DEB)
      cat(file = stdout(), "dataLoadStim\n")
    
109 110 111 112 113 114 115 116 117 118 119 120
    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
121 122
  # UI for loading csv with cell IDs for trajectory removal
  output$uiFileLoadTrajRem = renderUI({
123 124
    if (DEB)
      cat(file = stdout(), 'UI uiFileLoadTrajRem\n')
dmattek's avatar
Added:  
dmattek committed
125 126 127 128
    
    if(input$chBtrajRem) 
      fileInput(
        'inFileLoadTrajRem',
129
        'Select file and press "Load Data"',
dmattek's avatar
Added:  
dmattek committed
130 131 132 133 134
        accept = c('text/csv', 'text/comma-separated-values,text/plain')
      )
  })
  
  output$uiButLoadTrajRem = renderUI({
135 136
    if (DEB)
      cat(file = stdout(), 'UI uiButLoadTrajRem\n')
dmattek's avatar
Added:  
dmattek committed
137 138 139 140 141
    
    if(input$chBtrajRem)
      actionButton("inButLoadTrajRem", "Load Data")
  })

142 143
  # UI for loading csv with stimulation pattern
  output$uiFileLoadStim = renderUI({
144 145
    if (DEB)
      cat(file = stdout(), 'UI uiFileLoadStim\n')
dmattek's avatar
Added:  
dmattek committed
146
    
147 148 149
    if(input$chBstim) 
      fileInput(
        'inFileLoadStim',
150
        'Select file and press "Load Data"',
151 152 153 154 155
        accept = c('text/csv', 'text/comma-separated-values,text/plain')
      )
  })
  
  output$uiButLoadStim = renderUI({
156 157
    if (DEB)
      cat(file = stdout(), 'UI uiButLoadStim\n')
dmattek's avatar
Added:  
dmattek committed
158
    
159 160
    if(input$chBstim)
      actionButton("inButLoadStim", "Load Data")
dmattek's avatar
Added:  
dmattek committed
161 162
  })
  
163

dmattek's avatar
dmattek committed
164
  
dmattek's avatar
dmattek committed
165
  # UI-side-panel-column-selection ----
dmattek's avatar
dmattek committed
166
  output$varSelTrackLabel = renderUI({
167 168 169
    if (DEB)
      cat(file = stdout(), 'UI varSelTrackLabel\n')
    
dmattek's avatar
dmattek committed
170
    locCols = getDataNucCols()
171
    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
172 173 174
    
    selectInput(
      'inSelTrackLabel',
dmattek's avatar
dmattek committed
175
      'Select ID:',
dmattek's avatar
dmattek committed
176 177 178 179 180 181 182
      locCols,
      width = '100%',
      selected = locColSel
    )
  })
  
  output$varSelTime = renderUI({
183 184 185
    if (DEB)
      cat(file = stdout(), 'UI varSelTime\n')
    
dmattek's avatar
dmattek committed
186
    locCols = getDataNucCols()
187
    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
188 189 190
    
    selectInput(
      'inSelTime',
191
      'Select time column:',
dmattek's avatar
dmattek committed
192 193 194 195 196
      locCols,
      width = '100%',
      selected = locColSel
    )
  })
197 198

  output$varSelTimeFreq = renderUI({
199 200
    if (DEB)
      cat(file = stdout(), 'UI varSelTimeFreq\n')
201
    
202 203 204
    if (input$chBtrajInter) {
      numericInput(
        'inSelTimeFreq',
205
        'Frequency of time units:',
206 207 208 209 210 211
        min = 1,
        step = 1,
        width = '100%',
        value = 1
      )
    }
212
  })
dmattek's avatar
dmattek committed
213
  
dmattek's avatar
dmattek committed
214
  # This is the main field to select plot facet grouping
dmattek's avatar
dmattek committed
215
  # It's typically a column with the entire experimental description,
dmattek's avatar
dmattek committed
216 217
  # 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
218
  output$varSelGroup = renderUI({
219 220
    if (DEB)
      cat(file = stdout(), 'UI varSelGroup\n')
dmattek's avatar
dmattek committed
221
    
dmattek's avatar
dmattek committed
222 223 224 225 226
    if (input$chBgroup) {
      
      locCols = getDataNucCols()
      
      if (!is.null(locCols)) {
227
        locColSel = locCols[grep('(G|g)roup|(S|s)tim|(S|s)timulation|(S|s)ite', locCols)[1]]
228 229

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

dmattek's avatar
dmattek committed
272 273
      selectInput(
        'inSelMeas1',
274
        'Select 1st meas.:',
dmattek's avatar
dmattek committed
275 276 277 278 279 280 281 282 283
        locCols,
        width = '100%',
        selected = locColSel
      )
    }
  })
  
  
  output$varSelMeas2 = renderUI({
284 285 286
    if (DEB)
      cat(file = stdout(), 'UI varSelMeas2\n')
    
dmattek's avatar
dmattek committed
287 288 289 290
    locCols = getDataNucCols()
    
    if (!is.null(locCols) &&
        !(input$inSelMath %in% c('', '1 / '))) {
Maciej Dobrzynski's avatar
Maciej Dobrzynski committed
291
      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
292

dmattek's avatar
dmattek committed
293 294
      selectInput(
        'inSelMeas2',
295
        'Select 2nd meas.',
dmattek's avatar
dmattek committed
296 297 298 299 300 301 302
        locCols,
        width = '100%',
        selected = locColSel
      )
    }
  })
  
dmattek's avatar
dmattek committed
303
  # UI-side-panel-trim x-axis (time) ----
dmattek's avatar
dmattek committed
304
  output$uiSlTimeTrim = renderUI({
305 306
    if (DEB)
      cat(file = stdout(), 'UI uiSlTimeTrim\n')
dmattek's avatar
dmattek committed
307 308 309 310 311 312 313 314 315 316 317 318 319 320 321 322 323 324 325 326 327
    
    if (input$chBtimeTrim) {
      locTpts  = getDataTpts()
      
      if(is.null(locTpts))
        return(NULL)
      
      locRTmin = min(locTpts)
      locRTmax = max(locTpts)
      
      sliderInput(
        'slTimeTrim',
        label = 'Time range to include',
        min = locRTmin,
        max = locRTmax,
        value = c(locRTmin, locRTmax),
        step = 1
      )
      
    }
  })
dmattek's avatar
dmattek committed
328
  
dmattek's avatar
dmattek committed
329
  # UI-side-panel-normalization ----
dmattek's avatar
dmattek committed
330
  output$uiChBnorm = renderUI({
331 332
    if (DEB)
      cat(file = stdout(), 'UI uiChBnorm\n')
dmattek's avatar
dmattek committed
333 334 335 336 337 338 339 340 341 342 343
    
    if (input$chBnorm) {
      radioButtons(
        'rBnormMeth',
        label = 'Select method',
        choices = list('fold-change' = 'mean', 'z-score' = 'z.score')
      )
    }
  })
  
  output$uiSlNorm = renderUI({
344 345
    if (DEB)
      cat(file = stdout(), 'UI uiSlNorm\n')
dmattek's avatar
dmattek committed
346 347 348 349 350 351 352 353 354 355 356 357 358 359 360
    
    if (input$chBnorm) {
      locTpts  = getDataTpts()
      
      if(is.null(locTpts))
        return(NULL)
      
      locRTmin = min(locTpts)
      locRTmax = max(locTpts)
      
      sliderInput(
        'slNormRtMinMax',
        label = 'Time range for norm.',
        min = locRTmin,
        max = locRTmax,
dmattek's avatar
dmattek committed
361 362
        value = c(locRTmin, 0.1 * locRTmax), 
        step = 1
dmattek's avatar
dmattek committed
363 364 365 366 367
      )
    }
  })
  
  output$uiChBnormRobust = renderUI({
368 369
    if (DEB)
      cat(file = stdout(), 'UI uiChBnormRobust\n')
dmattek's avatar
dmattek committed
370 371 372 373 374 375 376 377 378
    
    if (input$chBnorm) {
      checkboxInput('chBnormRobust',
                    label = 'Robust stats',
                    FALSE)
    }
  })
  
  output$uiChBnormGroup = renderUI({
379 380
    if (DEB)
      cat(file = stdout(), 'UI uiChBnormGroup\n')
dmattek's avatar
dmattek committed
381 382 383
    
    if (input$chBnorm) {
      radioButtons('chBnormGroup',
dmattek's avatar
Mod:  
dmattek committed
384
                   label = 'Normalisation grouping',
385
                   choices = list('Entire dataset' = 'none', 'Per facet' = 'group', 'Per trajectory' = 'id'))
dmattek's avatar
dmattek committed
386 387 388 389
    }
  })
  
  
dmattek's avatar
dmattek committed
390
  
dmattek's avatar
dmattek committed
391

dmattek's avatar
dmattek committed
392
  # Processing-data ----
dmattek's avatar
dmattek committed
393
  
394 395 396 397 398 399 400 401 402 403 404 405
  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
    
406
    # Don't wrap around if(DEB)
407 408 409 410 411 412 413 414 415 416 417 418 419 420 421 422 423 424 425 426 427 428 429 430 431 432 433 434 435 436 437 438 439 440 441 442
    cat(
      "dataInBoth\ninGen1: ",
      locInGen1,
      "   prev=",
      isolate(counter$dataGen1),
      "\ninDataNuc: ",
      locInLoadNuc,
      "   prev=",
      isolate(counter$dataLoadNuc),
      # "\ninDataStim: ",
      # locInLoadStim,
      # "   prev=",
      # isolate(counter$dataLoadStim),
      "\n"
    )
    
    # isolate the checks of counter reactiveValues
    # as we set the values in this same reactive
    if (locInGen1 != isolate(counter$dataGen1)) {
      cat("dataInBoth if inDataGen1\n")
      dm = dataGen1()
      # no need to isolate updating the counter reactive values!
      counter$dataGen1 <- locInGen1
    } else if (locInLoadNuc != isolate(counter$dataLoadNuc)) {
      cat("dataInBoth if inDataLoadNuc\n")
      dm = dataLoadNuc()
      # no need to isolate updating the counter reactive values!
      counter$dataLoadNuc <- locInLoadNuc
    } else {
      cat("dataInBoth else\n")
      dm = NULL
    }
    return(dm)
  })
  
  # return column names of the main dt
dmattek's avatar
dmattek committed
443
  getDataNucCols <- reactive({
444 445 446
    if (DEB)
      cat(file = stdout(), 'getDataNucCols: in\n')
    
447 448 449 450 451 452 453 454 455 456
    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({
457 458 459
    if (DEB)
      cat(file = stdout(), 'dataMod\n')
    
460 461
    loc.dt = dataInBoth()
    
dmattek's avatar
dmattek committed
462
    if (is.null(loc.dt))
463 464
      return(NULL)
    
465
    if (input$chBtrackUni) {
466
      # create unique track ID based on columns specified in input$inSelSite field and combine with input$inSelTrackLabel
467
      loc.dt[, (COLIDUNI) := do.call(paste, c(.SD, sep = "_")), .SDcols = c(input$inSelSite, input$inSelTrackLabel) ]
dmattek's avatar
Added:  
dmattek committed
468
    } else {
469
      # stay with track ID provided in the loaded dataset; has to be unique
470
      loc.dt[, (COLIDUNI) := get(input$inSelTrackLabel)]
dmattek's avatar
Added:  
dmattek committed
471 472
    }
    
dmattek's avatar
dmattek committed
473
    
dmattek's avatar
Added:  
dmattek committed
474 475 476
    # remove trajectories based on uploaded csv

    if (input$chBtrajRem) {
477 478
      if (DEB)
        cat(file = stdout(), 'dataMod: trajRem not NULL\n')
dmattek's avatar
Added:  
dmattek committed
479 480
      
      loc.dt.rem = dataLoadTrajRem()
dmattek's avatar
dmattek committed
481
      loc.dt = loc.dt[!(trackObjectsLabelUni %in% loc.dt.rem[[1]])]
dmattek's avatar
Added:  
dmattek committed
482 483
    }
    
484 485 486
    return(loc.dt)
  })
  
dmattek's avatar
dmattek committed
487 488 489
  # return all unique track object labels (created in dataMod)
  # This will be used to display in UI for trajectory highlighting
  getDataTrackObjLabUni <- reactive({
490 491 492
    if (DEB)
      cat(file = stdout(), 'getDataTrackObjLabUni\n')
    
dmattek's avatar
dmattek committed
493
    loc.dt = dataMod()
494
    
dmattek's avatar
dmattek committed
495 496 497 498
    if (is.null(loc.dt))
      return(NULL)
    else
      return(unique(loc.dt$trackObjectsLabelUni))
499 500
  })
  
dmattek's avatar
Mod:  
dmattek committed
501
  
dmattek's avatar
dmattek committed
502 503 504
  # 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
505
  getDataTpts <- reactive({
506 507 508
    if (DEB)
      cat(file = stdout(), 'getDataTpts\n')
    
dmattek's avatar
dmattek committed
509
    loc.dt = dataMod()
510
    
dmattek's avatar
dmattek committed
511 512 513 514
    if (is.null(loc.dt))
      return(NULL)
    else
      return(unique(loc.dt[[input$inSelTime]]))
515 516
  })
  
dmattek's avatar
dmattek committed
517
  
518 519 520
  
  # prepare data for plotting time courses
  # returns dt with these columns:
dmattek's avatar
dmattek committed
521
  #    realtime - selected from input
dmattek's avatar
dmattek committed
522
  #    y        - measurement selected from input
dmattek's avatar
dmattek committed
523
  #               (can be a single column or result of an operation on two cols)
524 525
  #    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
526 527
  #    group    - grouping variable for facetting from input
  #    mid.in   - column with trajectory selection status from the input file or
528 529 530 531
  #               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
532
  data4trajPlot <- reactive({
533 534
    if (DEB)
      cat(file = stdout(), 'data4trajPlot\n')
535 536
    
    loc.dt = dataMod()
dmattek's avatar
dmattek committed
537
    if (is.null(loc.dt))
538 539
      return(NULL)
    
540
    # create expression for 'y' column based on measurements and math operations selected in UI
dmattek's avatar
dmattek committed
541
    if (input$inSelMath == '')
542 543 544 545 546 547
      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)
    
548
    # create expression for 'group' column
549 550
    # creates a merged column based on other columns from input
    # used for grouping of plot facets
dmattek's avatar
dmattek committed
551 552 553 554 555 556 557 558 559 560 561
    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')"
    }
562
    
dmattek's avatar
dmattek committed
563 564

    # column name with time
565 566
    loc.s.rt = input$inSelTime
    
dmattek's avatar
dmattek committed
567 568
    # Assign tracks selected for highlighting in UI
    loc.tracks.highlight = input$inSelHighlight
569
    locButHighlight = input$chBhighlightTraj
dmattek's avatar
dmattek committed
570
    
dmattek's avatar
Added:  
dmattek committed
571 572
    
    # Find column names with position
573
    loc.s.pos.x = names(loc.dt)[grep('(L|l)ocation.*X|(P|p)os.x|(P|p)osx', names(loc.dt))[1]]
574
    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
575
    
576 577
    if (DEB)
      cat('Position columns: ', loc.s.pos.x, loc.s.pos.y, '\n')
578 579
    
    if (!is.na(loc.s.pos.x) & !is.na(loc.s.pos.y))
dmattek's avatar
Added:  
dmattek committed
580 581 582 583
      locPos = TRUE
    else
      locPos = FALSE
    
584 585 586 587
    
    # Find column names with ObjectNumber
    # This is different from TrackObject_Label and is handy to keep
    # because labels on segmented images are typically ObjectNumber
588 589 590 591 592 593
    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
594
      loc.s.objnum = loc.s.objnum[1]
595
      locObjNum = TRUE
dmattek's avatar
dmattek committed
596
    }
597 598
    
    
599 600
    # if dataset contains column mid.in with trajectory filtering status,
    # then, include it in plotting
601 602 603 604 605 606 607 608 609 610 611 612 613 614 615 616 617 618 619 620 621 622 623 624 625 626 627 628 629 630 631 632 633 634 635 636 637 638 639 640 641 642
    if (sum(names(loc.dt) %in% 'mid.in') > 0)
      locMidIn = TRUE
    else
      locMidIn = FALSE
    
    ## Build expression for selecting columns from loc.dt
    # Core columns
    s.colexpr = paste0('.(y = ', loc.s.y,
                       ', id = trackObjectsLabelUni', 
                       ', group = ', loc.s.gr,
                       ', realtime = ', loc.s.rt)
    
    # account for the presence of 'mid.in' column in uploaded data
    if(locMidIn)
      s.colexpr = paste0(s.colexpr, 
                         ', mid.in = mid.in')
    
    # include position x,y columns in uploaded data
    if(locPos)
      s.colexpr = paste0(s.colexpr, 
                         ', pos.x = ', loc.s.pos.x,
                         ', pos.y = ', loc.s.pos.y)
    
    # include ObjectNumber column
    if(locObjNum)
      s.colexpr = paste0(s.colexpr, 
                         ', obj.num = ', loc.s.objnum)
    
    # 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))]
    
    
    # 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)
        loc.out[, mid.in := ifelse(id %in% loc.tracks.highlight, 'SELECTED', mid.in)]
      else
dmattek's avatar
Mod:  
dmattek committed
643
        # add a column with status of track selection
644
        loc.out[, mid.in := ifelse(id %in% loc.tracks.highlight, 'SELECTED', 'NOT SEL')]
645
    }
646
      
dmattek's avatar
dmattek committed
647

648
    ## Interpolate missing data and NA data points
649
    # From: https://stackoverflow.com/questions/28073752/r-how-to-add-rows-for-missing-values-for-unique-group-sequences
650
    # Tracks are interpolated only within first and last time points of every track id
651 652
    # 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.
653
    # Therefore, we cannot rely on that info to get time frequency; user must provide this number!
654
    
655
    setkeyv(loc.out, c(COLGR, COLID, COLRT))
656

657 658
    if (input$chBtrajInter) {
      # here we fill missing data with NA's
dmattek's avatar
dmattek committed
659
      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'))]
660 661
      
      # x-check: print all rows with NA's
662 663 664 665
      if (DEB) {
        cat(file = stdout(), 'Rows with NAs:\n')
        print(loc.out[rowSums(is.na(loc.out)) > 0, ])
      }
666 667 668 669
      
      # NA's may be already present in the dataset'.
      # Interpolate (linear) them with na.interpolate as well
      if(locPos)
dmattek's avatar
dmattek committed
670
        s.cols = c(COLY, COLPOSX, COLPOSY)
671
      else
dmattek's avatar
dmattek committed
672
        s.cols = c(COLY)
673
      
674 675 676 677 678 679 680 681 682 683 684 685
      # 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]]))

        loc.out[, (col) := na.interpolation(get(col)), by = c(COLID)]        
      }
      
      # loc.out[, (s.cols) := lapply(.SD, na.interpolation), by = c(COLID), .SDcols = s.cols]
686 687 688 689 690 691 692 693 694 695 696 697 698 699
      
      
      # !!! 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
700
    
701
    ## Trim x-axis (time)
dmattek's avatar
dmattek committed
702
    if(input$chBtimeTrim) {
dmattek's avatar
dmattek committed
703
      loc.out = loc.out[get(COLRT) >= input$slTimeTrim[[1]] & get(COLRT) <= input$slTimeTrim[[2]] ]
dmattek's avatar
dmattek committed
704
    }
dmattek's avatar
dmattek committed
705
    
706
    ## Normalization
dmattek's avatar
dmattek committed
707
    # F-n normTraj adds additional column with .norm suffix
dmattek's avatar
dmattek committed
708
    if (input$chBnorm) {
dmattek's avatar
dmattek committed
709
      loc.out = LOCnormTraj(
dmattek's avatar
dmattek committed
710
        in.dt = loc.out,
dmattek's avatar
dmattek committed
711 712
        in.meas.col = COLY,
        in.rt.col = COLRT,
dmattek's avatar
dmattek committed
713 714 715 716 717 718 719
        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
720 721
      # 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
722 723
      
      loc.out[, c(COLY) := NULL]
dmattek's avatar
dmattek committed
724
      setnames(loc.out, 'y.norm', COLY)
dmattek's avatar
dmattek committed
725 726 727
    }
    
    return(loc.out)
dmattek's avatar
dmattek committed
728 729
  })
  
dmattek's avatar
dmattek committed
730 731 732 733 734 735
  
  # prepare data for clustering
  # return a matrix with:
  # cells as columns
  # time points as rows
  data4clust <- reactive({
736 737
    if (DEB)  
      cat(file = stdout(), 'data4clust\n')
dmattek's avatar
dmattek committed
738
    
dmattek's avatar
dmattek committed
739
    loc.dt = data4trajPlotNoOut()
dmattek's avatar
dmattek committed
740 741 742
    if (is.null(loc.dt))
      return(NULL)
    
dmattek's avatar
Added:  
dmattek committed
743
    #print(loc.dt)
744
    loc.out = dcast(loc.dt, as.formula(paste0(COLID, "~", COLRT)), value.var = COLY)
dmattek's avatar
Added:  
dmattek committed
745
    #print(loc.out)
746
    loc.rownames = loc.out[[COLID]]
dmattek's avatar
dmattek committed
747
    
dmattek's avatar
Mod:  
dmattek committed
748
    
dmattek's avatar
dmattek committed
749 750
    loc.out = as.matrix(loc.out[, -1])
    rownames(loc.out) = loc.rownames
dmattek's avatar
Added:  
dmattek committed
751
    
752 753
    # This might be removed entirely because all NA treatment happens in data4trajPlot
    # Clustering should work with NAs present. These might result from data itself or from missing time point rows that were turned into NAs when dcast-ing from long format.
dmattek's avatar
Added:  
dmattek committed
754 755 756 757
    # Remove NA's
    # na.interpolation from package imputeTS works with multidimensional data
    # but imputation is performed for each column independently
    # The matrix for clustering contains time series in rows, hence transposing it twice
758
    # loc.out = t(na.interpolation(t(loc.out)))
dmattek's avatar
Added:  
dmattek committed
759
    
dmattek's avatar
dmattek committed
760
    return(loc.out)
dmattek's avatar
Mod:  
dmattek committed
761
  }) 
dmattek's avatar
dmattek committed
762
  
dmattek's avatar
dmattek committed
763
  
764 765 766
  # prepare data with stimulation pattern
  # this dataset is displayed underneath of trajectory plot (modules/trajPlot.R) as geom_segment
  data4stimPlot <- reactive({
767 768
    if (DEB)  
      cat(file = stdout(), 'data4stimPlot\n')
769 770
    
    if (input$chBstim) {
771 772
      if (DEB)  
        cat(file = stdout(), 'data4stimPlot: stim not NULL\n')
773 774 775 776
      
      loc.dt.stim = dataLoadStim()
      return(loc.dt.stim)
    } else {
777 778 779
      if (DEB)  
        cat(file = stdout(), 'data4stimPlot: stim is NULL\n')
      
780 781 782 783
      return(NULL)
    }
  })
  
dmattek's avatar
Added:  
dmattek committed
784 785 786
  # download data as prepared for plotting
  # after all modification
  output$downloadDataClean <- downloadHandler(
dmattek's avatar
dmattek committed
787
    filename = FCSVTCCLEAN,
dmattek's avatar
Added:  
dmattek committed
788
    content = function(file) {
dmattek's avatar
dmattek committed
789
      write.csv(data4trajPlotNoOut(), file, row.names = FALSE)
dmattek's avatar
Added:  
dmattek committed
790 791 792
    }
  )
  
dmattek's avatar
dmattek committed
793 794 795
  # Plotting-trajectories ----

  # UI for selecting trajectories
796
  # The output data table of data4trajPlot is modified based on inSelHighlight field
dmattek's avatar
dmattek committed
797
  output$varSelHighlight = renderUI({
798 799
    if (DEB)  
      cat(file = stdout(), 'UI varSelHighlight\n')
dmattek's avatar
dmattek committed
800
    
dmattek's avatar
dmattek committed
801 802 803
    locBut = input$chBhighlightTraj
    if (!locBut)
      return(NULL)
dmattek's avatar
dmattek committed
804
    
dmattek's avatar
dmattek committed
805
    loc.v = getDataTrackObjLabUni()
dmattek's avatar
dmattek committed
806
    if (!is.null(loc.v)) {
807
      selectInput(
dmattek's avatar
dmattek committed
808 809 810
        'inSelHighlight',
        'Select one or more rajectories:',
        loc.v,
811
        width = '100%',
dmattek's avatar
dmattek committed
812
        multiple = TRUE
813
      )
dmattek's avatar
dmattek committed
814 815 816
    }
  })
  
dmattek's avatar
dmattek committed
817 818 819
  # Taking out outliers 
  data4trajPlotNoOut = callModule(modSelOutliers, 'returnOutlierIDs', data4trajPlot)
  
dmattek's avatar
dmattek committed
820 821
  # Trajectory plotting - ribbon
  callModule(modTrajRibbonPlot, 'modTrajRibbon', 
dmattek's avatar
dmattek committed
822
             in.data = data4trajPlotNoOut,
dmattek's avatar
dmattek committed
823
             in.data.stim = data4stimPlot,
dmattek's avatar
dmattek committed
824
             in.fname = function() return(FPDFTCMEAN))
dmattek's avatar
dmattek committed
825
  
dmattek's avatar
dmattek committed
826
  # Trajectory plotting - individual
dmattek's avatar
dmattek committed
827
  callModule(modTrajPlot, 'modTrajPlot', 
dmattek's avatar
dmattek committed
828
             in.data = data4trajPlotNoOut, 
dmattek's avatar
dmattek committed
829
             in.data.stim = data4stimPlot,
dmattek's avatar
dmattek committed
830
             in.fname = function() {return(FPDFTCSINGLE)})
dmattek's avatar
dmattek committed
831
  
832 833 834 835 836
  # Trajectory plotting - PSD
  callModule(modPSDPlot, 'modPSDPlot',
             in.data = data4trajPlotNoOut,
             in.fname = function() {return(FPDFTCPSD)})
  
dmattek's avatar
dmattek committed
837 838
  
  # Tabs ----
839
  ###### AUC calculation and plotting
dmattek's avatar
dmattek committed
840
  callModule(modAUCplot, 'tabAUC', data4trajPlotNoOut, in.fname = function() return(FPDFBOXAUC))
dmattek's avatar
Added:  
dmattek committed
841
  
dmattek's avatar
Added:  
dmattek committed
842
  ###### Box-plot
dmattek's avatar
dmattek committed
843
  callModule(tabBoxPlot, 'tabBoxPlot', data4trajPlotNoOut, in.fname = function() return(FPDFBOXTP))
dmattek's avatar
dmattek committed
844
  
dmattek's avatar
dmattek committed
845
  ###### Scatter plot
dmattek's avatar
dmattek committed
846
  callModule(tabScatterPlot, 'tabScatter', data4trajPlotNoOut, in.fname = function() return(FPDFSCATTER))
dmattek's avatar
dmattek committed
847
  
dmattek's avatar
dmattek committed
848
  ##### Hierarchical clustering
dmattek's avatar
dmattek committed
849
  callModule(clustHier, 'tabClHier', data4clust, data4trajPlotNoOut, data4stimPlot)
dmattek's avatar
dmattek committed
850 851
  
  ##### Sparse hierarchical clustering using sparcl
dmattek's avatar
dmattek committed
852
  callModule(clustHierSpar, 'tabClHierSpar', data4clust, data4trajPlotNoOut, data4stimPlot)
dmattek's avatar
dmattek committed
853

dmattek's avatar
Mod:  
dmattek committed
854
  
dmattek's avatar
dmattek committed
855
})