server.R 27 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
26
library(dtw) # for dynamic time warping
library(imputeTS) # for interpolating NAs
Maciej Dobrzynski's avatar
Maciej Dobrzynski committed
27
28
29
30
library(robust) # for robust linear regression
library(MASS)
library(pracma) # for trapz

dmattek's avatar
dmattek committed
31

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

dmattek's avatar
dmattek committed
36
37
38
# colour of loader spinner (shinycssloaders)
options(spinner.color="#00A8AA")

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

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

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

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

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

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

  output$varSelTimeFreq = renderUI({
205
    if (DEB)
206
      cat(file = stdout(), 'server:varSelTimeFreq\n')
207
    
208
209
210
    if (input$chBtrajInter) {
      numericInput(
        'inSelTimeFreq',
dmattek's avatar
dmattek committed
211
        'Interval between two time points:',
212
213
214
215
216
217
        min = 1,
        step = 1,
        width = '100%',
        value = 1
      )
    }
218
  })
dmattek's avatar
dmattek committed
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',
238
          'Select columns for plot grouping:',
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',
260
        'Select grouping columns to add to track label:',
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
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',
280
        'Select 1st meas.:',
dmattek's avatar
dmattek committed
281
282
283
284
285
286
287
288
289
        locCols,
        width = '100%',
        selected = locColSel
      )
    }
  })
  
  
  output$varSelMeas2 = renderUI({
290
    if (DEB)
291
      cat(file = stdout(), 'server:varSelMeas2\n')
292
    
dmattek's avatar
dmattek committed
293
294
295
296
    locCols = getDataNucCols()
    
    if (!is.null(locCols) &&
        !(input$inSelMath %in% c('', '1 / '))) {
Maciej Dobrzynski's avatar
Maciej Dobrzynski committed
297
      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
298

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

dmattek's avatar
dmattek committed
420
  # Processing-data ----
dmattek's avatar
dmattek committed
421
  
422
423
424
425
426
427
428
429
430
431
432
433
  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
    
434
    # Don't wrap around if(DEB) !!!
435
    cat(
436
      "server:dataInBoth\n   inGen1: ",
437
      locInGen1,
438
      "      prev=",
439
      isolate(counter$dataGen1),
440
      "\n   inDataNuc: ",
441
442
443
444
445
446
447
448
449
450
      locInLoadNuc,
      "   prev=",
      isolate(counter$dataLoadNuc),
      # "\ninDataStim: ",
      # locInLoadStim,
      # "   prev=",
      # isolate(counter$dataLoadStim),
      "\n"
    )
    
451
    # isolate the checks of the counter reactiveValues
452
453
    # as we set the values in this same reactive
    if (locInGen1 != isolate(counter$dataGen1)) {
454
      cat("server:dataInBoth if inDataGen1\n")
455
456
457
458
      dm = dataGen1()
      # no need to isolate updating the counter reactive values!
      counter$dataGen1 <- locInGen1
    } else if (locInLoadNuc != isolate(counter$dataLoadNuc)) {
459
      cat("server:dataInBoth if inDataLoadNuc\n")
460
      dm = dataLoadNuc()
461
462
463
464
465
466
467
468
469
470
471
472
473
      
      # convert to long format if radio box set to "wide"
      # the input data in long format should contain:
      # - the first row with a header: ID, 1, 2, 3...
      # - consecutive rows with time series, where columns are time points
      if (input$inRbutLongWide == 1) {
        # long to wide
        dm = melt(dm, id.vars = names(dm)[1], variable.name = COLRT, value.name = COLY)

        # convert column names with time points to a number
        dm[, (COLRT) := as.numeric(levels(get(COLRT)))[get(COLRT)]]
      }
      
474
475
476
      # no need to isolate updating the counter reactive values!
      counter$dataLoadNuc <- locInLoadNuc
    } else {
477
      cat("server:dataInBoth else\n")
478
479
480
481
482
483
      dm = NULL
    }
    return(dm)
  })
  
  # return column names of the main dt
dmattek's avatar
dmattek committed
484
  getDataNucCols <- reactive({
485
    if (DEB)
486
      cat(file = stdout(), 'server:getDataNucCols: in\n')
487
    
488
489
490
491
492
493
494
495
496
497
    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({
498
    if (DEB)
499
      cat(file = stdout(), 'server:dataMod\n')
500
    
501
502
    loc.dt = dataInBoth()
    
dmattek's avatar
dmattek committed
503
    if (is.null(loc.dt))
504
505
      return(NULL)
    
506
    if (input$chBtrackUni) {
507
      # create unique track ID based on columns specified in input$inSelSite field and combine with input$inSelTrackLabel
508
      loc.dt[, (COLIDUNI) := do.call(paste, c(.SD, sep = "_")), .SDcols = c(input$inSelSite, input$inSelTrackLabel) ]
dmattek's avatar
Added:    
dmattek committed
509
    } else {
510
      # stay with track ID provided in the loaded dataset; has to be unique
511
      loc.dt[, (COLIDUNI) := get(input$inSelTrackLabel)]
dmattek's avatar
Added:    
dmattek committed
512
513
    }
    
dmattek's avatar
dmattek committed
514
    
dmattek's avatar
Added:    
dmattek committed
515
516
517
    # remove trajectories based on uploaded csv

    if (input$chBtrajRem) {
518
      if (DEB)
519
        cat(file = stdout(), 'server:dataMod: trajRem not NULL\n')
dmattek's avatar
Added:    
dmattek committed
520
521
      
      loc.dt.rem = dataLoadTrajRem()
dmattek's avatar
dmattek committed
522
      loc.dt = loc.dt[!(trackObjectsLabelUni %in% loc.dt.rem[[1]])]
dmattek's avatar
Added:    
dmattek committed
523
524
    }
    
525
526
527
    return(loc.dt)
  })
  
dmattek's avatar
dmattek committed
528
529
530
  # return all unique track object labels (created in dataMod)
  # This will be used to display in UI for trajectory highlighting
  getDataTrackObjLabUni <- reactive({
531
    if (DEB)
532
      cat(file = stdout(), 'server:getDataTrackObjLabUni\n')
533
    
dmattek's avatar
dmattek committed
534
    loc.dt = dataMod()
535
    
dmattek's avatar
dmattek committed
536
537
538
539
    if (is.null(loc.dt))
      return(NULL)
    else
      return(unique(loc.dt$trackObjectsLabelUni))
540
541
  })
  
dmattek's avatar
Mod:    
dmattek committed
542
  
dmattek's avatar
dmattek committed
543
544
545
  # 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
546
  getDataTpts <- reactive({
547
    if (DEB)
548
      cat(file = stdout(), 'server:getDataTpts\n')
549
    
dmattek's avatar
dmattek committed
550
    loc.dt = dataMod()
551
    
dmattek's avatar
dmattek committed
552
553
554
555
    if (is.null(loc.dt))
      return(NULL)
    else
      return(unique(loc.dt[[input$inSelTime]]))
556
557
  })
  
dmattek's avatar
dmattek committed
558
  
559
560
561
  
  # prepare data for plotting time courses
  # returns dt with these columns:
dmattek's avatar
dmattek committed
562
  #    realtime - selected from input
dmattek's avatar
dmattek committed
563
  #    y        - measurement selected from input
dmattek's avatar
dmattek committed
564
  #               (can be a single column or result of an operation on two cols)
565
566
  #    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
567
568
  #    group    - grouping variable for facetting from input
  #    mid.in   - column with trajectory selection status from the input file or
569
570
571
572
  #               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
573
  data4trajPlot <- reactive({
574
    if (DEB)
575
      cat(file = stdout(), 'server:data4trajPlot\n')
576
577
    
    loc.dt = dataMod()
dmattek's avatar
dmattek committed
578
    if (is.null(loc.dt))
579
580
      return(NULL)
    
581
    # create expression for 'y' column based on measurements and math operations selected in UI
dmattek's avatar
dmattek committed
582
    if (input$inSelMath == '')
583
584
585
586
587
588
      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)
    
589
    # create expression for 'group' column
590
591
    # creates a merged column based on other columns from input
    # used for grouping of plot facets
dmattek's avatar
dmattek committed
592
593
594
595
596
597
598
599
600
601
602
    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')"
    }
603
    
dmattek's avatar
dmattek committed
604
605

    # column name with time
606
607
    loc.s.rt = input$inSelTime
    
dmattek's avatar
dmattek committed
608
609
    # Assign tracks selected for highlighting in UI
    loc.tracks.highlight = input$inSelHighlight
610
    locButHighlight = input$chBhighlightTraj
dmattek's avatar
dmattek committed
611
    
dmattek's avatar
Added:    
dmattek committed
612
613
    
    # Find column names with position
614
    loc.s.pos.x = names(loc.dt)[grep('(L|l)ocation.*X|(P|p)os.x|(P|p)osx', names(loc.dt))[1]]
615
    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
616
    
617
    if (DEB)
618
      cat('server:data4trajPlot:\n   Position columns: ', loc.s.pos.x, loc.s.pos.y, '\n')
619
620
    
    if (!is.na(loc.s.pos.x) & !is.na(loc.s.pos.y))
dmattek's avatar
Added:    
dmattek committed
621
622
623
624
      locPos = TRUE
    else
      locPos = FALSE
    
625
626
627
628
    
    # Find column names with ObjectNumber
    # This is different from TrackObject_Label and is handy to keep
    # because labels on segmented images are typically ObjectNumber
629
630
631
632
633
634
    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
635
      loc.s.objnum = loc.s.objnum[1]
636
      locObjNum = TRUE
dmattek's avatar
dmattek committed
637
    }
638
639
    
    
640
641
    # if dataset contains column mid.in with trajectory filtering status,
    # then, include it in plotting
dmattek's avatar
dmattek committed
642
    if (sum(names(loc.dt) %in% COLIN) > 0)
643
644
645
646
647
648
      locMidIn = TRUE
    else
      locMidIn = FALSE
    
    ## Build expression for selecting columns from loc.dt
    # Core columns
dmattek's avatar
dmattek committed
649
650
651
652
    s.colexpr = paste0('.(',  COLY, ' = ', loc.s.y,
                       ', ', COLID, ' = ', COLIDUNI, 
                       ', ', COLGR, ' = ', loc.s.gr,
                       ', ', COLRT, ' = ', loc.s.rt)
653
654
    
    # account for the presence of 'mid.in' column in uploaded data
dmattek's avatar
dmattek committed
655
    # future: choose this column in UI
656
657
    if(locMidIn)
      s.colexpr = paste0(s.colexpr, 
dmattek's avatar
dmattek committed
658
                         ',', COLIN, ' = ', COLIN)
659
660
661
662
    
    # include position x,y columns in uploaded data
    if(locPos)
      s.colexpr = paste0(s.colexpr, 
dmattek's avatar
dmattek committed
663
664
                         ', ', COLPOSX, '= ', loc.s.pos.x,
                         ', ', COLPOSY, '= ', loc.s.pos.y)
665
666
667
668
    
    # include ObjectNumber column
    if(locObjNum)
      s.colexpr = paste0(s.colexpr, 
dmattek's avatar
dmattek committed
669
                         ', ', COLOBJN, ' = ', loc.s.objnum)
670
671
672
673
674
675
676
    
    # 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))]
    
677
678
679
680
681
    # 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))]
    
682
683
684
685
686
687
    
    # 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
688
        loc.out[, mid.in := ifelse(get(COLID) %in% loc.tracks.highlight, 'SELECTED', get(COLIN))]
689
      else
dmattek's avatar
Mod:    
dmattek committed
690
        # add a column with status of track selection
dmattek's avatar
dmattek committed
691
        loc.out[, mid.in := ifelse(get(COLID) %in% loc.tracks.highlight, 'SELECTED', 'NOT SEL')]
692
    }
693
      
dmattek's avatar
dmattek committed
694

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

704
705
    if (input$chBtrajInter) {
      # here we fill missing data with NA's
dmattek's avatar
dmattek committed
706
      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'))]
707
708
      
      # x-check: print all rows with NA's
709
      if (DEB) {
710
        cat(file = stdout(), 'server:data4trajPlot: Rows with NAs:\n')
711
712
        print(loc.out[rowSums(is.na(loc.out)) > 0, ])
      }
713
714
715
716
      
      # NA's may be already present in the dataset'.
      # Interpolate (linear) them with na.interpolate as well
      if(locPos)
dmattek's avatar
dmattek committed
717
        s.cols = c(COLY, COLPOSX, COLPOSY)
718
      else
dmattek's avatar
dmattek committed
719
        s.cols = c(COLY)
720
      
721
722
723
724
725
726
727
728
729
730
731
732
      # 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]
733
734
735
736
737
738
739
740
741
742
743
744
745
746
      
      
      # !!! 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
747
    
748
    ## Trim x-axis (time)
dmattek's avatar
dmattek committed
749
    if(input$chBtimeTrim) {
dmattek's avatar
dmattek committed
750
      loc.out = loc.out[get(COLRT) >= input$slTimeTrim[[1]] & get(COLRT) <= input$slTimeTrim[[2]] ]
dmattek's avatar
dmattek committed
751
    }
dmattek's avatar
dmattek committed
752
    
753
    ## Normalization
dmattek's avatar
dmattek committed
754
    # F-n normTraj adds additional column with .norm suffix
dmattek's avatar
dmattek committed
755
    if (input$chBnorm) {
dmattek's avatar
dmattek committed
756
      loc.out = LOCnormTraj(
dmattek's avatar
dmattek committed
757
        in.dt = loc.out,
dmattek's avatar
dmattek committed
758
759
        in.meas.col = COLY,
        in.rt.col = COLRT,
dmattek's avatar
dmattek committed
760
761
762
763
764
765
766
        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
767
768
      # 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
769
770
      
      loc.out[, c(COLY) := NULL]
dmattek's avatar
dmattek committed
771
      setnames(loc.out, 'y.norm', COLY)
dmattek's avatar
dmattek committed
772
773
774
    }
    
    return(loc.out)
dmattek's avatar
dmattek committed
775
776
  })
  
dmattek's avatar
dmattek committed
777
778
779
780
781
782
  
  # prepare data for clustering
  # return a matrix with:
  # cells as columns
  # time points as rows
  data4clust <- reactive({
783
    if (DEB)  
784
      cat(file = stdout(), 'server:data4clust\n')
dmattek's avatar
dmattek committed
785
    
dmattek's avatar
dmattek committed
786
    loc.dt = data4trajPlotNoOut()
dmattek's avatar
dmattek committed
787
788
789
    if (is.null(loc.dt))
      return(NULL)
    
790
791
792
793
    # 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
794
    
795
796
    # store row names for later
    loc.rownames = loc.dt.wide[[COLID]]
dmattek's avatar
Mod:    
dmattek committed
797
    
798
799
    # omit first column that contains row names
    loc.m.out = as.matrix(loc.dt.wide[, -1])
dmattek's avatar
Added:    
dmattek committed
800
    
801
802
    # assign row names to the matrix
    rownames(loc.m.out) = loc.rownames
dmattek's avatar
Added:    
dmattek committed
803
    
804
    return(loc.m.out)
dmattek's avatar
Mod:    
dmattek committed
805
  }) 
dmattek's avatar
dmattek committed
806
  
dmattek's avatar
dmattek committed
807
  
808
809
810
  # prepare data with stimulation pattern
  # this dataset is displayed underneath of trajectory plot (modules/trajPlot.R) as geom_segment
  data4stimPlot <- reactive({
811
    if (DEB)  
812
      cat(file = stdout(), 'server:data4stimPlot\n')
813
814
    
    if (input$chBstim) {
815
      if (DEB)  
816
        cat(file = stdout(), 'server:data4stimPlot: stim not NULL\n')
817
818
819
820
      
      loc.dt.stim = dataLoadStim()
      return(loc.dt.stim)
    } else {
821
      if (DEB)  
822
        cat(file = stdout(), 'server:data4stimPlot: stim is NULL\n')
823
      
824
825
826
827
      return(NULL)
    }
  })
  
dmattek's avatar
dmattek committed
828
829
830
831
832
833
834
835
836
837
  # 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
838
839
840
  # download data as prepared for plotting
  # after all modification
  output$downloadDataClean <- downloadHandler(
dmattek's avatar
dmattek committed
841
    filename = FCSVTCCLEAN,
dmattek's avatar
Added:    
dmattek committed
842
    content = function(file) {
dmattek's avatar
dmattek committed
843
      write.csv(data4trajPlotNoOut(), file, row.names = FALSE)
dmattek's avatar
Added:    
dmattek committed
844
845
846
    }
  )
  
dmattek's avatar
dmattek committed
847
848
849
  # Plotting-trajectories ----

  # UI for selecting trajectories
850
  # The output data table of data4trajPlot is modified based on inSelHighlight field
dmattek's avatar
dmattek committed
851
  output$varSelHighlight = renderUI({
852
    if (DEB)  
853
      cat(file = stdout(), 'server:varSelHighlight\n')
dmattek's avatar
dmattek committed
854
    
dmattek's avatar
dmattek committed
855
856
857
    locBut = input$chBhighlightTraj
    if (!locBut)
      return(NULL)
dmattek's avatar
dmattek committed
858
    
dmattek's avatar
dmattek committed
859
    loc.v = getDataTrackObjLabUni()
dmattek's avatar
dmattek committed
860
    if (!is.null(loc.v)) {
861
      selectInput(
dmattek's avatar
dmattek committed
862
        'inSelHighlight',
863
        'Select one or more trajectories:',
dmattek's avatar
dmattek committed
864
        loc.v,
865
        width = '100%',
dmattek's avatar
dmattek committed
866
        multiple = TRUE
867
      )
dmattek's avatar
dmattek committed
868
869
870
    }
  })
  
dmattek's avatar
dmattek committed
871
872
873
  # Taking out outliers 
  data4trajPlotNoOut = callModule(modSelOutliers, 'returnOutlierIDs', data4trajPlot)
  
dmattek's avatar
dmattek committed
874
875
  # Trajectory plotting - ribbon
  callModule(modTrajRibbonPlot, 'modTrajRibbon', 
dmattek's avatar
dmattek committed
876
             in.data = data4trajPlotNoOut,
dmattek's avatar
dmattek committed
877
             in.data.stim = data4stimPlot,
dmattek's avatar
dmattek committed
878
             in.fname = function() return(FPDFTCMEAN))
dmattek's avatar
dmattek committed
879
  
dmattek's avatar
dmattek committed
880
  # Trajectory plotting - individual
dmattek's avatar
dmattek committed
881
  callModule(modTrajPlot, 'modTrajPlot', 
dmattek's avatar
dmattek committed
882
             in.data = data4trajPlotNoOut, 
dmattek's avatar
dmattek committed
883
             in.data.stim = data4stimPlot,
dmattek's avatar
dmattek committed
884
885
             in.fname = function() {return(FPDFTCSINGLE)},
             in.ylab = createYaxisLabel)
dmattek's avatar
dmattek committed
886
  
887
888
889
890
891
  # Trajectory plotting - PSD
  callModule(modPSDPlot, 'modPSDPlot',
             in.data = data4trajPlotNoOut,
             in.fname = function() {return(FPDFTCPSD)})
  
dmattek's avatar
dmattek committed