server.R 25.5 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
    if (DEB)
62
      cat("server: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
    if (DEB)
70
      cat("server:dataLoadNuc\n")
71

72
73
74
75
76
77
78
    locFilePath = input$inFileLoadNuc$datapath
    
    counter$dataLoadNuc <- input$inButLoadNuc - 1
    
    if (is.null(locFilePath) || locFilePath == '')
      return(NULL)
    else {
79
      return(fread(locFilePath, strip.white = T))
80
81
82
    }
  })
  
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
    if (DEB)
91
      cat(file = stdout(), "server:dataLoadTrajRem\n")
92
    
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
    if (DEB)
107
      cat(file = stdout(), "server:dataLoadStim\n")
108
    
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
    if (DEB)
124
      cat(file = stdout(), 'server: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
    if (DEB)
136
      cat(file = stdout(), 'server: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
    if (DEB)
145
      cat(file = stdout(), 'server: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
    if (DEB)
157
      cat(file = stdout(), 'server: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
    if (DEB)
168
      cat(file = stdout(), 'server:varSelTrackLabel\n')
169
    
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
    if (DEB)
184
      cat(file = stdout(), 'server:varSelTime\n')
185
    
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
    if (DEB)
200
      cat(file = stdout(), 'server: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
    if (DEB)
220
      cat(file = stdout(), 'server: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
    if (DEB)
246
      cat(file = stdout(), 'server: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
    if (DEB)
266
      cat(file = stdout(), 'server: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
    if (DEB)
285
      cat(file = stdout(), 'server:varSelMeas2\n')
286
    
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
    if (DEB)
306
      cat(file = stdout(), 'server: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
    if (DEB)
332
      cat(file = stdout(), 'server: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
    if (DEB)
345
      cat(file = stdout(), 'server: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
    if (DEB)
369
      cat(file = stdout(), 'server: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
    if (DEB)
380
      cat(file = stdout(), 'server: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
    cat(
408
      "server:dataInBoth\n   inGen1: ",
409
      locInGen1,
410
      "      prev=",
411
      isolate(counter$dataGen1),
412
      "\n   inDataNuc: ",
413
414
415
416
417
418
419
420
421
422
      locInLoadNuc,
      "   prev=",
      isolate(counter$dataLoadNuc),
      # "\ninDataStim: ",
      # locInLoadStim,
      # "   prev=",
      # isolate(counter$dataLoadStim),
      "\n"
    )
    
423
    # isolate the checks of the counter reactiveValues
424
425
    # as we set the values in this same reactive
    if (locInGen1 != isolate(counter$dataGen1)) {
426
      cat("server:dataInBoth if inDataGen1\n")
427
428
429
430
      dm = dataGen1()
      # no need to isolate updating the counter reactive values!
      counter$dataGen1 <- locInGen1
    } else if (locInLoadNuc != isolate(counter$dataLoadNuc)) {
431
      cat("server:dataInBoth if inDataLoadNuc\n")
432
      dm = dataLoadNuc()
433
434
435
436
437
438
439
440
441
442
443
444
445
      
      # 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)]]
      }
      
446
447
448
      # no need to isolate updating the counter reactive values!
      counter$dataLoadNuc <- locInLoadNuc
    } else {
449
      cat("server:dataInBoth else\n")
450
451
452
453
454
455
      dm = NULL
    }
    return(dm)
  })
  
  # return column names of the main dt
dmattek's avatar
dmattek committed
456
  getDataNucCols <- reactive({
457
    if (DEB)
458
      cat(file = stdout(), 'server:getDataNucCols: in\n')
459
    
460
461
462
463
464
465
466
467
468
469
    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({
470
    if (DEB)
471
      cat(file = stdout(), 'server:dataMod\n')
472
    
473
474
    loc.dt = dataInBoth()
    
dmattek's avatar
dmattek committed
475
    if (is.null(loc.dt))
476
477
      return(NULL)
    
478
    if (input$chBtrackUni) {
479
      # create unique track ID based on columns specified in input$inSelSite field and combine with input$inSelTrackLabel
480
      loc.dt[, (COLIDUNI) := do.call(paste, c(.SD, sep = "_")), .SDcols = c(input$inSelSite, input$inSelTrackLabel) ]
dmattek's avatar
Added:    
dmattek committed
481
    } else {
482
      # stay with track ID provided in the loaded dataset; has to be unique
483
      loc.dt[, (COLIDUNI) := get(input$inSelTrackLabel)]
dmattek's avatar
Added:    
dmattek committed
484
485
    }
    
dmattek's avatar
dmattek committed
486
    
dmattek's avatar
Added:    
dmattek committed
487
488
489
    # remove trajectories based on uploaded csv

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

    # column name with time
578
579
    loc.s.rt = input$inSelTime
    
dmattek's avatar
dmattek committed
580
581
    # Assign tracks selected for highlighting in UI
    loc.tracks.highlight = input$inSelHighlight
582
    locButHighlight = input$chBhighlightTraj
dmattek's avatar
dmattek committed
583
    
dmattek's avatar
Added:    
dmattek committed
584
585
    
    # Find column names with position
586
    loc.s.pos.x = names(loc.dt)[grep('(L|l)ocation.*X|(P|p)os.x|(P|p)osx', names(loc.dt))[1]]
587
    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
588
    
589
    if (DEB)
590
      cat('server:data4trajPlot:\n   Position columns: ', loc.s.pos.x, loc.s.pos.y, '\n')
591
592
    
    if (!is.na(loc.s.pos.x) & !is.na(loc.s.pos.y))
dmattek's avatar
Added:    
dmattek committed
593
594
595
596
      locPos = TRUE
    else
      locPos = FALSE
    
597
598
599
600
    
    # Find column names with ObjectNumber
    # This is different from TrackObject_Label and is handy to keep
    # because labels on segmented images are typically ObjectNumber
601
602
603
604
605
606
    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
607
      loc.s.objnum = loc.s.objnum[1]
608
      locObjNum = TRUE
dmattek's avatar
dmattek committed
609
    }
610
611
    
    
612
613
    # if dataset contains column mid.in with trajectory filtering status,
    # then, include it in plotting
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
643
644
645
646
647
    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))]
    
648
649
650
651
652
    # 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))]
    
653
654
655
656
657
658
659
660
    
    # 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
661
        # add a column with status of track selection
662
        loc.out[, mid.in := ifelse(id %in% loc.tracks.highlight, 'SELECTED', 'NOT SEL')]
663
    }
664
      
dmattek's avatar
dmattek committed
665

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

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

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

dmattek's avatar
Mod:    
dmattek committed
869
  
dmattek's avatar
dmattek committed
870
})