dashboard.Rmd

---
title: "Runboard"
output:
  flexdashboard::flex_dashboard:
    css: styles.css
    theme: united
    orientation: rows
    social: menu
runtime: shiny
---

```{r setup, include=FALSE}
library(flexdashboard)
library(dplyr)
library(leaflet)
library(tidyr)
library(dygraphs)
library(xts)
library(stringr)
library(elevatr)
library(rgdal)
library(ggplot2)
library(devtools)
library(lubridate)
library(fit) # install_github('octomike/fit')
library(zoo)
library(caret)
library(doFuture)

options(shiny.reactlog = FALSE)


# helpers
source("helpers/data.R")
source("helpers/plots.R")
source("helpers/control.R")

# global variables
rv <- reactiveValues(hrMethod=NULL, hrdpCache=NULL, elevationCache=NULL,
                     mfit=NULL, lmfit=NULL, lmspan=1.0)


methodLabels <- c('Automatic', 'Lactate Threshold', 'Maximum HR')
hlCallbackJs <- 'function(event, x, points, row, seriesName){
  Shiny.setInputValue("currentLocation", row);
  for (div of document.getElementsByClassName("dygraphs")) {
    var obj = window.HTMLWidgets.getInstance(div);
    obj.dygraph.setSelection(row);
  }
}'
```


Sidebar {.sidebar}
==================

### Open activity:

```{r}
fileInput('fitFile', label='Only .fit files are supported')
```


### Heart Rate Zones

```{r}
selectInput('hrmethod', 'HR based on', choices=methodLabels)
renderUI({
  req(rv$hrMethod)
  sliderInput('hrthresh', 'HR in bpm', value=0, min=150, max=220)
})
```

### Elevation Correction

If your device does not have a barometric sensor, try to correct GPS data with a [external radar data](https://cran.r-project.org/web/packages/elevatr/vignettes/introduction_to_elevatr.html).


```{r}
checkboxInput('fixElevation', label='Fix Elevation', value=FALSE)
```

```{r}
# hide inputs with css and use from dygraph JS callbacks
numericInput('currentLocation', '', value=1)
dateInput('rangeStart', '', NULL)
dateInput('rangeStop', '', NULL)
```


```{r}

# rawdata contains $records timeseries and some metadata
rawdata <- reactive({
  req(input$fitFile)
  rv$elevationCache <- NULL # drop the cache
  read.fit(input$fitFile$datapath)
})

# data cleaning and fixing
fitdata <- reactive({
    raw <- rawdata()
    validate(need('heart_rate' %in% colnames(raw$record), 'Fit file does not contain any heart rate data :('))
    if ( str_detect(input$hrmethod, 'Automatic') ) {
      if ( 'threshold_heart_rate' %in% names(raw$zones_target) ) {
        updateSelectInput(session, 'hrmethod', label='HR based on', selected=methodLabels[1])
        updateSliderInput(session, 'hrthresh',label='Lactate HR in bpm (from fit file)',
                          min=150, max=220, value=raw$zones_target$threshold_heart_rate)
        rv$hrMethod <- 'lactate'
      }
      else if ( 'max_heart_rate' %in% names(raw$zones_target) ) {
        updateSelectInput(session, 'hrmethod', label='HR based on', selected=methodLabels[1])
        updateSliderInput(session, 'hrthresh', label='Max HR in bpm (from fit file)',
                          min=150, max=220, value=raw$zones_target$max_heart_rate)
        rv$hrMethod <- 'max'
      }
    }
    else if ( str_detect(input$hrmethod, 'Lactate') ) {
      updateSliderInput(session, 'hrthresh', label='Lactate HR in bpm',
                        min=150, max=220, value=180)
      updateSelectInput(session, 'hrmethod', label='HR based on', selected=methodLabels[2])
      rv$hrMethod <- 'lactate'
    }
    else {
      updateSliderInput(session, 'hrthresh', label='Max HR in bpm',
                        min=150, max=220, value=180)
      updateSelectInput(session, 'hrmethod', label='HR based on', selected=methodLabels[3])
      rv$hrMethod <- 'max'
    }
    df <- raw$record %>%
      run_rename_select %>%
      run_fix_resolution %>%
      run_fix_units %>%
      run_fix_elevation(skip=!input$fixElevation, cache=rv$elevationCache) %>%
      run_add_pace %>%
      run_add_features

    # cache elevation fix
    if (input$fixElevation && is.null(rv$elevationCache)) {
      rv$elevationCache <- df$altitude
    }
    df
})

# limit data to range from dygraph
rangedata <- reactive({
  fitdata() %>%
    run_filter_fitrange(start=input$rangeStart, stop=input$rangeStop) %>%
    run_add_hrzone(breaks=HR()$thresh, labels=HR()$label[1:5])
})

# reduce resolution for a quicker map
mapdata <- reactive({
  rd <- rangedata()
  tres <- rd %>%
            mutate(td=as.numeric(timestamp - lag(timestamp), units='secs')) %>%
            select('td') %>%
            drop_na() %>%
            summarise(mean(td)) %>%
            pull()

  if (tres <= 2){
    rd %>%
      filter(!(row_number() %% max(floor(tres), 4))) # TODO consider making this an input
  } else {
    rd
  }
})

# prepare data for HRDP estimation
analysisdata <- reactive({
  req(input$modelType)
  rd <- rangedata() %>%
      select(c('slope', 'pace', 'hr', 'hrzone', 'hrchange', 'speed')) %>%
      filter(pace < 8) %>%
      filter(hr >= HR()$thresh[2])
  tryCatch({
    d <- rd %>%
      filter(speed >= quantile(speed, probs=c(input$lowerQuantile))) %>%
      group_by(hr) %>%
      filter((abs(speed - median(speed)) * as.numeric(input$removeOutlier)
              <= 2 * sd(speed))) %>%
      mutate(med=median(speed), avg=mean(speed),
             lwr=median(speed)-sd(speed), upr=median(speed) + sd(speed)) %>%
      ungroup
  }, error=function(cond){
    message(cond)
  })
  if (grepl('SVM', input$modelType)) {
    message('fitting SVM linear')
    rv$lmfit <- svm_lin(d, slopes=input$slopes, hrs=input$hrs, C=0.05)
    message(sprintf("error %f\n", rv$lmfit$error))
    message(sprintf("tune: %s\n", rv$lmfit$bestTune))
    message('fitting radial SVM')
    rv$mfit <- svm_radial(d, slopes=input$slopes, hrs=input$hrs,
                          C=2^input$C, sigma=2^input$sigma)
    message(sprintf("error %f\n", rv$mfit$error))
    message(sprintf("tune: %s\n", rv$mfit$bestTune))
  } else {
    message('fitting LOESS/LM models')
    lmfit <- lm(hr ~ speed, data = d)
    lmfit$x <- d$speed
    rv$lmfit <- lmfit
    rv$mfit <- loess(hr ~ speed, data = d, span=rv$lmspan)
  }
  d
})


# depends on maximum hr in range
# used to restrict the extent of HR-Zone-plot
HR <- reactive({
  req(input$hrthresh, input$fitFile, rv$hrMethod)
  print(paste('hrMethod: ', rv$hrMethod))
  thresh <- input$hrthresh
  print(paste('thresh: ', thresh))
  # https://www.researchgate.net/figure/The-5-Zone-Intensity-Scale-Used-by-the-Norwegian-Olympic-Federation-and-the-3-Zone-and_tbl1_259652575
  breaks <- thresh * c(0.72, 0.82, 0.87, 0.92, 1)
  if ( rv$hrMethod == 'lactate' ) {
    breaks <- breaks / 0.92
  }
  breaks <- round(breaks)
  breaks <- c(min(breaks[1]-1, min(fitdata()$hr)), breaks)
  print(paste('breaks: ', breaks))
  # the return object
  data.frame(color=c('#c4ced4', '#b0dcf7', '#b9dfa0', '#fbd89e', '#faa8b1', '#ffffff'),
             label=c('Z1', 'Z2', 'Z3', 'Z4', 'Z5', 'Zmax'),
             thresh=breaks)
})

```

Summary
=======

Row
---

### Pace {.value-box}

```{r}
avgPace <- reactive({
  data <- rangedata()
  d <- ( tail(data$distance, n=1) - head(data$distance, n=1) ) / 1000
  t <- as.numeric(
    difftime(tail(data$timestamp, n=1), data$timestamp[1], unit='mins')
  )
  pace <- t/d
  list(string=sprintf('%s min/km', format_pace(pace)), numeric=pace)
})

renderValueBox({
  value <- avgPace()
  valueBox(value = value$string, icon = "fa-tachometer-alt", color = 'success')
})
```

### Time {.value-box}

```{r}
runTime <- reactive({
  data <- rangedata() %>% slice(1, n()) %>%  select('timestamp') %>% pull
  diff <- as.numeric(data[2] - data[1], units='secs')
  diff <- round(diff, digits=2) # seconds_to_period carries floating point imprecision
  list(string=tolower(seconds_to_period(diff)), numeric=diff)
})

renderValueBox({
  value <- runTime()
  valueBox(value = value$string, icon = "fa-stopwatch", color = 'success')
})
```

### Distance {.value-box}

```{r}
distance <- reactive({
  data <- rangedata()
  d <- ( tail(data$distance, n=1) - head(data$distance, n=1) ) / 1000
  list(string=sprintf('%0.2f km', d), numeric=d)
})

renderValueBox({
  value = distance()
  valueBox(value = value$string, icon = "fa-road", color = 'success' )
})
```

### Elevation Slope (imprecise) {.value-box}

```{r}
elSlope <- reactive({
  slopes <- rangedata()$slope
  slope <- floor(sum(slopes[slopes>0]))
  list(string=sprintf('%d m', slope), numeric=slope)
})

renderValueBox({
  value <- elSlope()
  valueBox(value = value$string, icon = "fa-signal", color = 'success')
})
```

### Average HR {.value-box}

```{r}
avgHR <- reactive({
  avghr <- round(mean(rangedata()$hr, na.rm=T))
  list(string=sprintf('%d bpm', avghr), numeric=avghr)
})

renderValueBox({
  value <- avgHR()
  valueBox(value = value$string, icon = "fa-heartbeat", color = 'danger')
})
```

Row
---

### Map {data-width=50%}

```{r}

output$map <- renderLeaflet({
  leaflet() %>%
    setView(lat=0, lng=0, zoom=2) %>%
    addProviderTiles(providers$CartoDB.Positron)
})
leafletOutput('map')

# clear, when we load new data
observe({
  req(input$fitFile)
  leafletProxy('map') %>%
    clearMarkers()
})

# set current marker when updated from charts or when we load new data
observe({
  runnerIcon <- makeIcon('res/running-solid.svg',
                         iconWidth='24', iconHeight='24')
  current <- fitdata() %>% slice(input$currentLocation)
  leafletProxy('map') %>%
    clearGroup('current') %>%
    addMarkers(lng=current$long, lat=current$lat, group='current', icon=runnerIcon)
})

# draw current range
observe({
  data <- mapdata() %>% select(c('lat', 'long', 'hrzone')) %>% drop_na() %>% rename(lng='long')

  leafletProxy('map') %>%
    clearMarkers() %>%
    fitBounds(min(data$lng), min(data$lat), max(data$lng), max(data$lat)) %>%
    addCircleMarkers(data=data %>% filter(hrzone==HR()$label[1]), fillColor=HR()$color[1],
                     stroke=FALSE, radius=4, fillOpacity=1)  %>%
    addCircleMarkers(data=data %>% filter(hrzone==HR()$label[2]), fillColor=HR()$color[2],
                     stroke=FALSE, radius=4, fillOpacity=1)  %>%
    addCircleMarkers(data=data %>% filter(hrzone==HR()$label[3]), fillColor=HR()$color[3],
                     stroke=FALSE, radius=4, fillOpacity=1)  %>%
    addCircleMarkers(data=data %>% filter(hrzone==HR()$label[4]), fillColor=HR()$color[4],
                     stroke=FALSE, radius=4, fillOpacity=1)  %>%
    addCircleMarkers(data=data %>% filter(hrzone==HR()$label[5]), fillColor=HR()$color[5],
                     stroke=FALSE, radius=4, fillOpacity=1)
})
```

### Charts {data-width=50%}

```{r}
fillCol(
  flex = c(5,7),
  renderDygraph({
    ZonePlot(fitdata, HR(), hlCallbackJs)
  }),
  renderDygraph({
    RacePlot(fitdata, hlCallbackJs)
  })
)
```

Analysis
========

Row {data-height=200}
--------------------

### Trim Data
```{r}
sliderInput('lowerQuantile', value=0.25, min=0, max=1.0, step=0.05,
            'Drop lower quantile speed data')
checkboxInput('removeOutlier', label='Remove Outlier', value=TRUE)
```

### Model
```{r}
radioButtons('modelType', 'Type of Model', choices = c('Loess/LM', 'SVM Regression (experimental)'), selected = 'Loess/LM')
```

### Features
```{r}
sliderInput('slopes', value=6, min=0, max=9, step=1,
            'log2(window size) of slope backlog')
sliderInput('hrs', value=8, min=0, max=9, step=1,
            'log2(window size) of HR backlog')
```

### SVM Params
```{r}
sliderInput('C', value=3, min=-5, max=15, step=0.5,
            'log2(C) / cost')
sliderInput('sigma', value=-2, min=-15, max=5, step=0.5,
            'log2(sigma) / RBF parameter')
```


Row {.tabset}
---

### HR-Pace Index

```{r}
renderPlot({
  HRPacePlot(analysisdata(), rv$lmfit, rv$mfit, HR(), rv)
})
```

### Gradient Histogram

```{r}
renderPlot({
  data <- analysisdata()

  ggplot(data, aes(abs(slope))) +
    geom_histogram()
})
```

Row {data-height=80}
--------------------

### Lactate Threshold (modified D-Max) {.value-box}

```{r}
hrdp <- reactive({
  data <- analysisdata()
  tryCatch(
    {
    hrdpval <- calc_hrdp(data, rv$lmfit, rv$mfit)
    print(hrdpval)

    rv$hrdpCache <- hrdpval
    list(string=sprintf('%s min/km | %.0f bpm', format_pace(hrdpval$pace), round(hrdpval$hr)),
         numeric=hrdpval$speed)
    }, error = function(cond) {
      rv$hrdpCache <- NULL
      message('Could not calculate hrdp')
      message(cond$message)
      list(string=cond$message)
    })
})

renderValueBox({
  value <- hrdp()
  valueBox(value = value$string, icon = "fa-tachometer-alt", color = 'primary')
})
```