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
})