Merge pull request #933 from wadpac/issue931_OOwindow

Issue931 onset-to-onset timewindow

NEWS.md

@@ -7,6 +7,9 @@
 - Part 2: Arguments hrs.del.start and hrs.del.end when combined with strategy = 3 and strategy = 5 now count
 relative to start and end of the most active time window as identified.  #905
 
+- Part 5: Added optioned "OO" to argument timewindow, which defines windows from
+sleep Onset to sleep Onset #931
+
 # CHANGES IN GGIR VERSION 2.10-4
 
 - Part 4: Now better able to handle nights without sustained inactivity bouts (rest) #911

R/g.part5.R

@@ -276,8 +276,11 @@ g.part5 = function(datadir = c(), metadatadir = c(), f0=c(), f1=c(),
             # Add first waking up time, if it is missing:
             ts = g.part5.addfirstwake(ts, summarysleep = summarysleep_tmp2, nightsi, sleeplog, ID,
                                       Nepochsinhour, SPTE_end)
+            # Convert time column from iso8601 to POSIX regardless of whether it is aggregated
+            # to ensure the format is consistent
+            ts$time = iso8601chartime2POSIX(ts$time,tz = params_general[["desiredtz"]])
             if (params_general[["part5_agg2_60seconds"]] == TRUE) { # Optionally aggregate to 1 minute epoch:
-              ts$time_num = floor(as.numeric(iso8601chartime2POSIX(ts$time,tz = params_general[["desiredtz"]])) / 60) * 60
+              ts$time_num = floor(as.numeric(ts$time) / 60) * 60
 
               # only include angle if angle is present
               angleColName = ifelse("angle" %in% names(ts), yes = "angle", no = NULL)
@@ -410,6 +413,7 @@ g.part5 = function(datadir = c(), metadatadir = c(), f0=c(), f1=c(),
                   #-------------------------------
                   # ignore all nights in 'inights' before the first waking up and after the last waking up
                   FM = which(diff(ts$diur) == -1)
+                  SO = which(diff(ts$diur) == 1)
                   # now 0.5+6+0.5 midnights and 7 days
                   for (timewindowi in params_output[["timewindow"]]) {
                     nightsi = nightsi_bu
@@ -419,6 +423,11 @@ g.part5 = function(datadir = c(), metadatadir = c(), f0=c(), f1=c(),
                         # ignore first and last midnight because we did not do sleep detection on it
                         nightsi = nightsi[nightsi > FM[1] & nightsi < FM[length(FM)]]
                       }
+                    } else if (timewindowi == "OO") {
+                      if (length(SO) > 0) {
+                        # ignore data before the first sleep onset and after the last sleep onset
+                        nightsi = nightsi[nightsi > SO[1] & nightsi < SO[length(SO)]]
+                      }
                     } else {
                       # if first night is missing then nights needs to align with diur
                       startend_sleep = which(abs(diff(ts$diur)) == 1)
@@ -428,8 +437,10 @@ g.part5 = function(datadir = c(), metadatadir = c(), f0=c(), f1=c(),
                     }
                     if (timewindowi == "MM") {
                       Nwindows = length(nightsi) + 1
-                    } else {
+                    } else if (timewindowi == "WW") {
                       Nwindows = length(which(diff(ts$diur) == -1))
+                    } else if (timewindowi == "OO") {
+                      Nwindows = length(which(diff(ts$diur) == 1))
                     }
                     indjump = 1
                     qqq_backup = c()
@@ -442,7 +453,6 @@ g.part5 = function(datadir = c(), metadatadir = c(), f0=c(), f1=c(),
                     }
                     for (wi in 1:Nwindows) { #loop through 7 windows (+1 to include the data after last awakening)
                       # Define indices of start and end of the day window (e.g. midnight-midnight, or waking-up or wakingup
-
                       defdays = g.part5.definedays(nightsi, wi, indjump,
                                                    nightsi_bu, epochSize = ws3new, qqq_backup, ts, 
                                                    timewindowi, Nwindows, qwindow = params_247[["qwindow"]],

R/g.part5.addfirstwake.R

@@ -1,4 +1,4 @@
-g.part5.addfirstwake =function(ts, summarysleep, nightsi, sleeplog, ID, 
+g.part5.addfirstwake = function(ts, summarysleep, nightsi, sleeplog, ID, 
                                Nepochsinhour, SPTE_end) {
   # Note related to if first and last night were ignored in part 4:
   # - diur lacks the first and last night at this point in the code.
@@ -72,7 +72,7 @@ g.part5.addfirstwake =function(ts, summarysleep, nightsi, sleeplog, ID,
     }
     if (is.na(wake_night1_index)) wake_night1_index = 0
     if (wake_night1_index < firstwake & wake_night1_index > 1 & (wake_night1_index-1) > nightsi[1]) {
-      ts$diur[nightsi[1]:(wake_night1_index-1)] = 1
+      ts$diur[1:(wake_night1_index-1)] = 1
     } else {
       # Person slept only during the afternoon on day 2
       # And there is no sleep data available for the first night
@@ -82,8 +82,8 @@ g.part5.addfirstwake =function(ts, summarysleep, nightsi, sleeplog, ID,
       # and merging of the sleep variables is still consistent with
       # the other recording.
       dummywake = max(firstonset - round(Nepochsinhour/12), nightsi[1] + round(Nepochsinhour * 6))
-      ts$diur[nightsi[1]:dummywake] = 1 
-      ts$nonwear[nightsi[1]:firstonset] = 1
+      ts$diur[1:dummywake] = 1 
+      ts$nonwear[1:firstonset] = 1
     }
   }
   return(ts)

R/g.part5.definedays.R

@@ -74,20 +74,22 @@ g.part5.definedays = function(nightsi, wi, indjump, nightsi_bu,
         }
       }
     } else {
+      qqq = c(NA, NA)
       if (length(qqq_backup) > 1) {
-        # if there is remaining time after previous day...
-        if (Nts > qqq_backup[2]) {
+        # If there is remaining time after previous day...
+        # but only do this if there is less than 24 hours.
+        # This is necessary if sleep is ignored for last night.
+        # In that case the last two calendar days should be ignored
+        # as no sleep onset will be available.
+        if (Nts - qqq_backup[2] < 24 * (60 / epochSize) * 60) {
           qqq = c(qqq_backup[2] + 1, Nts)
         }
-      } else {
-        # else, do not analyse this day
-        qqq = c(NA, NA)
       }
     }
     # in MM, also define segments of the day based on qwindow
     if (!is.na(qqq[1]) & !is.na(qqq[2])) {
       fullQqq = qqq[1]:qqq[2]
-      lastepoch = substr(ts$time[qqq[2]], 12, 19)
+      lastepoch = format(ts$time[qqq[2]],  "%H:%M:%S")
       qnames = NULL
       if (is.data.frame(qwindow)) {
         date_of_interest = substr(ts$time[qqq[1]], 1, 10)
@@ -137,25 +139,27 @@ g.part5.definedays = function(nightsi, wi, indjump, nightsi_bu,
       }
       names(segments) = segments_timing
     }
-  } else if (timewindowi == "WW") {
-    if (wi <= (Nwindows - 1)) { # all full wake to wake days
-      qqq[1] = which(diff(ts$diur) == -1)[wi] + 1
-      qqq[2] = which(diff(ts$diur) == -1)[wi + 1]
+  } else if (timewindowi == "WW" || timewindowi == "OO") {
+    windowEdge = ifelse(timewindowi == "WW", yes = -1, no = 1)
+    if (wi <= (Nwindows - 1)) { # all full windows
+      qqq[1] = which(diff(ts$diur) == windowEdge)[wi] + 1
+      qqq[2] = which(diff(ts$diur) == windowEdge)[wi + 1]
     } else {
-      # time after last reliable waking up (this can be more than 24 hours)
-      # ignore this day, because if the night was ignored for sleep analysis
+      # time after last reliable waking up or onset (this can be more than 24 hours)
+      # ignore this window, because if the night was ignored for sleep analysis
       # then the description of the day in part 5 including that night is
       # not informative.
       qqq = c(NA, NA)
     }
     # build up segments
     if (!is.na(qqq[1]) & !is.na(qqq[2])) {
       segments = list(qqq)
-      start = substr(ts$time[qqq[1]], 12, 19)
-      end = substr(ts$time[qqq[2]], 12, 19)
+      start = format(ts$time[qqq[1]], "%H:%M:%S")
+      end = format(ts$time[qqq[2]], "%H:%M:%S")
       names(segments) = paste(start, end, sep = "-")
-      segments_names = "WW"
+      segments_names = timewindowi
     }
+
   }
   return(invisible(list(qqq = qqq, qqq_backup = qqq_backup, 
                         segments = segments, segments_names = segments_names)))

R/g.part5.onsetwaketiming.R

@@ -2,22 +2,23 @@ g.part5.onsetwaketiming = function(qqq, ts, min, sec, hour, timewindowi) {
   onset = wake = 0
   skiponset = TRUE; skipwake = TRUE
   # Onset index
-  if (timewindowi == "WW") {
-    onseti = c(qqq[1]:qqq[2])[which(diff(ts$diur[qqq[1]:(qqq[2]-1)]) == 1)+1]
-    if (length(onseti) > 1) {
-      onseti = onseti[length(onseti)] # in the case if MM use last onset
-    }
+  if (timewindowi == "OO") {
+    # For OO onset is by definition the start and end of the window
+    onseti = qqq[1] + 1
   } else {
-    onseti = c(qqq[1]:qqq[2])[which(diff(ts$diur[qqq[1]:(qqq[2]-1)]) == 1)+1]
+    # For WW and MM onset needs to be search in the window
+    onseti = c(qqq[1]:qqq[2])[which(diff(ts$diur[qqq[1]:(qqq[2] - 1)]) == 1)] + 1
     if (length(onseti) > 1) {
       onseti = onseti[length(onseti)] # in the case if MM use last onset
     }
   }
   # Wake index
   if (timewindowi == "WW") {
-    wakei = qqq[2]+1
+    # For WW wake is by definition the start and end of the window
+    wakei = qqq[2] + 1
   } else {
-    wakei = c(qqq[1]:qqq[2])[which(diff(ts$diur[qqq[1]:(qqq[2]-1)]) == -1)+1]
+    # For OO and MM wake needs to be search in the window
+    wakei = c(qqq[1]:qqq[2])[which(diff(ts$diur[qqq[1]:(qqq[2] - 1)]) == -1) + 1]
     if (length(wakei) > 1) wakei = wakei[1] # in the case if MM use first wake-up time
   }
   # Onset time

R/g.part5_analyseSegment.R

@@ -48,8 +48,8 @@ g.part5_analyseSegment = function(indexlog, timeList, levelList,
   # The following is to avoid issue with merging sleep variables from part 4
   # This code extract them the time series (ts) object create in g.part5
   # Note that this means that for MM windows there can be multiple or no wake or onsets
-  date = as.Date(ts$time[segStart + 1])
-  if (add_one_day_to_next_date == TRUE & timewindowi == "WW") { # see below for explanation
+  date = as.Date(ts$time[segStart + 1], tz = params_general[["desiredtz"]])
+  if (add_one_day_to_next_date == TRUE & timewindowi %in% c("WW", "OO")) { # see below for explanation
     date = date + 1
     add_one_day_to_next_date = FALSE
   }
@@ -69,7 +69,13 @@ g.part5_analyseSegment = function(indexlog, timeList, levelList,
   skiponset = onsetwaketiming$skiponset; skipwake = onsetwaketiming$skipwake
   if (wake < 24 & timewindowi == "WW") {
     # waking up before midnight means that next WW window
-    # will start a day before the day we refer to when discussing it's SPT
+    # will start a day before the date we refer to when discussing it's SPT
+    # So, for next window we have to do date = date + 1
+    add_one_day_to_next_date = TRUE
+  }
+  if (onset > 24 & timewindowi == "OO") {
+    # onset after midnight means that next OO window
+    # will start a date after the date we refer to when discussing it;s SPT
     # So, for next window we have to do date = date + 1
     add_one_day_to_next_date = TRUE
   }
@@ -421,11 +427,15 @@ g.part5_analyseSegment = function(indexlog, timeList, levelList,
     # LIGHT, IF AVAILABLE
     if ("lightpeak" %in% colnames(ts) & length(params_247[["LUX_day_segments"]]) > 0) {
       # mean LUX
-      dsummary[si,fi] =  round(max(ts$lightpeak[sse[ts$diur[sse] == 0]], na.rm = TRUE), digits = 1)
-      dsummary[si,fi + 1] = round(mean(ts$lightpeak[sse[ts$diur[sse] == 0]], na.rm = TRUE), digits = 1)
-      dsummary[si,fi + 2] = round(mean(ts$lightpeak[sse[ts$diur[sse] == 1]], na.rm = TRUE), digits = 1)
-      dsummary[si,fi + 3] = round(mean(ts$lightpeak[sse[ts$diur[sse] == 0 & ts$ACC[sse] > TRMi]], na.rm = TRUE),
-                                  digits = 1)
+      if (length(which(ts$diur[sse] == 0)) > 0 & length(which(ts$diur[sse] == 1)) > 0) {
+        dsummary[si,fi] =  round(max(ts$lightpeak[sse[ts$diur[sse] == 0]], na.rm = TRUE), digits = 1)
+        dsummary[si,fi + 1] = round(mean(ts$lightpeak[sse[ts$diur[sse] == 0]], na.rm = TRUE), digits = 1)
+        dsummary[si,fi + 2] = round(mean(ts$lightpeak[sse[ts$diur[sse] == 1]], na.rm = TRUE), digits = 1)
+        dsummary[si,fi + 3] = round(mean(ts$lightpeak[sse[ts$diur[sse] == 0 & ts$ACC[sse] > TRMi]], na.rm = TRUE),
+                                    digits = 1)
+      } else {
+        dsummary[si,fi:(fi + 3)] = NA
+      }
       ds_names[fi:(fi + 3)] = c("LUX_max_day", "LUX_mean_day", "LUX_mean_spt", "LUX_mean_day_mvpa"); fi = fi + 4
       # time in LUX ranges
       Nluxt = length(params_247[["LUXthresholds"]])
@@ -440,7 +450,7 @@ g.part5_analyseSegment = function(indexlog, timeList, levelList,
         }
       }
       fi = fi + Nluxt
-      if (timewindowi == "WW") {
+      if (timewindowi %in% c("WW", "OO")) {
         # LUX per segment of the day
         luxperseg = g.part5.lux_persegment(ts, sse,
                                            LUX_day_segments = params_247[["LUX_day_segments"]],
@@ -487,7 +497,6 @@ g.part5_analyseSegment = function(indexlog, timeList, levelList,
                   columnIndex = fi)
   timeList = list(ts = ts,
                   epochSize = ws3new)
-
   invisible(list(
     indexlog = indexlog,
     ds_names = ds_names,

R/g.report.part5.R

@@ -45,7 +45,7 @@ g.report.part5 = function(metadatadir = c(), f0 = c(), f1 = c(), loglocation = c
       }
       # Note: Below we intentionally only sets a criteria on daytime, because for
       # the night time we only need start and end of the SPT window.
-      if (window == "WW") {
+      if (window == "WW" | window == "OO") {
         indices = which(x$nonwear_perc_day <= maxpernwday &
                           x$dur_spt_min > 0 & x$dur_day_min > 0 & include_window == TRUE)
       } else if (window == "MM") {
@@ -153,8 +153,8 @@ g.report.part5 = function(metadatadir = c(), f0 = c(), f1 = c(), loglocation = c
     # split results to different spreadsheets in order to minimize individual
     # filesize and to ease organising dataset
     uwi = as.character(unique(outputfinal$window))
-    if (!all(uwi %in% c("MM", "WW"))) {
-      uwi = c(uwi[uwi %in% c("MM", "WW")], "Segments")
+    if (!all(uwi %in% c("MM", "WW", "OO"))) {
+      uwi = c(uwi[uwi %in% c("MM", "WW", "OO")], "Segments")
     }
     uTRLi = as.character(unique(outputfinal$TRLi))
     uTRMi = as.character(unique(outputfinal$TRMi))
@@ -197,7 +197,7 @@ g.report.part5 = function(metadatadir = c(), f0 = c(), f1 = c(), loglocation = c
                            "-", uTRVi[h3], "-", usleepparam[h4]))
               }
               select_window = as.character(outputfinal$window) == uwi[j]
-              if (!(uwi[j] %in% c("MM", "WW"))) select_window = !(as.character(outputfinal$window) %in% c("MM", "WW"))
+              if (!(uwi[j] %in% c("MM", "WW", "OO"))) select_window = !(as.character(outputfinal$window) %in% c("MM", "WW", "OO"))
               seluwi = which(select_window &
                                as.character(outputfinal$TRLi) == uTRLi[h1] &
                                as.character(outputfinal$TRMi) == uTRMi[h2] &
@@ -210,7 +210,7 @@ g.report.part5 = function(metadatadir = c(), f0 = c(), f1 = c(), loglocation = c
                 CN = colnames(outputfinal)
                 outputfinal2 = outputfinal
                 colnames(outputfinal2) = CN
-                delcol = grep(pattern = "window|TRLi|TRMi|TRVi|sleepparam",
+                delcol = grep(pattern = "TRLi|TRMi|TRVi|sleepparam",
                               x = colnames(outputfinal2))
                 if (uwi[j] != "Segments") {
                   delcol = c(delcol, which(colnames(outputfinal2) == "window"))

man/GGIR.Rd

@@ -1278,7 +1278,7 @@ GGIR(mode = 1:5,
         Character (default = c("MM", "WW")).
         In \link{g.part5}: Timewindow over which summary statistics are derived.
         Value can be "MM" (midnight to midnight), "WW" (waking time to waking time),
-        or both c("MM","WW").}
+        "OO" (sleep onset to sleep onset), or any combination of them.}
 
       \item{save_ms5rawlevels}{
         Boolean (default = FALSE).

tests/testthat/test_chainof5parts.R

@@ -201,7 +201,7 @@ test_that("chainof5parts", {
           overwrite = TRUE, excludefirstlast = FALSE, do.parallel = do.parallel,
           frag.metrics = "all", save_ms5rawlevels = TRUE,
           part5_agg2_60seconds = TRUE, do.sibreport = TRUE, nap_model = "hip3yr",
-          iglevels = 1)
+          iglevels = 1, timewindow = c("MM", "WW", "OO"))
   sibreport_dirname = "output_test/meta/ms5.outraw/sib.reports"
   expect_true(dir.exists(sibreport_dirname))
   expect_true(file.exists(paste0(sibreport_dirname, "/sib_report_123A_testaccfile_T5A5.csv")))
@@ -214,9 +214,11 @@ test_that("chainof5parts", {
 
   expect_true(dir.exists(dirname))
   expect_true(file.exists(rn[1]))
-  expect_that(nrow(output),equals(4))
+  expect_that(nrow(output),equals(5)) # changed because OO window is exported
   expect_that(ncol(output),equals(197))
   expect_that(round(as.numeric(output$wakeup[2]), digits = 4), equals(36))
+  expect_that(as.numeric(output$dur_day_spt_min[4]), equals(1150)) # WW window duration
+  expect_that(as.numeric(output$dur_day_spt_min[5]), equals(1680)) # OO window duration
   dirname_raw = "output_test/meta/ms5.outraw/40_100_400"
   rn2 = dir(dirname_raw,full.names = TRUE, recursive = T)
   expect_true(file.exists(rn2[1]))
@@ -240,6 +242,9 @@ test_that("chainof5parts", {
   expect_true(file.exists("output_test/results/part4_nightsummary_sleep_cleaned.csv"))
   expect_true(file.exists("output_test/results/part4_summary_sleep_cleaned.csv"))
   expect_true(file.exists("output_test/results/file summary reports/Report_123A_testaccfile.csv.pdf"))
+  expect_true(file.exists("output_test/results/part5_daysummary_MM_L40M100V400_T5A5.csv"))
+  expect_true(file.exists("output_test/results/part5_daysummary_WW_L40M100V400_T5A5.csv"))
+  expect_true(file.exists("output_test/results/part5_daysummary_OO_L40M100V400_T5A5.csv"))
   dn = "output_test"
 
   #=======================

tests/testthat/test_part5_qwindow.R

@@ -6,10 +6,9 @@ test_that("g.part5.definedays considers qwindow to generate segments", {
   # set windows
   timestamp = seq.POSIXt(from = as.POSIXct("2022-01-01 00:00:00"), 
                          to = as.POSIXct("2022-01-04 23:59:00"), by = 60)
-  timestamp = POSIXtime2iso8601(timestamp, tz = "")
   ts = data.frame(timestamp = timestamp,
                   ACC = 0)
-  nightsi = grep("00:00:00", ts$timestamp)
+  nightsi = grep("00:00:00", format(ts$timestamp))
   qwindow = c(7, 18)
 
   definedays = g.part5.definedays(nightsi = nightsi, wi = 1, indjump = 1,