Merge branch 'master' into issue708_centralise_cosinor

DESCRIPTION

@@ -1,8 +1,8 @@
 Package: GGIR
 Type: Package
 Title: Raw Accelerometer Data Analysis
-Version: 2.10-3
-Date: 2023-09-30
+Version: 2.10-4
+Date: 2023-10-05
 Authors@R: c(person("Vincent T","van Hees",role=c("aut","cre"),
                     email="[email protected]"),
              person("Jairo H","Migueles",role="aut",

R/GGIR.R

@@ -68,13 +68,18 @@ GGIR = function(mode = 1:5, datadir = c(), outputdir = c(),
   if (is.null(f0) || f0 < 1) { # What file to start with?
     f0 = 1
   }
-  if (is.null(f1) || f1 < 1) {  # What file to end with?
+  if ((is.null(f1) || f1 < 1) && 1 %in% mode) {  # What file to end with?
+    # Do not modify f1 here when not attempting to process GGIR part 1
+    # we only expect datadir to exist when running part 1
     if (filelist == FALSE) {
       f1 <- length(dir(datadir, recursive = TRUE, ignore.case = TRUE, pattern = "[.](csv|bin|Rda|wa|cw|gt3)")) # modified by JH
     } else {
       f1 = length(datadir) #modified
     }
   }
+  if (is.null(f1)) {
+    f1 = 0
+  }
   # Establish which parts need to be processed:
   dopart1 = dopart2 = dopart3 = dopart4 = dopart5 = FALSE
   if (length(which(mode == 0)) > 0) {
@@ -176,36 +181,37 @@ GGIR = function(mode = 1:5, datadir = c(), outputdir = c(),
       stop("If you want to derive circadian rhythm indicators, please install package: ActCR.", call. = FALSE)
     }
   }
-
-  checkFormat = TRUE
-  if (all(dir.exists(datadir)) == TRUE) {
-    rawaccfiles = dir(datadir, full.names = TRUE)[f0:f1]
-  } else if (all(file.exists(datadir))) {
-    rawaccfiles = datadir[f0:f1]
-  } else {
-    checkFormat = FALSE
-  }
-
-  if (checkFormat == TRUE) {
-    is_GGIRread_installed = is.element('GGIRread', installed.packages()[,1])
-    is_read.gt3x_installed = is.element('read.gt3x', installed.packages()[,1])
-    # skip this check if GGIRread and read.gt3x are both available
-    if (is_GGIRread_installed == FALSE | is_read.gt3x_installed == FALSE) {
-      getExt = function(x) {
-        tmp = unlist(strsplit(x, "[.]"))
-        return(tmp[length(tmp)])
-      }
-      rawaccfiles_formats = unique(unlist(lapply(rawaccfiles, FUN = getExt)))
-      # axivity (cwa, wav), geneactive (bin), genea (bin):
-      if (any(grepl("cwa|wav|bin", rawaccfiles_formats))) {
-        if (is_GGIRread_installed == FALSE) {
-          stop("If you are working with axivity, geneactiv, or genea files, please install package: GGIRread.", call. = FALSE)
+  if (1 %in% mode) {
+    checkFormat = TRUE
+    if (all(dir.exists(datadir)) == TRUE) {
+      rawaccfiles = dir(datadir, full.names = TRUE)[f0:f1]
+    } else if (all(file.exists(datadir))) {
+      rawaccfiles = datadir[f0:f1]
+    } else {
+      checkFormat = FALSE
+    }
+
+    if (checkFormat == TRUE) {
+      is_GGIRread_installed = is.element('GGIRread', installed.packages()[,1])
+      is_read.gt3x_installed = is.element('read.gt3x', installed.packages()[,1])
+      # skip this check if GGIRread and read.gt3x are both available
+      if (is_GGIRread_installed == FALSE | is_read.gt3x_installed == FALSE) {
+        getExt = function(x) {
+          tmp = unlist(strsplit(x, "[.]"))
+          return(tmp[length(tmp)])
         }
-      }
-      # actigraph (gt3x)
-      if (any(grepl("gt3x", rawaccfiles_formats))) {
-        if (is_read.gt3x_installed == FALSE) {
-          stop(paste0("If you are working with actigraph files, please install package: read.gt3x.", call. = FALSE))
+        rawaccfiles_formats = unique(unlist(lapply(rawaccfiles, FUN = getExt)))
+        # axivity (cwa, wav), geneactive (bin), genea (bin):
+        if (any(grepl("cwa|wav|bin", rawaccfiles_formats))) {
+          if (is_GGIRread_installed == FALSE) {
+            stop("If you are working with axivity, geneactiv, or genea files, please install package: GGIRread.", call. = FALSE)
+          }
+        }
+        # actigraph (gt3x)
+        if (any(grepl("gt3x", rawaccfiles_formats))) {
+          if (is_read.gt3x_installed == FALSE) {
+            stop(paste0("If you are working with actigraph files, please install package: read.gt3x.", call. = FALSE))
+          }
         }
       }
     }
@@ -359,7 +365,9 @@ GGIR = function(mode = 1:5, datadir = c(), outputdir = c(),
       # if (N.files.ms2.out < f1) f1 = N.files.ms2.out
       if (length(f0) == 0) f0 = 1
       if (f1 == 0) f1 = N.files.ms2.out
-      if (length(params_247[["qwindow"]]) > 2 | is.character(params_247[["qwindow"]])) {
+      if (length(params_247[["qwindow"]]) > 2 |
+          is.character(params_247[["qwindow"]]) |
+        (length(params_247[["qwindow"]]) == 2 & !all(c(0, 24) %in% params_247[["qwindow"]]))) {
         store.long = TRUE
       } else {
         store.long = FALSE

R/g.getM5L5.R

@@ -5,10 +5,11 @@ g.getM5L5 = function(varnum, ws3, t0_LFMF, t1_LFMF, M5L5res, winhr,
   # placed in g.analyse because it is not related to M5L5 analyse (23-7-2019)
   meanVarnum = mean(varnum)
   do.M5L5 = meanVarnum > 0 # & length(varnum) > 1440*(60/ws3)
-  if (length(do.M5L5) == 0 | is.na(do.M5L5) == TRUE) do.M5L5 = FALSE
+  nwindow_f = (t1_LFMF-winhr) - t0_LFMF #number of windows for L5M5 analyses
+  if (length(do.M5L5) == 0 | is.na(do.M5L5) == TRUE | nwindow_f < 1) do.M5L5 = FALSE
+
   if (do.M5L5 == TRUE) { # only do the analysis if varnum has values other than zero
     reso = M5L5res #resolution
-    nwindow_f = (t1_LFMF-winhr) - t0_LFMF #number of windows for L5M5 analyses
     nwindow_f = nwindow_f * (60/reso)
     DAYrunav5 = matrix(NA,nwindow_f,1)
     first_hri = (t0_LFMF*(60/reso))

R/g.getstarttime.R

@@ -1,32 +1,10 @@
 g.getstarttime = function(datafile, P, header, mon, dformat, desiredtz, configtz = NULL) {
   #get input variables (relevant when read.myacc.csv is used)
   #------------------------------------------------------------
-  if (mon == MONITOR$GENEA && dformat == FORMAT$BIN) {
-    starttime = P$timestamps2[1]
-  } else if (mon == MONITOR$AXIVITY && dformat == FORMAT$CWA) {
+  if (mon == MONITOR$AXIVITY && dformat == FORMAT$CWA) {
     starttime = P$data[1,1]
     starttime = as.POSIXlt(starttime, tz = desiredtz, origin = "1970-01-01")
     starttime = POSIXtime2iso8601(starttime, tz = desiredtz)
-  } else if (dformat == FORMAT$WAV) {
-    starttime = c()
-    #It seems that Axivity does not store timestamp in a consistent position
-    # therefore, we need to search for it in the data:
-    starttime = format(header[which(rownames(header) == "ICMTzTime"),1])
-    rn = rownames(header)
-    vl = header$value
-    if (length(starttime) == 0) {
-      if (length(which(rn == "Start")) > 0) {
-        starttime = format(header$value[which(rn == "Start")])
-        #in one of the files starttime is hidden in rowname
-        if (length(starttime) == 0) starttime = rownames(header)[2]
-      }
-      #in one of the files start variable name is hidden in the values
-      if (length(which(vl == "Start")) > 0) {
-        starttime = header$value[2]
-      }
-    }
-    if (length(starttime) == 0) starttime = P$timestamp # initially used, but apparently its is corrupted sometimes, so I am now using ICMTzTime
-    if (length(P$timestamp) == 0) starttime = format(P$hvalues[which(P$hnames == "Start")])
   } else if (mon == MONITOR$GENEACTIV && dformat == FORMAT$BIN) {
     starttime = as.POSIXlt(P$data.out$time[1], tz = desiredtz, origin = "1970-01-01")
     starttime = POSIXtime2iso8601(starttime, tz = desiredtz)

R/get_starttime_weekday_meantemp_truncdata.R

@@ -10,16 +10,16 @@ get_starttime_weekday_meantemp_truncdata = function(temp.available, monc, dforma
     use.temp = FALSE
   }
   meantemp = c()
-  if (monc == 2 | (monc == 4 & dformat == 4) | monc == 5 | (monc == 0 & use.temp == TRUE)) {
-    if (monc == 2) {
+  if (monc == MONITOR$GENEACTIV || (monc == MONITOR$AXIVITY && dformat == FORMAT$CWA) || monc == MONITOR$MOVISENS || (monc == MONITOR$AD_HOC && use.temp == TRUE)) {
+    if (monc == MONITOR$GENEACTIV) {
       if ("temperature" %in% colnames(data)) {
         tempcolumn = which(colnames(data) == "temperature") #GGIRread
       } else {
         tempcolumn = 7
       }
     }
-    if (monc == 4 | monc == 0) tempcolumn = 5
-    if (monc == 5) tempcolumn = 4
+    if (monc == MONITOR$AXIVITY || monc == MONITOR$AD_HOC) tempcolumn = 5
+    if (monc == MONITOR$MOVISENS) tempcolumn = 4
     meantemp = mean(as.numeric(data[, tempcolumn]), na.rm = TRUE)
     if (is.na(meantemp) == T) { #mean(as.numeric(data[1:10,7]))
       warning("temperature is NA", call. = FALSE)
@@ -48,7 +48,7 @@ get_starttime_weekday_meantemp_truncdata = function(temp.available, monc, dforma
     timezone = attr(unclass(as.POSIXlt(starttime[1])), which = "tzone")
     starttimebefore = as.POSIXlt(starttime)
     # assuming that timestamps is good, but that timezone might be lost in conversion from string to POSIXct
-    if (dformat == 1) { #not sure whether this is required for csv-format (2)
+    if (dformat == FORMAT$BIN) { #not sure whether this is required for csv-format (2)
       if (length(which(timezone == "GMT")) > 0) {
         if (length(desiredtz) == 0) {
           warning("desiredtz not specified, local timezoneused as default", call. = FALSE)
@@ -65,16 +65,16 @@ get_starttime_weekday_meantemp_truncdata = function(temp.available, monc, dforma
     wdayname = weekdays[wday]
     #======================================================
     #assess how much data to delete till next 15 minute period
-    temp = unlist(strsplit(format(starttime)," "))
-    if (length(temp) > 1) {
-      starttime2 = as.numeric(unlist(strsplit(temp[2],":")))
-    } else if (length(temp) == 1 & !grepl("T", temp)) {
-      starttime2 = as.numeric(unlist(strsplit(temp[2],":")))
+    tmp = unlist(strsplit(format(starttime)," "))
+    if (length(tmp) > 1) {
+      starttime2 = as.numeric(unlist(strsplit(tmp[2],":")))
+    } else if (length(tmp) == 1 && !grepl("T", tmp)) {
+      starttime2 = as.numeric(unlist(strsplit(tmp[2],":")))
     } else {
       # first get char to POSIX
-      temp = iso8601chartime2POSIX(starttime, tz = desiredtz)
-      temp = unlist(strsplit(format(temp)," ")) # to keep it consistent with what we had
-      starttime2 = as.numeric(unlist(strsplit(format(temp[2]),":")))
+      tmp = iso8601chartime2POSIX(starttime, tz = desiredtz)
+      tmp = unlist(strsplit(format(tmp)," ")) # to keep it consistent with what we had
+      starttime2 = as.numeric(unlist(strsplit(format(tmp[2]),":")))
     }
     if (length(which(is.na(starttime2) ==  TRUE)) > 0 |
         length(starttime2) == 0) {
@@ -109,7 +109,7 @@ get_starttime_weekday_meantemp_truncdata = function(temp.available, monc, dforma
         remem2add24 = TRUE #remember to add 24 hours because this is now the wrong day
       }
     }
-    starttime3 = paste0(temp[1], " ", start_hr, ":", newmin, ":", newsec)
+    starttime3 = paste0(tmp[1], " ", start_hr, ":", newmin, ":", newsec)
     #create timestamp from string (now desiredtz is added)
     if (length(desiredtz) == 0) {
       warning("desiredtz not specified, local timezone used as default", call. = FALSE)

man/GGIR-package.Rd

@@ -21,8 +21,8 @@
   \tabular{ll}{
   Package: \tab GGIR\cr
   Type: \tab Package\cr
-  Version: \tab 2.10-2\cr
-  Date: \tab 2023-09-01\cr
+  Version: \tab 2.10-4\cr
+  Date: \tab 2023-10-05\cr
   License: \tab Apache License (== 2.0)\cr
   Discussion group: \tab https://groups.google.com/forum/#!forum/rpackageggir\cr
   }

tests/testthat/test_greadaccfile.R

@@ -22,8 +22,8 @@ test_that("g.readaccfile and g.inspectfile can read gt3x, cwa, and actigraph csv
   cat("\nActigraph .gt3x")
   # actigraph .gt3x
   Igt3x = g.inspectfile(gt3xfile, desiredtz = desiredtz)
-  expect_equal(Igt3x$monc, 3)
-  expect_equal(Igt3x$dformc, 6)
+  expect_equal(Igt3x$monc, MONITOR$ACTIGRAPH)
+  expect_equal(Igt3x$dformc, FORMAT$GT3X)
   expect_equal(Igt3x$sf, 30)
   EHV = g.extractheadervars(Igt3x)
   expect_equal(EHV$deviceSerialNumber, "MOS2E39180594_firmware_1.9.2")
@@ -34,8 +34,8 @@ test_that("g.readaccfile and g.inspectfile can read gt3x, cwa, and actigraph csv
   cat("\nAxivity .cwa")
   # axivity .cwa
   Icwa = g.inspectfile(cwafile, desiredtz = desiredtz, params_rawdata = params_rawdata)
-  expect_equal(Icwa$monc, 4)
-  expect_equal(Icwa$dformc, 4)
+  expect_equal(Icwa$monc, MONITOR$AXIVITY)
+  expect_equal(Icwa$dformc, FORMAT$CWA)
   expect_equal(Icwa$sf, 100)
   EHV = g.extractheadervars(Icwa)
   expect_equal(EHV$deviceSerialNumber,"39434")
@@ -47,8 +47,8 @@ test_that("g.readaccfile and g.inspectfile can read gt3x, cwa, and actigraph csv
   cat("\nGENEActiv .bin")
   # GENEActiv .bin
   IGA = g.inspectfile(GAfile, desiredtz = desiredtz)
-  expect_equal(IGA$monc,2)
-  expect_equal(IGA$dformc,1)
+  expect_equal(IGA$monc, MONITOR$GENEACTIV)
+  expect_equal(IGA$dformc, FORMAT$BIN)
   expect_equal(IGA$sf,85.7)
 
   EHV = g.extractheadervars(IGA)
@@ -59,9 +59,9 @@ test_that("g.readaccfile and g.inspectfile can read gt3x, cwa, and actigraph csv
 
 
   # test decimal separator recognition extraction
-  decn =  g.dotorcomma(cwafile,dformat = 4, mon = 4, desiredtz = desiredtz)
+  decn =  g.dotorcomma(cwafile,dformat = FORMAT$CWA, mon = MONITOR$AXIVITY, desiredtz = desiredtz)
   expect_equal(decn,".")
-  decn =  g.dotorcomma(GAfile, dformat = 1, mon = 2, desiredtz = desiredtz)
+  decn =  g.dotorcomma(GAfile, dformat = FORMAT$BIN, mon = MONITOR$GENEACTIV, desiredtz = desiredtz)
   expect_equal(decn,".")
   filequality = list(filecorrupt = FALSE, filetooshort = FALSE)
   dayborder = 0

vignettes/GGIR.Rmd

@@ -776,6 +776,12 @@ Part 4 generates the following output:
 -   part4_nightsummary_sleep_cleaned.csv
 -   QC/part4_nightsummary_sleep_full.csv
 
+The latter with '_full' in the name is intended to aid
+clarifying why some nights (if any) are excluded from the cleaned summary report.
+Although, nights where the accelerometer was not worn at all are excluded from this. So,
+if you have a 30 day recording where the accelerometer was not worn from day 7 onward then
+you will not find the last 22 nights in either csv-report.
+
 The csv. files contain the variables as shown below.
 
 | (Part of) variable name    | Description         |
@@ -1837,7 +1843,7 @@ also the intensity gradient. If argument `winhr` is a vector then
 descriptive values for LX and MX are derived per value in `winhr`.
 Within GGIR part 2 MXLX is calculated per calendar day and, if argument
 `qwindow` is specified, per segment of the day. Within GGIR part 5 MXLX
-is calculated per window, and if used in combination with the GENEActiv
+is calculated per window, and if used in combination with the GENEActiv or Axivity
 accelerometer brand LUX estimates per LX and MX are included.
 
 The MX metric describe here should not be confused by the MX metrics as proposed
@@ -2035,6 +2041,27 @@ for midnight. If midnight is not found the sleep detection is skipped.
 - Ideally two valid consecutive nights and the waking hours in between.
 
 
+## LUX sensor data processing
+
+Although GGIR focuses on accelerometer data a few brands come with LUX data.
+
+In part 1 GGIR calculates the peak lux per long epoch at a default resolution of 15 minutes,
+which can be modified with argument windowsizes. Peak light offers a more reliable 
+estimate of light exposure per time window compared with taking the average. Further,
+LUX is used in the auto-calibration.
+
+In GGIR part 2 we visualise the LUX values in the qc plot. In part 3 and 4 LUX is 
+not used for sleep classification because relation between light exposure and sleep is weak.
+
+In part 5 we calculate the mean and maximum of the peak LUX per epoch across all
+waking hours of the day. Here, the mean (peak per epoch) LUX would then indicate
+average light exposure per time segment, while max peak would indicate the maximum
+light exposure per day. Further, we calculate the max and mean peak LUX per most 
+active consecutive X hour of the day. This is intended to offer an alternative to LUX
+exposure during waking hours which relies on correct sleep classification. LUX
+exposure during M10 may be seen as an alternative if you are unsure whether you 
+can trust the sleep classification in your data set.
+
 # Other Resources
 
 -   The [GGIR package manual](https://CRAN.R-project.org/package=GGIR)