actigraph / agcounts

Code for the technical report on the ActiLife counts algorithm.
GNU General Public License v3.0
29 stars 13 forks source link

Discrepancy from ActiLife Output #19

Open muschellij2 opened 9 months ago

muschellij2 commented 9 months ago

I downloaded some open data from https://springernature.figshare.com/collections/Upper_limb_activity_of_twenty_myoelectric_prosthesis_users_and_twenty_healthy_anatomically_intact_adults_/4457855 in order to ensure the output from @bhelsel and the agcounts R package match up with ActiLife output. That package is using the code from this repo, adapted for R.

The GT3X files are available at https://springernature.figshare.com/articles/dataset/Unprocessed_raw_30Hz_acceleration_data_stored_as_gt3x/7946189?backTo=/collections/Upper_limb_activity_of_twenty_myoelectric_prosthesis_users_and_twenty_healthy_anatomically_intact_adults_/4457855 or raw CSVs: https://springernature.figshare.com/articles/dataset/Unprocessed_raw_30Hz_acceleration_data_stored_as_csv/7946186?backTo=/collections/Upper_limb_activity_of_twenty_myoelectric_prosthesis_users_and_twenty_healthy_anatomically_intact_adults_/4457855

I picked out the odd case and uploaded it to figshare.

Here is a reproducible example demonstrating the difference. I don't know why this occurs, but would be interested in your thoughts on where the discrepancy is.

library(readr) # install.packages("readr")
library(dplyr) # install.packages("dplyr")
#> 
#> Attaching package: 'dplyr'
#> The following objects are masked from 'package:stats':
#> 
#>     filter, lag
#> The following objects are masked from 'package:base':
#> 
#>     intersect, setdiff, setequal, union
library(agcounts) # remotes::install_github("bhelsel/agcounts")
library(SummarizedActigraphy) # remotes::install_github("muschellij2/SummarizedActigraphy")

options(digits.secs = 3)

Read in the data

raw_url = "https://figshare.com/ndownloader/files/42434487"
onesec_url = "https://figshare.com/ndownloader/files/42434481"
read_it = function(url) {
  file = tempfile(fileext = ".csv.gz")
  curl::curl_download(url, file)
  x = SummarizedActigraphy::read_acc_csv(file)
  readr::stop_for_problems(x)
  x
}
data = read_it(raw_url)
data$header
#> $serial
#> [1] "NEO1B41100262"
#> 
#> $start_time
#> [1] "14:00:00"
#> 
#> $start_date
#> [1] "18/07/2017"
#> 
#> $epoch
#> [1] "00:00:00"
#> 
#> $download_time
#> [1] "19:44:43"
#> 
#> $download_date
#> [1] "31/07/2017"
#> 
#> $battery_voltage
#> [1] "4"
#> 
#> $memory_address
#> [1] "0"
#> 
#> $mode
#> [1] "12"
#> 
#> $firmware
#> [1] "3.2.1"
#> 
#> $actilife_version
#> [1] "6.11.9"
#> 
#> $filter
#> [1] "Normal"
#> 
#> $date_format
#> [1] "dmy"
#> 
#> $sample_rate
#> [1] 30

csv = read_it(onesec_url)
#> Warning in extract_acc_header(file): NAs introduced by coercion
#> Warning in parse_acc_header(hdr): NAs introduced by coercion

#> Warning in parse_acc_header(hdr): NAs introduced by coercion
#> Warning in SummarizedActigraphy::read_acc_csv(file): X/Y/Z not in the data, but
#> only_xyz = TRUE, setting only_xyz = FALSE
#> Warning in get_sample_rate(df): Guessing sample_rate from the data
csv$header
#> $serial
#> [1] "NEO1B41100262"
#> 
#> $start_time
#> [1] "14:00:00"
#> 
#> $start_date
#> [1] "18/07/2017"
#> 
#> $epoch
#> [1] "00:00:01"
#> 
#> $download_time
#> [1] "19:44:43"
#> 
#> $download_date
#> [1] "31/07/2017"
#> 
#> $battery_voltage
#> [1] "4"
#> 
#> $memory_address
#> [1] "0"
#> 
#> $mode
#> [1] "61"
#> 
#> $firmware
#> [1] "3.2.1"
#> 
#> $actilife_version
#> [1] "6.11.9"
#> 
#> $filter
#> [1] "Normal"
#> 
#> $date_format
#> [1] "dmy"
#> 
#> $sample_rate
#> [1] NA

Downloading/Reading 1s Counts from AciLife

x = csv$data
x = x %>% 
  mutate(time = lubridate::dmy_hms(paste0(Date, " ", Time)))
x = x %>% 
  select(time, Axis1, Axis2, Axis3)
head(x)
#> # A tibble: 6 × 4
#>   time                    Axis1 Axis2 Axis3
#>   <dttm>                  <dbl> <dbl> <dbl>
#> 1 2017-07-18 14:00:00.000     0     0     0
#> 2 2017-07-18 14:00:01.000     0     0     0
#> 3 2017-07-18 14:00:02.000     0     0     0
#> 4 2017-07-18 14:00:03.000     0     0     0
#> 5 2017-07-18 14:00:04.000   102    89    55
#> 6 2017-07-18 14:00:05.000   131    76    92

y = data
freq = data$freq
y = tibble::as_tibble(y$data)
attr(y, "sample_rate") = freq
head(y)
#> # A tibble: 6 × 4
#>   time                        X     Y      Z
#>   <dttm>                  <dbl> <dbl>  <dbl>
#> 1 2017-07-18 14:00:00.000 0.314 0.968 -0.196
#> 2 2017-07-18 14:00:00.033 0.331 0.971 -0.22 
#> 3 2017-07-18 14:00:00.066 0.299 0.974 -0.229
#> 4 2017-07-18 14:00:00.099 0.296 0.974 -0.208
#> 5 2017-07-18 14:00:00.133 0.328 0.968 -0.196
#> 6 2017-07-18 14:00:00.166 0.317 0.956 -0.211

Calculate AC from Open Source Version

res = calculate_counts(raw = y, epoch = 1)
res = res %>% 
  select(time, Axis1, Axis2, Axis3)
res = tibble::as_tibble(res)
head(res)
#> # A tibble: 6 × 4
#>   time                    Axis1 Axis2 Axis3
#>   <dttm>                  <dbl> <dbl> <dbl>
#> 1 2017-07-18 14:00:00.000     0     0     0
#> 2 2017-07-18 14:00:01.000     0     0     0
#> 3 2017-07-18 14:00:02.000     0     0     0
#> 4 2017-07-18 14:00:03.000     0     0     0
#> 5 2017-07-18 14:00:04.000   102    89    55
#> 6 2017-07-18 14:00:05.000   131    76    92

the_check = all.equal(res, x)
if (!isTRUE(the_check)) {
  print(the_check)
}
#> [1] "Component \"Axis2\": Mean relative difference: 0.5087719"
#> [2] "Component \"Axis3\": Mean relative difference: 0.0754717"

ind = rowSums(res != x) > 0
res[ind,]
#> # A tibble: 3 × 4
#>   time                    Axis1 Axis2 Axis3
#>   <dttm>                  <dbl> <dbl> <dbl>
#> 1 2017-07-23 22:57:30.000     0    42    39
#> 2 2017-07-23 22:57:31.000     4     2    14
#> 3 2017-07-23 22:57:32.000     0    13     2
x[ind,]
#> # A tibble: 3 × 4
#>   time                    Axis1 Axis2 Axis3
#>   <dttm>                  <dbl> <dbl> <dbl>
#> 1 2017-07-23 22:57:30.000     0    17    41
#> 2 2017-07-23 22:57:31.000     4     4    12
#> 3 2017-07-23 22:57:32.000     0    11     2

y %>% 
  filter(between(time, lubridate::as_datetime("2017-07-23 22:57:30"), 
                 lubridate::as_datetime("2017-07-23 22:57:33"))) %>% 
  print(n = nrow(y))
#> # A tibble: 91 × 4
#>    time                         X      Y      Z
#>    <dttm>                   <dbl>  <dbl>  <dbl>
#>  1 2017-07-23 22:57:30.000  0.768  0.384  0.196
#>  2 2017-07-23 22:57:30.033  0.328  0.106  0.073
#>  3 2017-07-23 22:57:30.066  0.455  0.158  0.27 
#>  4 2017-07-23 22:57:30.099  0.311  0.446  0.504
#>  5 2017-07-23 22:57:30.133  1.20   0.839  1.27 
#>  6 2017-07-23 22:57:30.166  1.37   0.93   1.12 
#>  7 2017-07-23 22:57:30.200  0.977  0.674  0.405
#>  8 2017-07-23 22:57:30.233  1.07   0.525  0.39 
#>  9 2017-07-23 22:57:30.266  0      0      0    
#> 10 2017-07-23 22:57:30.299  0.252  0.015  0.205
#> 11 2017-07-23 22:57:30.333  0.428  0.317  0.543
#> 12 2017-07-23 22:57:30.366  1.35   0.76   1.09 
#> 13 2017-07-23 22:57:30.400  1.40   1.02   1.12 
#> 14 2017-07-23 22:57:30.433  1.07   0.771  0.384
#> 15 2017-07-23 22:57:30.466  1.08   0.566  0.457
#> 16 2017-07-23 22:57:30.500  0.161  0.088  0.065
#> 17 2017-07-23 22:57:30.533  0.337  0.097  0.235
#> 18 2017-07-23 22:57:30.566  0.261  0.399  0.698
#> 19 2017-07-23 22:57:30.599  0.903  0.672  0.757
#> 20 2017-07-23 22:57:30.633  1.47   0.918  1.17 
#> 21 2017-07-23 22:57:30.666  0.959  0.827  0.651
#> 22 2017-07-23 22:57:30.700  0.314  0.487  0.343
#> 23 2017-07-23 22:57:30.733  0.633  0.27   0.267
#> 24 2017-07-23 22:57:30.766  0.337  0.106  0.199
#> 25 2017-07-23 22:57:30.799  0.413  0.284  0.34 
#> 26 2017-07-23 22:57:30.833  0.748  0.654  0.76 
#> 27 2017-07-23 22:57:30.866  1.33   0.897  1.15 
#> 28 2017-07-23 22:57:30.900  1.01   0.724  0.572
#> 29 2017-07-23 22:57:30.933  0.856  0.572  0.317
#> 30 2017-07-23 22:57:30.966  0.997  0.581  0.51 
#> 31 2017-07-23 22:57:31.000  0.067  0.018  0.123
#> 32 2017-07-23 22:57:31.033  0.46   0.144  0.334
#> 33 2017-07-23 22:57:31.066  0.516  0.343  0.569
#> 34 2017-07-23 22:57:31.099  1.07   0.672  0.613
#> 35 2017-07-23 22:57:31.133  1.34   1.01   1.12 
#> 36 2017-07-23 22:57:31.166  0.944  0.818  0.677
#> 37 2017-07-23 22:57:31.200  0.575  0.628  0.326
#> 38 2017-07-23 22:57:31.233  1.01   0.44   0.434
#> 39 2017-07-23 22:57:31.266 -0.167 -0.097 -0.021
#> 40 2017-07-23 22:57:31.299  0.548  0.185  0.323
#> 41 2017-07-23 22:57:31.333  0.66   0.405  0.572
#> 42 2017-07-23 22:57:31.366  1.14   0.894  0.853
#> 43 2017-07-23 22:57:31.400  1.30   1.07   1    
#> 44 2017-07-23 22:57:31.433  1.10   0.824  0.543
#> 45 2017-07-23 22:57:31.466  0.739  0.554  0.384
#> 46 2017-07-23 22:57:31.500  0.305  0.232  0.358
#> 47 2017-07-23 22:57:31.533  0.015 -0.056  0.158
#> 48 2017-07-23 22:57:31.566  0.402  0.161  0.364
#> 49 2017-07-23 22:57:31.599  0.856  0.475  0.677
#> 50 2017-07-23 22:57:31.633  1.37   0.988  1.23 
#> 51 2017-07-23 22:57:31.666  1.18   0.991  0.806
#> 52 2017-07-23 22:57:31.700  0.71   0.61   0.323
#> 53 2017-07-23 22:57:31.733  0.903  0.405  0.258
#> 54 2017-07-23 22:57:31.766  0.39   0.123  0.223
#> 55 2017-07-23 22:57:31.799  0.375  0.117  0.408
#> 56 2017-07-23 22:57:31.833  0.457  0.475  0.554
#> 57 2017-07-23 22:57:31.866  1.07   0.809  0.891
#> 58 2017-07-23 22:57:31.900  1.24   0.924  0.985
#> 59 2017-07-23 22:57:31.933  1.11   0.73   0.66 
#> 60 2017-07-23 22:57:31.966  0.525  0.516  0.37 
#> 61 2017-07-23 22:57:32.000  0.557  0.346  0.434
#> 62 2017-07-23 22:57:32.033  0.155  0.07   0.226
#> 63 2017-07-23 22:57:32.066  0.399  0.191  0.293
#> 64 2017-07-23 22:57:32.099  0.46   0.501  0.493
#> 65 2017-07-23 22:57:32.133  1.15   0.789  0.979
#> 66 2017-07-23 22:57:32.166  1.36   0.909  1.08 
#> 67 2017-07-23 22:57:32.200  1.42   0.798  0.789
#> 68 2017-07-23 22:57:32.233  0.625  0.469  0.296
#> 69 2017-07-23 22:57:32.266  0.589  0.287  0.431
#> 70 2017-07-23 22:57:32.299  0.276  0.041  0.208
#> 71 2017-07-23 22:57:32.333  0.396  0.217  0.384
#> 72 2017-07-23 22:57:32.366  0.625  0.525  0.566
#> 73 2017-07-23 22:57:32.400  1.25   0.824  1.15 
#> 74 2017-07-23 22:57:32.433  1.24   0.871  0.865
#> 75 2017-07-23 22:57:32.466  0.809  0.654  0.196
#> 76 2017-07-23 22:57:32.500  0.947  0.501  0.317
#> 77 2017-07-23 22:57:32.533  0.246  0.026  0.147
#> 78 2017-07-23 22:57:32.566  0.416  0.141  0.39 
#> 79 2017-07-23 22:57:32.599  0.496  0.343  0.616
#> 80 2017-07-23 22:57:32.633  1.22   0.654  0.941
#> 81 2017-07-23 22:57:32.666  1.43   1.01   1.07 
#> 82 2017-07-23 22:57:32.700  1.14   0.883  0.534
#> 83 2017-07-23 22:57:32.733  0.839  0.604  0.393
#> 84 2017-07-23 22:57:32.766  0.457  0.24   0.355
#> 85 2017-07-23 22:57:32.799 -0.018 -0.076  0.196
#> 86 2017-07-23 22:57:32.833  0.331  0.238  0.46 
#> 87 2017-07-23 22:57:32.866  1.10   0.452  0.718
#> 88 2017-07-23 22:57:32.900  1.17   0.818  0.93 
#> 89 2017-07-23 22:57:32.933  1.05   0.9    0.657
#> 90 2017-07-23 22:57:32.966  1.23   0.733  0.572
#> 91 2017-07-23 22:57:33.000  0.754  0.51   0.587

Fixing the zeros

Here we “fix” any row with all zeroes with a last observation carried forward, to mimic what we think what happens in idle sleep mode in ActiLife. Is that correct?

y2 = fix_zeros(y)

We see the there are no zero rows

y2 %>% 
  filter(between(time, lubridate::as_datetime("2017-07-23 22:57:30"), 
                 lubridate::as_datetime("2017-07-23 22:57:33")))  %>% 
  print(n = nrow(y2))
#> # A tibble: 91 × 4
#>    time                         X      Y      Z
#>    <dttm>                   <dbl>  <dbl>  <dbl>
#>  1 2017-07-23 22:57:30.000  0.768  0.384  0.196
#>  2 2017-07-23 22:57:30.033  0.328  0.106  0.073
#>  3 2017-07-23 22:57:30.066  0.455  0.158  0.27 
#>  4 2017-07-23 22:57:30.099  0.311  0.446  0.504
#>  5 2017-07-23 22:57:30.133  1.20   0.839  1.27 
#>  6 2017-07-23 22:57:30.166  1.37   0.93   1.12 
#>  7 2017-07-23 22:57:30.200  0.977  0.674  0.405
#>  8 2017-07-23 22:57:30.233  1.07   0.525  0.39 
#>  9 2017-07-23 22:57:30.266  1.07   0.525  0.39 
#> 10 2017-07-23 22:57:30.299  0.252  0.015  0.205
#> 11 2017-07-23 22:57:30.333  0.428  0.317  0.543
#> 12 2017-07-23 22:57:30.366  1.35   0.76   1.09 
#> 13 2017-07-23 22:57:30.400  1.40   1.02   1.12 
#> 14 2017-07-23 22:57:30.433  1.07   0.771  0.384
#> 15 2017-07-23 22:57:30.466  1.08   0.566  0.457
#> 16 2017-07-23 22:57:30.500  0.161  0.088  0.065
#> 17 2017-07-23 22:57:30.533  0.337  0.097  0.235
#> 18 2017-07-23 22:57:30.566  0.261  0.399  0.698
#> 19 2017-07-23 22:57:30.599  0.903  0.672  0.757
#> 20 2017-07-23 22:57:30.633  1.47   0.918  1.17 
#> 21 2017-07-23 22:57:30.666  0.959  0.827  0.651
#> 22 2017-07-23 22:57:30.700  0.314  0.487  0.343
#> 23 2017-07-23 22:57:30.733  0.633  0.27   0.267
#> 24 2017-07-23 22:57:30.766  0.337  0.106  0.199
#> 25 2017-07-23 22:57:30.799  0.413  0.284  0.34 
#> 26 2017-07-23 22:57:30.833  0.748  0.654  0.76 
#> 27 2017-07-23 22:57:30.866  1.33   0.897  1.15 
#> 28 2017-07-23 22:57:30.900  1.01   0.724  0.572
#> 29 2017-07-23 22:57:30.933  0.856  0.572  0.317
#> 30 2017-07-23 22:57:30.966  0.997  0.581  0.51 
#> 31 2017-07-23 22:57:31.000  0.067  0.018  0.123
#> 32 2017-07-23 22:57:31.033  0.46   0.144  0.334
#> 33 2017-07-23 22:57:31.066  0.516  0.343  0.569
#> 34 2017-07-23 22:57:31.099  1.07   0.672  0.613
#> 35 2017-07-23 22:57:31.133  1.34   1.01   1.12 
#> 36 2017-07-23 22:57:31.166  0.944  0.818  0.677
#> 37 2017-07-23 22:57:31.200  0.575  0.628  0.326
#> 38 2017-07-23 22:57:31.233  1.01   0.44   0.434
#> 39 2017-07-23 22:57:31.266 -0.167 -0.097 -0.021
#> 40 2017-07-23 22:57:31.299  0.548  0.185  0.323
#> 41 2017-07-23 22:57:31.333  0.66   0.405  0.572
#> 42 2017-07-23 22:57:31.366  1.14   0.894  0.853
#> 43 2017-07-23 22:57:31.400  1.30   1.07   1    
#> 44 2017-07-23 22:57:31.433  1.10   0.824  0.543
#> 45 2017-07-23 22:57:31.466  0.739  0.554  0.384
#> 46 2017-07-23 22:57:31.500  0.305  0.232  0.358
#> 47 2017-07-23 22:57:31.533  0.015 -0.056  0.158
#> 48 2017-07-23 22:57:31.566  0.402  0.161  0.364
#> 49 2017-07-23 22:57:31.599  0.856  0.475  0.677
#> 50 2017-07-23 22:57:31.633  1.37   0.988  1.23 
#> 51 2017-07-23 22:57:31.666  1.18   0.991  0.806
#> 52 2017-07-23 22:57:31.700  0.71   0.61   0.323
#> 53 2017-07-23 22:57:31.733  0.903  0.405  0.258
#> 54 2017-07-23 22:57:31.766  0.39   0.123  0.223
#> 55 2017-07-23 22:57:31.799  0.375  0.117  0.408
#> 56 2017-07-23 22:57:31.833  0.457  0.475  0.554
#> 57 2017-07-23 22:57:31.866  1.07   0.809  0.891
#> 58 2017-07-23 22:57:31.900  1.24   0.924  0.985
#> 59 2017-07-23 22:57:31.933  1.11   0.73   0.66 
#> 60 2017-07-23 22:57:31.966  0.525  0.516  0.37 
#> 61 2017-07-23 22:57:32.000  0.557  0.346  0.434
#> 62 2017-07-23 22:57:32.033  0.155  0.07   0.226
#> 63 2017-07-23 22:57:32.066  0.399  0.191  0.293
#> 64 2017-07-23 22:57:32.099  0.46   0.501  0.493
#> 65 2017-07-23 22:57:32.133  1.15   0.789  0.979
#> 66 2017-07-23 22:57:32.166  1.36   0.909  1.08 
#> 67 2017-07-23 22:57:32.200  1.42   0.798  0.789
#> 68 2017-07-23 22:57:32.233  0.625  0.469  0.296
#> 69 2017-07-23 22:57:32.266  0.589  0.287  0.431
#> 70 2017-07-23 22:57:32.299  0.276  0.041  0.208
#> 71 2017-07-23 22:57:32.333  0.396  0.217  0.384
#> 72 2017-07-23 22:57:32.366  0.625  0.525  0.566
#> 73 2017-07-23 22:57:32.400  1.25   0.824  1.15 
#> 74 2017-07-23 22:57:32.433  1.24   0.871  0.865
#> 75 2017-07-23 22:57:32.466  0.809  0.654  0.196
#> 76 2017-07-23 22:57:32.500  0.947  0.501  0.317
#> 77 2017-07-23 22:57:32.533  0.246  0.026  0.147
#> 78 2017-07-23 22:57:32.566  0.416  0.141  0.39 
#> 79 2017-07-23 22:57:32.599  0.496  0.343  0.616
#> 80 2017-07-23 22:57:32.633  1.22   0.654  0.941
#> 81 2017-07-23 22:57:32.666  1.43   1.01   1.07 
#> 82 2017-07-23 22:57:32.700  1.14   0.883  0.534
#> 83 2017-07-23 22:57:32.733  0.839  0.604  0.393
#> 84 2017-07-23 22:57:32.766  0.457  0.24   0.355
#> 85 2017-07-23 22:57:32.799 -0.018 -0.076  0.196
#> 86 2017-07-23 22:57:32.833  0.331  0.238  0.46 
#> 87 2017-07-23 22:57:32.866  1.10   0.452  0.718
#> 88 2017-07-23 22:57:32.900  1.17   0.818  0.93 
#> 89 2017-07-23 22:57:32.933  1.05   0.9    0.657
#> 90 2017-07-23 22:57:32.966  1.23   0.733  0.572
#> 91 2017-07-23 22:57:33.000  0.754  0.51   0.587

Calculating AC on LOCF Data

res2 = calculate_counts(raw = y2, epoch = 1)
res2 = res2 %>% 
  select(time, Axis1, Axis2, Axis3)
res2 = tibble::as_tibble(res2)
the_check = all.equal(res2, x)
if (!isTRUE(the_check)) {
  print(the_check)
}
#> [1] "Component \"Axis2\": Mean relative difference: 0.5087719"
#> [2] "Component \"Axis3\": Mean relative difference: 0.0754717"

Still does not match up

ind = rowSums(res2 != x) > 0
res2[ind,]
#> # A tibble: 3 × 4
#>   time                    Axis1 Axis2 Axis3
#>   <dttm>                  <dbl> <dbl> <dbl>
#> 1 2017-07-23 22:57:30.000     0    42    39
#> 2 2017-07-23 22:57:31.000     4     2    14
#> 3 2017-07-23 22:57:32.000     0    13     2
x[ind,]
#> # A tibble: 3 × 4
#>   time                    Axis1 Axis2 Axis3
#>   <dttm>                  <dbl> <dbl> <dbl>
#> 1 2017-07-23 22:57:30.000     0    17    41
#> 2 2017-07-23 22:57:31.000     4     4    12
#> 3 2017-07-23 22:57:32.000     0    11     2
res[ind,]
#> # A tibble: 3 × 4
#>   time                    Axis1 Axis2 Axis3
#>   <dttm>                  <dbl> <dbl> <dbl>
#> 1 2017-07-23 22:57:30.000     0    42    39
#> 2 2017-07-23 22:57:31.000     4     2    14
#> 3 2017-07-23 22:57:32.000     0    13     2

Created on 2023-09-22 with reprex v2.0.2

Session info ``` r sessioninfo::session_info() #> ─ Session info ─────────────────────────────────────────────────────────────── #> setting value #> version R version 4.2.1 (2022-06-23) #> os macOS Big Sur ... 10.16 #> system x86_64, darwin17.0 #> ui X11 #> language (EN) #> collate en_US.UTF-8 #> ctype en_US.UTF-8 #> tz America/New_York #> date 2023-09-22 #> pandoc 3.1.5 @ /usr/local/bin/ (via rmarkdown) #> #> ─ Packages ─────────────────────────────────────────────────────────────────── #> package * version date (UTC) lib source #> agcounts * 0.6.4 2023-07-26 [1] local #> anytime 0.3.9 2020-08-27 [1] CRAN (R 4.2.0) #> bit 4.0.5 2022-11-15 [1] CRAN (R 4.2.0) #> bit64 4.0.5 2020-08-30 [1] CRAN (R 4.2.0) #> blob 1.2.4 2023-03-17 [1] CRAN (R 4.2.0) #> bslib 0.5.0 2023-06-09 [1] CRAN (R 4.2.0) #> cachem 1.0.8 2023-05-01 [1] CRAN (R 4.2.0) #> cli 3.6.1 2023-03-23 [1] CRAN (R 4.2.0) #> colorspace 2.1-0 2023-01-23 [1] CRAN (R 4.2.0) #> crayon 1.5.2 2022-09-29 [1] CRAN (R 4.2.0) #> curl 5.0.1 2023-06-07 [1] CRAN (R 4.2.0) #> data.table 1.14.9 2023-08-21 [1] local #> DBI 1.1.3 2022-06-18 [1] CRAN (R 4.2.0) #> digest 0.6.33 2023-07-07 [1] CRAN (R 4.2.0) #> dplyr * 1.1.2 2023-04-20 [1] CRAN (R 4.2.0) #> ellipsis 0.3.2 2021-04-29 [1] CRAN (R 4.2.0) #> evaluate 0.21 2023-05-05 [1] CRAN (R 4.2.0) #> fansi 1.0.4 2023-01-22 [1] CRAN (R 4.2.0) #> fastmap 1.1.1 2023-02-24 [1] CRAN (R 4.2.0) #> fs 1.6.3 2023-07-20 [1] CRAN (R 4.2.0) #> generics 0.1.3 2022-07-05 [1] CRAN (R 4.2.0) #> GGIR 2.9-0 2023-03-24 [1] CRAN (R 4.2.0) #> ggplot2 3.4.2 2023-04-03 [1] CRAN (R 4.2.0) #> glue 1.6.2 2022-02-24 [1] CRAN (R 4.2.0) #> gsignal 0.3-5 2022-05-15 [1] CRAN (R 4.2.0) #> gtable 0.3.3 2023-03-21 [1] CRAN (R 4.2.0) #> hms 1.1.3 2023-03-21 [1] CRAN (R 4.2.0) #> htmltools 0.5.5 2023-03-23 [1] CRAN (R 4.2.0) #> htmlwidgets 1.6.2 2023-03-17 [1] CRAN (R 4.2.0) #> httpuv 1.6.11 2023-05-11 [1] CRAN (R 4.2.0) #> jquerylib 0.1.4 2021-04-26 [1] CRAN (R 4.2.0) #> jsonlite 1.8.7 2023-06-29 [1] CRAN (R 4.2.0) #> knitr 1.43 2023-05-25 [1] CRAN (R 4.2.0) #> later 1.3.1 2023-05-02 [1] CRAN (R 4.2.0) #> lattice 0.21-8 2023-04-05 [1] CRAN (R 4.2.0) #> lifecycle 1.0.3 2022-10-07 [1] CRAN (R 4.2.0) #> lubridate 1.9.2 2023-02-10 [1] CRAN (R 4.2.0) #> magrittr 2.0.3 2022-03-30 [1] CRAN (R 4.2.0) #> Matrix 1.6-0 2023-07-08 [1] CRAN (R 4.2.1) #> memoise 2.0.1 2021-11-26 [1] CRAN (R 4.2.0) #> mime 0.12 2021-09-28 [1] CRAN (R 4.2.0) #> munsell 0.5.0 2018-06-12 [1] CRAN (R 4.2.0) #> pillar 1.9.0 2023-03-22 [1] CRAN (R 4.2.0) #> pkgconfig 2.0.3 2019-09-22 [1] CRAN (R 4.2.0) #> png 0.1-8 2022-11-29 [1] CRAN (R 4.2.0) #> pracma 2.4.2 2022-09-22 [1] CRAN (R 4.2.0) #> promises 1.2.0.1 2021-02-11 [1] CRAN (R 4.2.0) #> purrr 1.0.1 2023-01-10 [1] CRAN (R 4.2.0) #> R.cache 0.16.0 2022-07-21 [1] CRAN (R 4.2.0) #> R.methodsS3 1.8.2 2022-06-13 [1] CRAN (R 4.2.0) #> R.oo 1.25.0 2022-06-12 [1] CRAN (R 4.2.0) #> R.utils 2.12.2 2022-11-11 [1] CRAN (R 4.2.0) #> R6 2.5.1 2021-08-19 [1] CRAN (R 4.2.0) #> Rcpp 1.0.11 2023-07-06 [1] CRAN (R 4.2.0) #> reactable 0.4.4 2023-03-12 [1] CRAN (R 4.2.0) #> read.gt3x 1.2.0 2022-06-30 [1] CRAN (R 4.2.0) #> readr * 2.1.4 2023-02-10 [1] CRAN (R 4.2.0) #> reprex 2.0.2 2022-08-17 [1] CRAN (R 4.2.0) #> reticulate 1.30 2023-06-09 [1] CRAN (R 4.2.0) #> rlang 1.1.1 2023-04-28 [1] CRAN (R 4.2.0) #> rmarkdown 2.23 2023-07-01 [1] CRAN (R 4.2.0) #> RSQLite 2.3.1 2023-04-03 [1] CRAN (R 4.2.0) #> rstudioapi 0.15.0 2023-07-07 [1] CRAN (R 4.2.0) #> sass 0.4.7 2023-07-15 [1] CRAN (R 4.2.0) #> scales 1.2.1 2022-08-20 [1] CRAN (R 4.2.0) #> sessioninfo 1.2.2 2021-12-06 [1] CRAN (R 4.2.0) #> shiny 1.7.4.1 2023-07-06 [1] CRAN (R 4.2.0) #> stringi 1.7.12 2023-01-11 [1] CRAN (R 4.2.0) #> stringr 1.5.0 2022-12-02 [1] CRAN (R 4.2.0) #> styler 1.10.1 2023-06-05 [1] CRAN (R 4.2.0) #> SummarizedActigraphy * 0.6.0 2023-08-18 [1] Bioconductor #> tibble 3.2.1 2023-03-20 [1] CRAN (R 4.2.0) #> tidyselect 1.2.0 2022-10-10 [1] CRAN (R 4.2.0) #> timechange 0.2.0 2023-01-11 [1] CRAN (R 4.2.0) #> tsibble 1.1.3 2022-10-09 [1] CRAN (R 4.2.0) #> tzdb 0.4.0 2023-05-12 [1] CRAN (R 4.2.0) #> utf8 1.2.3 2023-01-31 [1] CRAN (R 4.2.0) #> vctrs 0.6.3 2023-06-14 [1] CRAN (R 4.2.0) #> vroom 1.6.3 2023-04-28 [1] CRAN (R 4.2.0) #> withr 2.5.0 2022-03-03 [1] CRAN (R 4.2.0) #> xfun 0.39 2023-04-20 [1] CRAN (R 4.2.0) #> xtable 1.8-4 2019-04-21 [1] CRAN (R 4.2.0) #> yaml 2.3.7 2023-01-23 [1] CRAN (R 4.2.0) #> zoo 1.8-12 2023-04-13 [1] CRAN (R 4.2.0) #> #> [1] /Library/Frameworks/R.framework/Versions/4.2/Resources/library #> #> ─ Python configuration ─────────────────────────────────────────────────────── #> python: /Users/johnmuschelli/miniconda3/bin/python3 #> libpython: /Users/johnmuschelli/miniconda3/lib/libpython3.11.dylib #> pythonhome: /Users/johnmuschelli/miniconda3:/Users/johnmuschelli/miniconda3 #> version: 3.11.4 (main, Jul 5 2023, 08:41:25) [Clang 14.0.6 ] #> numpy: /Users/johnmuschelli/miniconda3/lib/python3.11/site-packages/numpy #> numpy_version: 1.25.2 #> pygt3x: /Users/johnmuschelli/miniconda3/lib/python3.11/site-packages/pygt3x #> #> NOTE: Python version was forced by RETICULATE_PYTHON #> #> ────────────────────────────────────────────────────────────────────────────── ```
rouzbeh commented 9 months ago

Hi @muschellij2, is there a reason why you are not reporting this on the R package itself?

The correct way of detecting idle sleep mode is documented here: https://github.com/actigraph/GT3X-File-Format#log-record-types . This way is implemented in our package, pygt3x, but I am not aware of other implementations.

I guess the first step is to compare the data that is read (export to csv from actilife and compare to the input of R's agcounts).

muschellij2 commented 9 months ago

I am not so comfortable with how you perform reproducible examples easily in python but I had this on hand but can make an example demonstrating this result in python. They have also released their gt3x files if that’s what you recommend starting from compared to actigraph CSV output

On Mon, Sep 25, 2023 at 6:22 AM Ali Neishabouri @.***> wrote:

Hi @muschellij2 https://github.com/muschellij2, is there a reason why you are not reporting this on the R package itself?

The correct way of detecting idle sleep mode is documented here: https://github.com/actigraph/GT3X-File-Format#log-record-types . This way is implemented in our package, pygt3x, but I am not aware of other implementations.

I guess the first step is to compare the data that is read (export to csv from actilife and compare to the input of R's agcounts).

— Reply to this email directly, view it on GitHub https://github.com/actigraph/agcounts/issues/19#issuecomment-1733381884, or unsubscribe https://github.com/notifications/unsubscribe-auth/AAIGPLVRLP44TSN47UG4ELLX4FLOXANCNFSM6AAAAAA5DTVVHQ . You are receiving this because you were mentioned.Message ID: @.***>

muschellij2 commented 7 months ago

hey @rouzbeh, I realized the issue has to do with the last observation carried forward due to idle sleep mode. See below. The code is still R code, but the agcounter package (https://github.com/muschellij2/agcounter) simply calls your Python code to get the output. So something is going on there (may be due to filtering). I wonder if we'd get the same or similar results for dropping rows with all 0.

library(readr) # install.packages("readr")
library(dplyr) # install.packages("dplyr")
#> 
#> Attaching package: 'dplyr'
#> The following objects are masked from 'package:stats':
#> 
#>     filter, lag
#> The following objects are masked from 'package:base':
#> 
#>     intersect, setdiff, setequal, union
library(agcounter) # remotes::install_github("bhelsel/agcounts")
library(SummarizedActigraphy) # remotes::install_github("muschellij2/SummarizedActigraphy")

options(digits.secs = 3)

raw_url = "https://figshare.com/ndownloader/files/42434487"
onesec_url = "https://figshare.com/ndownloader/files/42434481"
read_it = function(url) {
  file = tempfile(fileext = ".csv.gz")
  curl::curl_download(url, file)
  x = SummarizedActigraphy::read_acc_csv(file)
  readr::stop_for_problems(x)
  x
}
data = read_it(raw_url)
data$header
#> $serial
#> [1] "NEO1B41100262"
#> 
#> $start_time
#> [1] "14:00:00"
#> 
#> $start_date
#> [1] "18/07/2017"
#> 
#> $epoch
#> [1] "00:00:00"
#> 
#> $download_time
#> [1] "19:44:43"
#> 
#> $download_date
#> [1] "31/07/2017"
#> 
#> $battery_voltage
#> [1] "4"
#> 
#> $memory_address
#> [1] "0"
#> 
#> $mode
#> [1] "12"
#> 
#> $firmware
#> [1] "3.2.1"
#> 
#> $actilife_version
#> [1] "6.11.9"
#> 
#> $filter
#> [1] "Normal"
#> 
#> $date_format
#> [1] "dmy"
#> 
#> $sample_rate
#> [1] 30

csv = read_it(onesec_url)
#> Warning in extract_acc_header(file): NAs introduced by coercion
#> Warning in parse_acc_header(hdr): NAs introduced by coercion

#> Warning in parse_acc_header(hdr): NAs introduced by coercion
#> Warning in SummarizedActigraphy::read_acc_csv(file): X/Y/Z not in the data, but
#> only_xyz = TRUE, setting only_xyz = FALSE
#> Warning in get_sample_rate(df): Guessing sample_rate from the data
csv$header
#> $serial
#> [1] "NEO1B41100262"
#> 
#> $start_time
#> [1] "14:00:00"
#> 
#> $start_date
#> [1] "18/07/2017"
#> 
#> $epoch
#> [1] "00:00:01"
#> 
#> $download_time
#> [1] "19:44:43"
#> 
#> $download_date
#> [1] "31/07/2017"
#> 
#> $battery_voltage
#> [1] "4"
#> 
#> $memory_address
#> [1] "0"
#> 
#> $mode
#> [1] "61"
#> 
#> $firmware
#> [1] "3.2.1"
#> 
#> $actilife_version
#> [1] "6.11.9"
#> 
#> $filter
#> [1] "Normal"
#> 
#> $date_format
#> [1] "dmy"
#> 
#> $sample_rate
#> [1] NA

x = csv$data
x = x %>% 
  mutate(time = lubridate::dmy_hms(paste0(Date, " ", Time)))
x = x %>% 
  select(time, Axis1, Axis2, Axis3)
head(x)
#> # A tibble: 6 × 4
#>   time                    Axis1 Axis2 Axis3
#>   <dttm>                  <dbl> <dbl> <dbl>
#> 1 2017-07-18 14:00:00.000     0     0     0
#> 2 2017-07-18 14:00:01.000     0     0     0
#> 3 2017-07-18 14:00:02.000     0     0     0
#> 4 2017-07-18 14:00:03.000     0     0     0
#> 5 2017-07-18 14:00:04.000   102    89    55
#> 6 2017-07-18 14:00:05.000   131    76    92

y = data
freq = data$freq
y = tibble::as_tibble(y$data)
attr(y, "sample_rate") = freq
head(y)
#> # A tibble: 6 × 4
#>   time                        X     Y      Z
#>   <dttm>                  <dbl> <dbl>  <dbl>
#> 1 2017-07-18 14:00:00.000 0.314 0.968 -0.196
#> 2 2017-07-18 14:00:00.033 0.331 0.971 -0.22 
#> 3 2017-07-18 14:00:00.066 0.299 0.974 -0.229
#> 4 2017-07-18 14:00:00.099 0.296 0.974 -0.208
#> 5 2017-07-18 14:00:00.133 0.328 0.968 -0.196
#> 6 2017-07-18 14:00:00.166 0.317 0.956 -0.211

This is calling the python files directly:

agcounter::extract_counts
#> function (raw, epoch_in_seconds = 1L, sample_rate = NULL, verbose = TRUE) 
#> {
#>     f = get_ag_functions()
#>     epoch_in_seconds = as.integer(epoch_in_seconds)
#>     verbose = as.integer(verbose)
#>     if (verbose > 0) {
#>         message("Resampling Data")
#>     }
#>     raw = f$`_resample`(raw = raw, frequency = sample_rate, verbose = verbose > 
#>         1)
#>     gc()
#>     if (verbose > 0) {
#>         message("Filtering Data")
#>     }
#>     raw = f$`_bpf_filter`(downsample_data = raw, verbose = verbose > 
#>         1)
#>     if (verbose > 0) {
#>         message("Trimming Data")
#>     }
#>     raw = f$`_trim_data`(bpf_data = raw, lfe_select = FALSE, 
#>         verbose = verbose > 1)
#>     if (verbose > 0) {
#>         message("Resampling 10Hz")
#>     }
#>     raw = f$`_resample_10hz`(trim_data = raw, verbose = verbose > 
#>         1)
#>     if (verbose > 0) {
#>         message("Getting counts")
#>     }
#>     raw = f$`_sum_counts`(downsample_10hz = raw, epoch_seconds = epoch_in_seconds, 
#>         verbose = verbose > 1)
#>     raw = t(raw)
#>     raw
#> }
#> <bytecode: 0x7fe01f062478>
#> <environment: namespace:agcounter>

res = get_counts(df = y, epoch_in_seconds = 1L)
#> lpf_data not needed
res = res %>% 
  select(time, Axis1 = Y, Axis2 = X, Axis3 = Z)
res = tibble::as_tibble(res)
head(res)
#> # A tibble: 6 × 4
#>   time                    Axis1 Axis2 Axis3
#>   <dttm>                  <dbl> <dbl> <dbl>
#> 1 2017-07-18 14:00:00.000     0     0     0
#> 2 2017-07-18 14:00:01.000     0     0     0
#> 3 2017-07-18 14:00:02.000     0     0     0
#> 4 2017-07-18 14:00:03.000     0     0     0
#> 5 2017-07-18 14:00:04.000   102    89    55
#> 6 2017-07-18 14:00:05.000   131    76    92

the_check = all.equal(res, x)
if (!isTRUE(the_check)) {
  print(the_check)
}

ind = rowSums(res != x) > 0
res[ind,]
#> # A tibble: 0 × 4
#> # ℹ 4 variables: time <dttm>, Axis1 <dbl>, Axis2 <dbl>, Axis3 <dbl>
x[ind,]
#> # A tibble: 0 × 4
#> # ℹ 4 variables: time <dttm>, Axis1 <dbl>, Axis2 <dbl>, Axis3 <dbl>

Python

res = get_counts_py(df = y, epoch_in_seconds = 1L)
#> lpf_data not needed
res = res %>% 
  select(time, Axis1 = Y, Axis2 = X, Axis3 = Z)
res = tibble::as_tibble(res)
head(res)
#> # A tibble: 6 × 4
#>   time                    Axis1 Axis2 Axis3
#>   <dttm>                  <dbl> <dbl> <dbl>
#> 1 2017-07-18 14:00:00.000     0     0     0
#> 2 2017-07-18 14:00:01.000     0     0     0
#> 3 2017-07-18 14:00:02.000     0     0     0
#> 4 2017-07-18 14:00:03.000     0     0     0
#> 5 2017-07-18 14:00:04.000   102    89    55
#> 6 2017-07-18 14:00:05.000   131    76    92

the_check = all.equal(res, x)
if (!isTRUE(the_check)) {
  print(the_check)
}

ind = rowSums(res != x) > 0
res[ind,]
#> # A tibble: 0 × 4
#> # ℹ 4 variables: time <dttm>, Axis1 <dbl>, Axis2 <dbl>, Axis3 <dbl>
x[ind,]
#> # A tibble: 0 × 4
#> # ℹ 4 variables: time <dttm>, Axis1 <dbl>, Axis2 <dbl>, Axis3 <dbl>

y %>% 
  filter(between(time, lubridate::as_datetime("2017-07-23 22:57:30"), 
                 lubridate::as_datetime("2017-07-23 22:57:33"))) %>% 
  print(n = nrow(y))
#> # A tibble: 91 × 4
#>    time                         X      Y      Z
#>    <dttm>                   <dbl>  <dbl>  <dbl>
#>  1 2017-07-23 22:57:30.000  0.768  0.384  0.196
#>  2 2017-07-23 22:57:30.033  0.328  0.106  0.073
#>  3 2017-07-23 22:57:30.066  0.455  0.158  0.27 
#>  4 2017-07-23 22:57:30.099  0.311  0.446  0.504
#>  5 2017-07-23 22:57:30.133  1.20   0.839  1.27 
#>  6 2017-07-23 22:57:30.166  1.37   0.93   1.12 
#>  7 2017-07-23 22:57:30.200  0.977  0.674  0.405
#>  8 2017-07-23 22:57:30.233  1.07   0.525  0.39 
#>  9 2017-07-23 22:57:30.266  0      0      0    
#> 10 2017-07-23 22:57:30.299  0.252  0.015  0.205
#> 11 2017-07-23 22:57:30.333  0.428  0.317  0.543
#> 12 2017-07-23 22:57:30.366  1.35   0.76   1.09 
#> 13 2017-07-23 22:57:30.400  1.40   1.02   1.12 
#> 14 2017-07-23 22:57:30.433  1.07   0.771  0.384
#> 15 2017-07-23 22:57:30.466  1.08   0.566  0.457
#> 16 2017-07-23 22:57:30.500  0.161  0.088  0.065
#> 17 2017-07-23 22:57:30.533  0.337  0.097  0.235
#> 18 2017-07-23 22:57:30.566  0.261  0.399  0.698
#> 19 2017-07-23 22:57:30.599  0.903  0.672  0.757
#> 20 2017-07-23 22:57:30.633  1.47   0.918  1.17 
#> 21 2017-07-23 22:57:30.666  0.959  0.827  0.651
#> 22 2017-07-23 22:57:30.700  0.314  0.487  0.343
#> 23 2017-07-23 22:57:30.733  0.633  0.27   0.267
#> 24 2017-07-23 22:57:30.766  0.337  0.106  0.199
#> 25 2017-07-23 22:57:30.799  0.413  0.284  0.34 
#> 26 2017-07-23 22:57:30.833  0.748  0.654  0.76 
#> 27 2017-07-23 22:57:30.866  1.33   0.897  1.15 
#> 28 2017-07-23 22:57:30.900  1.01   0.724  0.572
#> 29 2017-07-23 22:57:30.933  0.856  0.572  0.317
#> 30 2017-07-23 22:57:30.966  0.997  0.581  0.51 
#> 31 2017-07-23 22:57:31.000  0.067  0.018  0.123
#> 32 2017-07-23 22:57:31.033  0.46   0.144  0.334
#> 33 2017-07-23 22:57:31.066  0.516  0.343  0.569
#> 34 2017-07-23 22:57:31.099  1.07   0.672  0.613
#> 35 2017-07-23 22:57:31.133  1.34   1.01   1.12 
#> 36 2017-07-23 22:57:31.166  0.944  0.818  0.677
#> 37 2017-07-23 22:57:31.200  0.575  0.628  0.326
#> 38 2017-07-23 22:57:31.233  1.01   0.44   0.434
#> 39 2017-07-23 22:57:31.266 -0.167 -0.097 -0.021
#> 40 2017-07-23 22:57:31.299  0.548  0.185  0.323
#> 41 2017-07-23 22:57:31.333  0.66   0.405  0.572
#> 42 2017-07-23 22:57:31.366  1.14   0.894  0.853
#> 43 2017-07-23 22:57:31.400  1.30   1.07   1    
#> 44 2017-07-23 22:57:31.433  1.10   0.824  0.543
#> 45 2017-07-23 22:57:31.466  0.739  0.554  0.384
#> 46 2017-07-23 22:57:31.500  0.305  0.232  0.358
#> 47 2017-07-23 22:57:31.533  0.015 -0.056  0.158
#> 48 2017-07-23 22:57:31.566  0.402  0.161  0.364
#> 49 2017-07-23 22:57:31.599  0.856  0.475  0.677
#> 50 2017-07-23 22:57:31.633  1.37   0.988  1.23 
#> 51 2017-07-23 22:57:31.666  1.18   0.991  0.806
#> 52 2017-07-23 22:57:31.700  0.71   0.61   0.323
#> 53 2017-07-23 22:57:31.733  0.903  0.405  0.258
#> 54 2017-07-23 22:57:31.766  0.39   0.123  0.223
#> 55 2017-07-23 22:57:31.799  0.375  0.117  0.408
#> 56 2017-07-23 22:57:31.833  0.457  0.475  0.554
#> 57 2017-07-23 22:57:31.866  1.07   0.809  0.891
#> 58 2017-07-23 22:57:31.900  1.24   0.924  0.985
#> 59 2017-07-23 22:57:31.933  1.11   0.73   0.66 
#> 60 2017-07-23 22:57:31.966  0.525  0.516  0.37 
#> 61 2017-07-23 22:57:32.000  0.557  0.346  0.434
#> 62 2017-07-23 22:57:32.033  0.155  0.07   0.226
#> 63 2017-07-23 22:57:32.066  0.399  0.191  0.293
#> 64 2017-07-23 22:57:32.099  0.46   0.501  0.493
#> 65 2017-07-23 22:57:32.133  1.15   0.789  0.979
#> 66 2017-07-23 22:57:32.166  1.36   0.909  1.08 
#> 67 2017-07-23 22:57:32.200  1.42   0.798  0.789
#> 68 2017-07-23 22:57:32.233  0.625  0.469  0.296
#> 69 2017-07-23 22:57:32.266  0.589  0.287  0.431
#> 70 2017-07-23 22:57:32.299  0.276  0.041  0.208
#> 71 2017-07-23 22:57:32.333  0.396  0.217  0.384
#> 72 2017-07-23 22:57:32.366  0.625  0.525  0.566
#> 73 2017-07-23 22:57:32.400  1.25   0.824  1.15 
#> 74 2017-07-23 22:57:32.433  1.24   0.871  0.865
#> 75 2017-07-23 22:57:32.466  0.809  0.654  0.196
#> 76 2017-07-23 22:57:32.500  0.947  0.501  0.317
#> 77 2017-07-23 22:57:32.533  0.246  0.026  0.147
#> 78 2017-07-23 22:57:32.566  0.416  0.141  0.39 
#> 79 2017-07-23 22:57:32.599  0.496  0.343  0.616
#> 80 2017-07-23 22:57:32.633  1.22   0.654  0.941
#> 81 2017-07-23 22:57:32.666  1.43   1.01   1.07 
#> 82 2017-07-23 22:57:32.700  1.14   0.883  0.534
#> 83 2017-07-23 22:57:32.733  0.839  0.604  0.393
#> 84 2017-07-23 22:57:32.766  0.457  0.24   0.355
#> 85 2017-07-23 22:57:32.799 -0.018 -0.076  0.196
#> 86 2017-07-23 22:57:32.833  0.331  0.238  0.46 
#> 87 2017-07-23 22:57:32.866  1.10   0.452  0.718
#> 88 2017-07-23 22:57:32.900  1.17   0.818  0.93 
#> 89 2017-07-23 22:57:32.933  1.05   0.9    0.657
#> 90 2017-07-23 22:57:32.966  1.23   0.733  0.572
#> 91 2017-07-23 22:57:33.000  0.754  0.51   0.587

Errors happen after fixing zeros

y2 = fix_zeros(y)

y2 %>% 
  filter(between(time, lubridate::as_datetime("2017-07-23 22:57:30"), 
                 lubridate::as_datetime("2017-07-23 22:57:33")))  %>% 
  print(n = nrow(y2))
#> # A tibble: 91 × 4
#>    time                         X      Y      Z
#>    <dttm>                   <dbl>  <dbl>  <dbl>
#>  1 2017-07-23 22:57:30.000  0.768  0.384  0.196
#>  2 2017-07-23 22:57:30.033  0.328  0.106  0.073
#>  3 2017-07-23 22:57:30.066  0.455  0.158  0.27 
#>  4 2017-07-23 22:57:30.099  0.311  0.446  0.504
#>  5 2017-07-23 22:57:30.133  1.20   0.839  1.27 
#>  6 2017-07-23 22:57:30.166  1.37   0.93   1.12 
#>  7 2017-07-23 22:57:30.200  0.977  0.674  0.405
#>  8 2017-07-23 22:57:30.233  1.07   0.525  0.39 
#>  9 2017-07-23 22:57:30.266  1.07   0.525  0.39 
#> 10 2017-07-23 22:57:30.299  0.252  0.015  0.205
#> 11 2017-07-23 22:57:30.333  0.428  0.317  0.543
#> 12 2017-07-23 22:57:30.366  1.35   0.76   1.09 
#> 13 2017-07-23 22:57:30.400  1.40   1.02   1.12 
#> 14 2017-07-23 22:57:30.433  1.07   0.771  0.384
#> 15 2017-07-23 22:57:30.466  1.08   0.566  0.457
#> 16 2017-07-23 22:57:30.500  0.161  0.088  0.065
#> 17 2017-07-23 22:57:30.533  0.337  0.097  0.235
#> 18 2017-07-23 22:57:30.566  0.261  0.399  0.698
#> 19 2017-07-23 22:57:30.599  0.903  0.672  0.757
#> 20 2017-07-23 22:57:30.633  1.47   0.918  1.17 
#> 21 2017-07-23 22:57:30.666  0.959  0.827  0.651
#> 22 2017-07-23 22:57:30.700  0.314  0.487  0.343
#> 23 2017-07-23 22:57:30.733  0.633  0.27   0.267
#> 24 2017-07-23 22:57:30.766  0.337  0.106  0.199
#> 25 2017-07-23 22:57:30.799  0.413  0.284  0.34 
#> 26 2017-07-23 22:57:30.833  0.748  0.654  0.76 
#> 27 2017-07-23 22:57:30.866  1.33   0.897  1.15 
#> 28 2017-07-23 22:57:30.900  1.01   0.724  0.572
#> 29 2017-07-23 22:57:30.933  0.856  0.572  0.317
#> 30 2017-07-23 22:57:30.966  0.997  0.581  0.51 
#> 31 2017-07-23 22:57:31.000  0.067  0.018  0.123
#> 32 2017-07-23 22:57:31.033  0.46   0.144  0.334
#> 33 2017-07-23 22:57:31.066  0.516  0.343  0.569
#> 34 2017-07-23 22:57:31.099  1.07   0.672  0.613
#> 35 2017-07-23 22:57:31.133  1.34   1.01   1.12 
#> 36 2017-07-23 22:57:31.166  0.944  0.818  0.677
#> 37 2017-07-23 22:57:31.200  0.575  0.628  0.326
#> 38 2017-07-23 22:57:31.233  1.01   0.44   0.434
#> 39 2017-07-23 22:57:31.266 -0.167 -0.097 -0.021
#> 40 2017-07-23 22:57:31.299  0.548  0.185  0.323
#> 41 2017-07-23 22:57:31.333  0.66   0.405  0.572
#> 42 2017-07-23 22:57:31.366  1.14   0.894  0.853
#> 43 2017-07-23 22:57:31.400  1.30   1.07   1    
#> 44 2017-07-23 22:57:31.433  1.10   0.824  0.543
#> 45 2017-07-23 22:57:31.466  0.739  0.554  0.384
#> 46 2017-07-23 22:57:31.500  0.305  0.232  0.358
#> 47 2017-07-23 22:57:31.533  0.015 -0.056  0.158
#> 48 2017-07-23 22:57:31.566  0.402  0.161  0.364
#> 49 2017-07-23 22:57:31.599  0.856  0.475  0.677
#> 50 2017-07-23 22:57:31.633  1.37   0.988  1.23 
#> 51 2017-07-23 22:57:31.666  1.18   0.991  0.806
#> 52 2017-07-23 22:57:31.700  0.71   0.61   0.323
#> 53 2017-07-23 22:57:31.733  0.903  0.405  0.258
#> 54 2017-07-23 22:57:31.766  0.39   0.123  0.223
#> 55 2017-07-23 22:57:31.799  0.375  0.117  0.408
#> 56 2017-07-23 22:57:31.833  0.457  0.475  0.554
#> 57 2017-07-23 22:57:31.866  1.07   0.809  0.891
#> 58 2017-07-23 22:57:31.900  1.24   0.924  0.985
#> 59 2017-07-23 22:57:31.933  1.11   0.73   0.66 
#> 60 2017-07-23 22:57:31.966  0.525  0.516  0.37 
#> 61 2017-07-23 22:57:32.000  0.557  0.346  0.434
#> 62 2017-07-23 22:57:32.033  0.155  0.07   0.226
#> 63 2017-07-23 22:57:32.066  0.399  0.191  0.293
#> 64 2017-07-23 22:57:32.099  0.46   0.501  0.493
#> 65 2017-07-23 22:57:32.133  1.15   0.789  0.979
#> 66 2017-07-23 22:57:32.166  1.36   0.909  1.08 
#> 67 2017-07-23 22:57:32.200  1.42   0.798  0.789
#> 68 2017-07-23 22:57:32.233  0.625  0.469  0.296
#> 69 2017-07-23 22:57:32.266  0.589  0.287  0.431
#> 70 2017-07-23 22:57:32.299  0.276  0.041  0.208
#> 71 2017-07-23 22:57:32.333  0.396  0.217  0.384
#> 72 2017-07-23 22:57:32.366  0.625  0.525  0.566
#> 73 2017-07-23 22:57:32.400  1.25   0.824  1.15 
#> 74 2017-07-23 22:57:32.433  1.24   0.871  0.865
#> 75 2017-07-23 22:57:32.466  0.809  0.654  0.196
#> 76 2017-07-23 22:57:32.500  0.947  0.501  0.317
#> 77 2017-07-23 22:57:32.533  0.246  0.026  0.147
#> 78 2017-07-23 22:57:32.566  0.416  0.141  0.39 
#> 79 2017-07-23 22:57:32.599  0.496  0.343  0.616
#> 80 2017-07-23 22:57:32.633  1.22   0.654  0.941
#> 81 2017-07-23 22:57:32.666  1.43   1.01   1.07 
#> 82 2017-07-23 22:57:32.700  1.14   0.883  0.534
#> 83 2017-07-23 22:57:32.733  0.839  0.604  0.393
#> 84 2017-07-23 22:57:32.766  0.457  0.24   0.355
#> 85 2017-07-23 22:57:32.799 -0.018 -0.076  0.196
#> 86 2017-07-23 22:57:32.833  0.331  0.238  0.46 
#> 87 2017-07-23 22:57:32.866  1.10   0.452  0.718
#> 88 2017-07-23 22:57:32.900  1.17   0.818  0.93 
#> 89 2017-07-23 22:57:32.933  1.05   0.9    0.657
#> 90 2017-07-23 22:57:32.966  1.23   0.733  0.572
#> 91 2017-07-23 22:57:33.000  0.754  0.51   0.587

res2 = get_counts(df = y2, epoch_in_seconds = 1L)
#> lpf_data not needed
res2 = res2 %>% 
  select(time, Axis1 = Y, Axis2 = X, Axis3 = Z)
res2 = tibble::as_tibble(res2)
the_check = all.equal(res2, x)
if (!isTRUE(the_check)) {
  print(the_check)
}
#> [1] "Component \"Axis2\": Mean relative difference: 0.5087719"
#> [2] "Component \"Axis3\": Mean relative difference: 0.0754717"

ind = rowSums(res2 != x) > 0
res2[ind,]
#> # A tibble: 3 × 4
#>   time                    Axis1 Axis2 Axis3
#>   <dttm>                  <dbl> <dbl> <dbl>
#> 1 2017-07-23 22:57:30.000     0    42    39
#> 2 2017-07-23 22:57:31.000     4     2    14
#> 3 2017-07-23 22:57:32.000     0    13     2
x[ind,]
#> # A tibble: 3 × 4
#>   time                    Axis1 Axis2 Axis3
#>   <dttm>                  <dbl> <dbl> <dbl>
#> 1 2017-07-23 22:57:30.000     0    17    41
#> 2 2017-07-23 22:57:31.000     4     4    12
#> 3 2017-07-23 22:57:32.000     0    11     2
res[ind,]
#> # A tibble: 3 × 4
#>   time                    Axis1 Axis2 Axis3
#>   <dttm>                  <dbl> <dbl> <dbl>
#> 1 2017-07-23 22:57:30.000     0    17    41
#> 2 2017-07-23 22:57:31.000     4     4    12
#> 3 2017-07-23 22:57:32.000     0    11     2

Created on 2023-11-28 with reprex v2.0.2

Session info ``` r sessioninfo::session_info() #> ─ Session info ─────────────────────────────────────────────────────────────── #> setting value #> version R version 4.3.1 (2023-06-16) #> os macOS Monterey 12.6 #> system x86_64, darwin20 #> ui X11 #> language (EN) #> collate en_US.UTF-8 #> ctype en_US.UTF-8 #> tz America/New_York #> date 2023-11-28 #> pandoc 3.1.5 @ /usr/local/bin/ (via rmarkdown) #> #> ─ Packages ─────────────────────────────────────────────────────────────────── #> package * version date (UTC) lib source #> agcounter * 0.2.0.9000 2023-10-25 [1] Github (muschellij2/agcounter@c747987) #> anytime 0.3.9 2020-08-27 [1] CRAN (R 4.3.0) #> bit 4.0.5 2022-11-15 [1] CRAN (R 4.3.0) #> bit64 4.0.5 2020-08-30 [1] CRAN (R 4.3.0) #> cli 3.6.1 2023-03-23 [1] CRAN (R 4.3.0) #> crayon 1.5.2 2022-09-29 [1] CRAN (R 4.3.0) #> curl 5.1.0 2023-10-02 [1] CRAN (R 4.3.0) #> digest 0.6.33 2023-07-07 [1] CRAN (R 4.3.0) #> dplyr * 1.1.3 2023-09-03 [1] CRAN (R 4.3.0) #> ellipsis 0.3.2 2021-04-29 [1] CRAN (R 4.3.0) #> evaluate 0.22 2023-09-29 [1] CRAN (R 4.3.0) #> fansi 1.0.5 2023-10-08 [1] CRAN (R 4.3.0) #> fastmap 1.1.1 2023-02-24 [1] CRAN (R 4.3.0) #> fs 1.6.3 2023-07-20 [1] CRAN (R 4.3.0) #> generics 0.1.3 2022-07-05 [1] CRAN (R 4.3.0) #> glue 1.6.2 2022-02-24 [1] CRAN (R 4.3.0) #> hms 1.1.3 2023-03-21 [1] CRAN (R 4.3.0) #> htmltools 0.5.6.1 2023-10-06 [1] CRAN (R 4.3.0) #> jsonlite 1.8.7 2023-06-29 [1] CRAN (R 4.3.0) #> knitr 1.44 2023-09-11 [1] CRAN (R 4.3.0) #> lattice 0.22-5 2023-10-24 [1] CRAN (R 4.3.0) #> lifecycle 1.0.3 2022-10-07 [1] CRAN (R 4.3.0) #> lubridate 1.9.3 2023-09-27 [1] CRAN (R 4.3.0) #> magrittr 2.0.3 2022-03-30 [1] CRAN (R 4.3.0) #> Matrix 1.6-1.1 2023-09-18 [1] CRAN (R 4.3.0) #> pillar 1.9.0 2023-03-22 [1] CRAN (R 4.3.0) #> pkgconfig 2.0.3 2019-09-22 [1] CRAN (R 4.3.0) #> png 0.1-8 2022-11-29 [1] CRAN (R 4.3.0) #> purrr 1.0.2 2023-08-10 [1] CRAN (R 4.3.0) #> R.cache 0.16.0 2022-07-21 [1] CRAN (R 4.3.0) #> R.methodsS3 1.8.2 2022-06-13 [1] CRAN (R 4.3.0) #> R.oo 1.25.0 2022-06-12 [1] CRAN (R 4.3.0) #> R.utils 2.12.2 2022-11-11 [1] CRAN (R 4.3.0) #> R6 2.5.1 2021-08-19 [1] CRAN (R 4.3.0) #> Rcpp 1.0.11 2023-07-06 [1] CRAN (R 4.3.0) #> readr * 2.1.4 2023-02-10 [1] CRAN (R 4.3.0) #> reprex 2.0.2 2022-08-17 [1] CRAN (R 4.3.0) #> reticulate 1.34.0 2023-10-12 [1] CRAN (R 4.3.0) #> rlang 1.1.1 2023-04-28 [1] CRAN (R 4.3.0) #> rmarkdown 2.25 2023-09-18 [1] CRAN (R 4.3.0) #> rstudioapi 0.15.0 2023-07-07 [1] CRAN (R 4.3.0) #> sessioninfo 1.2.2 2021-12-06 [1] CRAN (R 4.3.0) #> styler 1.10.2 2023-08-29 [1] CRAN (R 4.3.0) #> SummarizedActigraphy * 0.6.1 2023-10-18 [1] Bioconductor #> tibble 3.2.1 2023-03-20 [1] CRAN (R 4.3.0) #> tidyselect 1.2.0 2022-10-10 [1] CRAN (R 4.3.0) #> timechange 0.2.0 2023-01-11 [1] CRAN (R 4.3.0) #> tsibble 1.1.3 2022-10-09 [1] CRAN (R 4.3.0) #> tzdb 0.4.0 2023-05-12 [1] CRAN (R 4.3.0) #> utf8 1.2.4 2023-10-22 [1] CRAN (R 4.3.0) #> vctrs 0.6.4 2023-10-12 [1] CRAN (R 4.3.0) #> vroom 1.6.4 2023-10-02 [1] CRAN (R 4.3.0) #> withr 2.5.1 2023-09-26 [1] CRAN (R 4.3.0) #> xfun 0.40 2023-08-09 [1] CRAN (R 4.3.0) #> yaml 2.3.7 2023-01-23 [1] CRAN (R 4.3.0) #> zoo 1.8-12 2023-04-13 [1] CRAN (R 4.3.0) #> #> [1] /Library/Frameworks/R.framework/Versions/4.3-x86_64/Resources/library #> #> ─ Python configuration ─────────────────────────────────────────────────────── #> python: /Users/johnmuschelli/miniconda3/bin/python3 #> libpython: /Users/johnmuschelli/miniconda3/lib/libpython3.11.dylib #> pythonhome: /Users/johnmuschelli/miniconda3:/Users/johnmuschelli/miniconda3 #> version: 3.11.4 (main, Jul 5 2023, 08:41:25) [Clang 14.0.6 ] #> numpy: /Users/johnmuschelli/miniconda3/lib/python3.11/site-packages/numpy #> numpy_version: 1.25.2 #> numpy: /Users/johnmuschelli/miniconda3/lib/python3.11/site-packages/numpy #> #> NOTE: Python version was forced by RETICULATE_PYTHON #> #> ────────────────────────────────────────────────────────────────────────────── ```
rouzbeh commented 6 months ago

So one of the packages is removing 0s, and the other isn't? That does indeed change the results.

I am not sure what the expected behavior is here. Removing rows of 0s as a proxy for idle sleep mode is widespread, but note that we can actually correctly remove idle sleep mode data by reading the corresponding flag from the GT3x file. Our python package pygt3x does this correctly.

muschellij2 commented 6 months ago

In some instances, people remove 0s. In other instances, and according to Actigraph (https://actigraphcorp.my.site.com/support/s/article/Idle-Sleep-Mode-Explained):

While in sleep mode, the last sampled accelerometer value is written into memory at the device's preset sample rate.

Which I believe is incorrect, but some software assumes (given this is Actigraph's publication) that the zeroes should be last observation carried forward (LOCF) to replicate the last sample value.

What is the behavior of pygt3x? I tried it with the first file and it failed as https://github.com/actigraph/pygt3x/issues/37

muschellij2 commented 6 months ago

As noted in https://github.com/actigraph/pygt3x/issues/38, how is the gap between seconds accounted for in the filtering? In the case posed in that issue, the sample rate is 30Hz, so https://github.com/actigraph/agcounts/blob/main/agcounts/extract.py#L51 would not do any resampling.

If resampling was done, I'm not sure if upsample_factor would pull the correct indices:

https://github.com/actigraph/agcounts/blob/fdf97c68536b1d255126ec8cadef5d7f347d022b/agcounts/extract.py#L65

It most likely would since the reading happens at the start of the second in pygt3x, but if there is any shift with gaps in time, this could be pulling wrong elements for resampling.

It would also seem that lfilter does not use the timestamp in the filtering, so I'm unsure how this would accounted for with gaps in time (as opposed to LOCF).

https://github.com/actigraph/agcounts/blob/fdf97c68536b1d255126ec8cadef5d7f347d022b/agcounts/extract.py#L123

rouzbeh commented 3 months ago

Hi @muschellij2 , have you managed to retry with the latest Pygt3x (that fixes the bugs you had found?).

This package (agcounts) does indeed assume there are no gaps in the data.

muschellij2 commented 2 months ago

I haven't been able to. I'm curious, do you have unit tests that check this? The data I linked above should be open to use for things such as these.

rouzbeh commented 2 months ago

We have a bunch of unit-tests (https://github.com/actigraph/agcounts/tree/main/tests), but dealing with gaps is out of the scope of agcounts. The behaviour of idle sleep mode when dealing with GT3x files is in the scope of pygt3x.

agcounts doesn't do LOCF or filling with 0. The function's interface is set so that we don't even know if there are gaps. It is the responsibility of the caller to deal with gaps before calling agcounts.