Exploring turbidity fluctuations

Source: vignettes/turbidity.Rmd

This is an example workflow of how to perform basic turbidity analysis with loadflux package and tidyverse ecosystem.

Hydrological Events

First, we need to demarcate hydrological events

data(djanturb)
df <- hydro_events(dataframe = djanturb,
                   q = discharge,
                   datetime = time,
                   window = 21)

df %>% 
  event_plot(q = discharge,
             datetime = time,
             he = he,
             ssc = ntu,
             y2label = "Turbidity")

Turbidity Index

Then we can calculate Turbidity Index TI for every hydrological event


TI_index <- df %>% 
  group_by(he) %>% 
  nest() %>% 
  mutate(TI = map_dbl(data, ~TI(.x, ntu, time))) %>% 
  select(-data) %>% 
  ungroup()

TI_index
#> # A tibble: 6 × 2
#>      he    TI
#>   <dbl> <dbl>
#> 1     1 0.193
#> 2     2 0.113
#> 3     3 0.147
#> 4     4 0.153
#> 5     5 0.108
#> 6     6 0.137

Hydrological events parameters

To summarize the hydrological events parameters an approach from features package can be used. For this purpose we need to transform our dataframe into tsibble object:

library(tsibble)

df_ts <- df %>% 
  as_tsibble(key = he,
             index = time)

df_ts
#> # A tsibble: 620 x 4 [10m] <Europe/Moscow>
#> # Key:       he [6]
#>       he time                discharge   ntu
#>    <dbl> <dttm>                  <dbl> <dbl>
#>  1     1 2016-06-16 09:30:00      1.23  751.
#>  2     1 2016-06-16 09:40:00      1.23  622.
#>  3     1 2016-06-16 10:00:00      1.23  964 
#>  4     1 2016-06-16 10:10:00      1.23  804.
#>  5     1 2016-06-16 10:30:00      1.25  927.
#>  6     1 2016-06-16 10:40:00      1.26  673.
#>  7     1 2016-06-16 11:00:00      1.27  631.
#>  8     1 2016-06-16 11:10:00      1.28  645.
#>  9     1 2016-06-16 11:30:00      1.28  665.
#> 10     1 2016-06-16 11:40:00      1.29  726.
#> # … with 610 more rows

Then we can calculate start, end and length of the every hydrological event:

library(feasts)
#> Loading required package: fabletools

df_ts %>% 
  features(time,
           feat_event)
#> # A tibble: 6 × 4
#>      he start               end                 length        
#>   <dbl> <dttm>              <dttm>              <drtn>        
#> 1     1 2016-06-16 09:30:00 2016-06-17 10:50:00 25.33333 hours
#> 2     2 2016-06-17 11:00:00 2016-06-18 09:50:00 22.83333 hours
#> 3     3 2016-06-18 10:00:00 2016-06-19 12:20:00 26.33333 hours
#> 4     4 2016-06-19 12:30:00 2016-06-20 11:00:00 22.50000 hours
#> 5     5 2016-06-20 12:00:00 2016-06-21 11:20:00 23.33333 hours
#> 6     6 2016-06-21 11:40:00 2016-06-22 10:20:00 22.66667 hours

Or with the help of brolgar and feasts packages we can calculate turbidity statistics, autocorrelation and spectral functions:

library(brolgar)
library(feasts)

df_ts %>% 
  features(ntu, feat_five_num)
#> # A tibble: 6 × 6
#>      he   min   q25   med   q75   max
#>   <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
#> 1     1 133.  193.   348.  493. 1018.
#> 2     2 126.  222.   306.  439. 1209.
#> 3     3  83.4 118.   154.  226.  520.
#> 4     4  62.4  86.2  114.  196.  451.
#> 5     5  51.8  63.3   76   106.  374.
#> 6     6  81.3 121.   164   262. 1364.

df_ts %>% 
  features(ntu, c(feat_spectral,
                  feat_acf))
#> # A tibble: 6 × 9
#>      he spectral_entropy  acf1 acf10 diff1_acf1 diff1_…¹ diff2…² diff2…³ seaso…⁴
#>   <dbl>            <dbl> <dbl> <dbl>      <dbl>    <dbl>   <dbl>   <dbl>   <dbl>
#> 1     1            0.477 0.735  3.42     -0.433   0.382   -0.598   0.645   0.495
#> 2     2            0.787 0.782  1.90     -0.127   0.0851  -0.533   0.359   0.331
#> 3     3            0.792 0.656  2.70     -0.436   0.217   -0.648   0.467   0.522
#> 4     4            0.517 0.769  4.26     -0.268   0.272   -0.477   0.371   0.637
#> 5     5            0.846 0.618  1.44     -0.201   0.107   -0.455   0.238   0.361
#> 6     6            0.647 0.875  4.24     -0.221   0.132   -0.600   0.425   0.645
#> # … with abbreviated variable names ¹​diff1_acf10, ²​diff2_acf1, ³​diff2_acf10,
#> #   ⁴​season_acf1