From 8069dbacc78ad3397048cf2fd3a57a9c7ea39f1f Mon Sep 17 00:00:00 2001 From: Vincent van Hees Date: Wed, 19 Oct 2022 15:34:44 +0200 Subject: [PATCH 01/18] progress on mapping overlap selfreported nonwear, nap and sibs #687 --- R/g.part5.R | 37 +++++++++++++++++++++++++++++++++++-- R/g.sibreport.R | 29 +++++++++++++++++------------ 2 files changed, 52 insertions(+), 14 deletions(-) diff --git a/R/g.part5.R b/R/g.part5.R index 01fd9ef37..63dfca499 100644 --- a/R/g.part5.R +++ b/R/g.part5.R @@ -329,10 +329,43 @@ g.part5 = function(datadir = c(), metadatadir = c(), f0=c(), f1=c(), } sibreport_fname = paste0(metadatadir,ms5.sibreport,"/sib_report_",fnames.ms3[i],"_",j,".csv") write.csv(x = sibreport, file = sibreport_fname, row.names = FALSE) + # nap/sib/nonwear overlap analysis - # TO DO - + # account for possibility that some of these categories do not exis + # identify overlapping and non-overlapping, (nap-sib, non-wear-sib, sib, nap, nonwear) + # calculate for all five categories number, total duration, mean duration + if (length(sibreport) > 0) { + sibs = which(sibreport$type == "sib") + if (length(sibs) > 0) { + + sibreport$start = as.POSIXlt(sibreport$start, tz = params_general[["desiredtz"]]) + sibreport$end = as.POSIXlt(sibreport$end, tz = params_general[["desiredtz"]]) + + classes = unique(sibreport$type) + selfreport = which(sibreport$type == "nonwear" | sibreport$type == "nap") + sibreport$overlapNonwear = 0 + sibreport$overlapNap = 0 + if (length(selfreport) > 0) { + for (si in sibs) { + for (sr in 1:length(selfreport)) { + if (sibreport$start[si] < sibreport$end[selfreport[sr]] & + sibreport$end[si] > sibreport$start[selfreport[sr]]) { + end_overlap = as.numeric(pmin(sibreport$end[si], sibreport$end[selfreport[sr]])) + start_overlap = as.numeric(pmax(sibreport$start[si], sibreport$start[selfreport[sr]])) + duration_overlap = end_overlap - start_overlap + duration_sib = as.numeric(sibreport$end[si]) - as.numeric(sibreport$start[si]) + if (sibreport$type[selfreport[sr]] == "nonwear") { + sibreport$overlapNonwear[si] = round(100 * (duration_overlap / duration_sib), digits = 1) + } else if (sibreport$type[selfreport[sr]] == "nap") { + sibreport$overlapNap[si] = round(100 * (duration_overlap / duration_sib), digits = 1) + } + } + } + } + } + } + } # nap detection if (params_general[["acc.metric"]] != "ENMO" | params_sleep[["HASIB.algo"]] != "vanHees2015") { diff --git a/R/g.sibreport.R b/R/g.sibreport.R index 87848db80..3733b6d51 100644 --- a/R/g.sibreport.R +++ b/R/g.sibreport.R @@ -5,7 +5,7 @@ g.sibreport = function(ts, ID, epochlength, logs_diaries=c(), desiredtz="") { Nsibs = length(sib_starts) sibreport = c() if (Nsibs > 0) { - sibreport = data.frame(ID= rep(ID, Nsibs), + sibreport = data.frame(ID = rep(ID, Nsibs), type = rep("sib", Nsibs), start = character(Nsibs), end = character(Nsibs), @@ -16,13 +16,13 @@ g.sibreport = function(ts, ID, epochlength, logs_diaries=c(), desiredtz="") { sibreport$start[sibi] = format(ts$time[dayind[sib_starts[sibi]]]) sibreport$end[sibi] = format(ts$time[dayind[sib_ends[sibi]]]) if (is.ISO8601(sibreport$start[sibi])) { - sibreport$start[sibi] = format(iso8601chartime2POSIX(sibreport$start[sibi],tz=desiredtz)) - sibreport$end[sibi] = format(iso8601chartime2POSIX(sibreport$end[sibi],tz=desiredtz)) + sibreport$start[sibi] = format(iso8601chartime2POSIX(sibreport$start[sibi], tz = desiredtz)) + sibreport$end[sibi] = format(iso8601chartime2POSIX(sibreport$end[sibi], tz = desiredtz)) } sibreport$duration[sibi] = ((sib_ends[sibi] - sib_starts[sibi]) + 1) / (60/epochlength) boutind = sib_starts[sibi]:sib_ends[sibi] - minute_before = (sib_starts[sibi]-(60/epochlength)):(sib_starts[sibi]-1) - minute_after = (sib_ends[sibi]+1):(sib_ends[sibi]+(60/epochlength)) + minute_before = (sib_starts[sibi] - (60/epochlength)):(sib_starts[sibi] - 1) + minute_after = (sib_ends[sibi] + 1):(sib_ends[sibi] + (60/epochlength)) if (min(minute_before) > 1) { sibreport$mean_acc_1min_before[sibi] = round(mean(ts$ACC[dayind[minute_before]]), digits = 3) } @@ -41,15 +41,19 @@ g.sibreport = function(ts, ID, epochlength, logs_diaries=c(), desiredtz="") { logreport = c() if (length(log) > 0) { relevant_rows = which(log$ID == ID) + if (length(relevant_rows) > 0) { log = log[relevant_rows,] # extract ID for (i in 1:nrow(log)) { # loop over lines (days) + tmp = log[i,] # convert into timestamps # only attempt if there are at least 2 timestamps to process nonempty = which(tmp[3:ncol(tmp)] != "" & tmp[3:ncol(tmp)] != "NA") if (length(nonempty) > 1) { date = as.Date(as.character(tmp[1,2]), format = dateformat) times = format(unlist(tmp[1,3:ncol(tmp)])) + times = grep(pattern = "NA", value = TRUE, invert = TRUE, x = times) + times = gsub(pattern = " ", replacement = "", x = times) times = times[which(times %in% c("", "NA") == FALSE)] # ignore entries without start and/or end time t_to_remove = c() @@ -66,13 +70,13 @@ g.sibreport = function(ts, ID, epochlength, logs_diaries=c(), desiredtz="") { if (length(times) > 1) { Nevents = floor(length(times) / 2) timestamps = sort(as.POSIXlt(paste0(date, " ", times), tz = desiredtz)) - logreport_tmp = data.frame(ID= rep(ID, Nevents), + logreport_tmp = data.frame(ID = rep(ID, Nevents), type = rep(logname, Nevents), start = rep("", Nevents), end = rep("", Nevents), stringsAsFactors = FALSE) for (bi in 1:Nevents) { - logreport_tmp$start[bi] = format(timestamps[(bi*2)-1]) - logreport_tmp$end[bi] = format(timestamps[(bi*2)]) + logreport_tmp$start[bi] = format(timestamps[(bi * 2) - 1]) + logreport_tmp$end[bi] = format(timestamps[(bi * 2)]) } } if (length(logreport) == 0) { @@ -86,20 +90,21 @@ g.sibreport = function(ts, ID, epochlength, logs_diaries=c(), desiredtz="") { } return(logreport) } - naplogreport = extract_logs(naplog, ID, logname="nap") - nonwearlogreport = extract_logs(nonwearlog, ID, logname="nonwear") + naplogreport = extract_logs(naplog, ID, logname = "nap") + nonwearlogreport = extract_logs(nonwearlog, ID, logname = "nonwear") logreport = sibreport # append all together in one output data.frame if (length(logreport) > 0 & length(naplogreport) > 0) { - logreport = merge(logreport, naplogreport, by=c("ID", "type", "start", "end"), all=TRUE) + logreport = merge(logreport, naplogreport, by = c("ID", "type", "start", "end"), all = TRUE) } else if (length(logreport) == 0 & length(naplogreport) > 0) { logreport = naplogreport } if (length(logreport) > 0 & length(nonwearlogreport) > 0) { - logreport = merge(logreport, nonwearlogreport, by=c("ID", "type", "start", "end"), all=TRUE) + logreport = merge(logreport, nonwearlogreport, by = c("ID", "type", "start", "end"), all = TRUE) } else if (length(logreport) == 0 & length(nonwearlogreport) > 0) { logreport = nonwearlogreport } + } else { logreport = sibreport } From 64c02fc97db32a4816fd43ceb04f5b2fe8fe73c8 Mon Sep 17 00:00:00 2001 From: Vincent van Hees Date: Wed, 19 Oct 2022 17:06:00 +0200 Subject: [PATCH 02/18] splitting nap classification from sibreport analyses, first draft of sibreport anlayses inside g.part5 for initial testing, next step is to move this to separate functions and add unit-tests #687 --- R/g.part5.R | 248 +++++++++++++++++++++++++++++++++++----------------- 1 file changed, 166 insertions(+), 82 deletions(-) diff --git a/R/g.part5.R b/R/g.part5.R index 63dfca499..061fbc0cb 100644 --- a/R/g.part5.R +++ b/R/g.part5.R @@ -62,7 +62,7 @@ g.part5 = function(datadir = c(), metadatadir = c(), f0=c(), f1=c(), #====================================================================== # compile lists of milestone data filenames fnames.ms3 = sort(dir(paste(metadatadir, "/meta/ms3.out", sep = ""))) - + fnames.ms5 = sort(dir(paste(metadatadir, "/meta/ms5.out", sep = ""))) # path to sleeplog milestonedata, if it exists: sleeplogRDA = paste(metadatadir, "/meta/sleeplog.RData", sep = "") @@ -314,7 +314,7 @@ g.part5 = function(datadir = c(), metadatadir = c(), f0=c(), f1=c(), #=============================================== # Use sib.report to classify naps, non-wear and integrate these in time series # does not depend on bout detection criteria or window definitions. - if (params_output[["do.sibreport"]] == TRUE & length(params_sleep[["nap_model"]]) > 0) { + if (params_output[["do.sibreport"]] == TRUE) { if (params_sleep[["sleeplogidnum"]] == TRUE) { IDtmp = as.numeric(ID) } else { @@ -329,90 +329,56 @@ g.part5 = function(datadir = c(), metadatadir = c(), f0=c(), f1=c(), } sibreport_fname = paste0(metadatadir,ms5.sibreport,"/sib_report_",fnames.ms3[i],"_",j,".csv") write.csv(x = sibreport, file = sibreport_fname, row.names = FALSE) - - # nap/sib/nonwear overlap analysis - # account for possibility that some of these categories do not exis - # identify overlapping and non-overlapping, (nap-sib, non-wear-sib, sib, nap, nonwear) - # calculate for all five categories number, total duration, mean duration - if (length(sibreport) > 0) { - sibs = which(sibreport$type == "sib") - if (length(sibs) > 0) { - - sibreport$start = as.POSIXlt(sibreport$start, tz = params_general[["desiredtz"]]) - sibreport$end = as.POSIXlt(sibreport$end, tz = params_general[["desiredtz"]]) - - classes = unique(sibreport$type) - selfreport = which(sibreport$type == "nonwear" | sibreport$type == "nap") - sibreport$overlapNonwear = 0 - sibreport$overlapNap = 0 - if (length(selfreport) > 0) { - for (si in sibs) { - for (sr in 1:length(selfreport)) { - if (sibreport$start[si] < sibreport$end[selfreport[sr]] & - sibreport$end[si] > sibreport$start[selfreport[sr]]) { - end_overlap = as.numeric(pmin(sibreport$end[si], sibreport$end[selfreport[sr]])) - start_overlap = as.numeric(pmax(sibreport$start[si], sibreport$start[selfreport[sr]])) - duration_overlap = end_overlap - start_overlap - duration_sib = as.numeric(sibreport$end[si]) - as.numeric(sibreport$start[si]) - if (sibreport$type[selfreport[sr]] == "nonwear") { - sibreport$overlapNonwear[si] = round(100 * (duration_overlap / duration_sib), digits = 1) - } else if (sibreport$type[selfreport[sr]] == "nap") { - sibreport$overlapNap[si] = round(100 * (duration_overlap / duration_sib), digits = 1) - } - } - } - } - } + if (length(params_sleep[["nap_model"]]) > 0) { + # nap detection + if (params_general[["acc.metric"]] != "ENMO" | + params_sleep[["HASIB.algo"]] != "vanHees2015") { + warning("\nNap classification currently assumes acc.metric = ENMO and HASIB.algo = vanHees2015, so output may not be meaningful") } - } - # nap detection - if (params_general[["acc.metric"]] != "ENMO" | - params_sleep[["HASIB.algo"]] != "vanHees2015") { - warning("\nNap classification currently assumes acc.metric = ENMO and HASIB.algo = vanHees2015, so output may not be meaningful") - } - naps_nonwear = g.part5.classifyNaps(sibreport = sibreport, - desiredtz = params_general[["desiredtz"]], - possible_nap_window = params_sleep[["possible_nap_window"]], - possible_nap_dur = params_sleep[["possible_nap_dur"]], - nap_model = params_sleep[["nap_model"]], - HASIB.algo = params_sleep[["HASIB.algo"]]) - # store in ts object, such that it is exported in as time series - ts$nap1_nonwear2 = 0 - # napsindices = which(naps_nonwear$probability_nap == 1) - # if (length(napsindices) > 0) { - if (length(naps_nonwear) > 0) { - for (nni in 1:nrow(naps_nonwear)) { - nnc_window = which(time_POSIX >= naps_nonwear$start[nni] & time_POSIX <= naps_nonwear$end[nni] & ts$diur == 0) - if (length(nnc_window) > 0) { - if (naps_nonwear$probability_nap[nni] == 1) { - ts$nap1_nonwear2[nnc_window] = 1 # nap - } else if (naps_nonwear$probability_nap[nni] == 0) { - ts$nap1_nonwear2[nnc_window] = 2 # nonwear + naps_nonwear = g.part5.classifyNaps(sibreport = sibreport, + desiredtz = params_general[["desiredtz"]], + possible_nap_window = params_sleep[["possible_nap_window"]], + possible_nap_dur = params_sleep[["possible_nap_dur"]], + nap_model = params_sleep[["nap_model"]], + HASIB.algo = params_sleep[["HASIB.algo"]]) + # store in ts object, such that it is exported in as time series + ts$nap1_nonwear2 = 0 + # napsindices = which(naps_nonwear$probability_nap == 1) + # if (length(napsindices) > 0) { + if (length(naps_nonwear) > 0) { + for (nni in 1:nrow(naps_nonwear)) { + nnc_window = which(time_POSIX >= naps_nonwear$start[nni] & time_POSIX <= naps_nonwear$end[nni] & ts$diur == 0) + if (length(nnc_window) > 0) { + if (naps_nonwear$probability_nap[nni] == 1) { + ts$nap1_nonwear2[nnc_window] = 1 # nap + } else if (naps_nonwear$probability_nap[nni] == 0) { + ts$nap1_nonwear2[nnc_window] = 2 # nonwear + } } } } - } - # impute non-naps episodes (non-wear) - nonwearindices = which(naps_nonwear$probability_nap == 0) - if (length(nonwearindices) > 0) { - for (nni in nonwearindices) { - nwwindow_start = which(time_POSIX >= naps_nonwear$start[nni] & time_POSIX <= naps_nonwear$end[nni] & ts$diur == 0) - if (length(nwwindow_start) > 0) { - Nepochsin24Hours = (60/ws3new) * 60 * 24 - if (nwwindow_start[1] > Nepochsin24Hours) { - nwwindow = nwwindow_start - Nepochsin24Hours # impute time series with preceding day - if (length(which(ts$nap1_nonwear2[nwwindow] == 2)) / length(nwwindow) > 0.5) { - # if there is also a lot of overlap with non-wear there then do next day - nwwindow = nwwindow_start + Nepochsin24Hours + # impute non-naps episodes (non-wear) + nonwearindices = which(naps_nonwear$probability_nap == 0) + if (length(nonwearindices) > 0) { + for (nni in nonwearindices) { + nwwindow_start = which(time_POSIX >= naps_nonwear$start[nni] & time_POSIX <= naps_nonwear$end[nni] & ts$diur == 0) + if (length(nwwindow_start) > 0) { + Nepochsin24Hours = (60/ws3new) * 60 * 24 + if (nwwindow_start[1] > Nepochsin24Hours) { + nwwindow = nwwindow_start - Nepochsin24Hours # impute time series with preceding day + if (length(which(ts$nap1_nonwear2[nwwindow] == 2)) / length(nwwindow) > 0.5) { + # if there is also a lot of overlap with non-wear there then do next day + nwwindow = nwwindow_start + Nepochsin24Hours + } + } else { + nwwindow = nwwindow_start + Nepochsin24Hours # if there is not preceding day use next day } - } else { - nwwindow = nwwindow_start + Nepochsin24Hours # if there is not preceding day use next day - } - if (max(nwwindow) <= nrow(ts)) { # only attempt imputation if possible - # check again that there is not a lot of overlap with non-wear - if (length(which(ts$nap1_nonwear2[nwwindow] == 2)) / length(nwwindow) > 0.5) { - ts$ACC[nwwindow_start] = ts$ACC[nwwindow] # impute + if (max(nwwindow) <= nrow(ts)) { # only attempt imputation if possible + # check again that there is not a lot of overlap with non-wear + if (length(which(ts$nap1_nonwear2[nwwindow] == 2)) / length(nwwindow) > 0.5) { + ts$ACC[nwwindow_start] = ts$ACC[nwwindow] # impute + } } } } @@ -887,6 +853,124 @@ g.part5 = function(datadir = c(), metadatadir = c(), f0=c(), f1=c(), fi = fi + length(luxperseg$values) } } + #======================================================= + # nap/sib/nonwear overlap analysis + #======================================================= + if (params_output[["do.sibreport"]] == TRUE) { + srep_tmp = sibreport[which(sibreport$start > min(ts$time[sse]) & + sibreport$end < max(ts$time[sse])),] + # account for possibility that some of these categories do not exis + # identify overlapping and non-overlapping, (nap-sib, non-wear-sib, sib, nap, nonwear) + # calculate for all five categories number, total duration, mean duration + if (nrow(srep_tmp) > 0) { + sibs = which(srep_tmp$type == "sib") + srep_tmp$overlapNonwear = 0 + srep_tmp$overlapNap = 0 + if (length(sibs) > 0) { + + srep_tmp$start = as.POSIXlt(srep_tmp$start, tz = params_general[["desiredtz"]]) + srep_tmp$end = as.POSIXlt(srep_tmp$end, tz = params_general[["desiredtz"]]) + + classes = unique(srep_tmp$type) + selfreport = which(srep_tmp$type == "nonwear" | srep_tmp$type == "nap") + + if (length(selfreport) > 0) { + for (si in sibs) { + for (sr in 1:length(selfreport)) { + if (srep_tmp$start[si] < srep_tmp$end[selfreport[sr]] & + srep_tmp$end[si] > srep_tmp$start[selfreport[sr]]) { + end_overlap = as.numeric(pmin(srep_tmp$end[si], srep_tmp$end[selfreport[sr]])) + start_overlap = as.numeric(pmax(srep_tmp$start[si], srep_tmp$start[selfreport[sr]])) + duration_overlap = end_overlap - start_overlap + duration_sib = as.numeric(srep_tmp$end[si]) - as.numeric(srep_tmp$start[si]) + if (srep_tmp$type[selfreport[sr]] == "nonwear") { + srep_tmp$overlapNonwear[si] = round(100 * (duration_overlap / duration_sib), digits = 1) + } else if (srep_tmp$type[selfreport[sr]] == "nap") { + srep_tmp$overlapNap[si] = round(100 * (duration_overlap / duration_sib), digits = 1) + } + } + } + } + } + } + sibs_indices = which(srep_tmp$type == "sib") + nap_indices = which(srep_tmp$type == "nap") + nonwear_indices = which(srep_tmp$type == "nonwear") + overlapNap_indices = which(srep_tmp$overlapNap != 0) + overlapNonwear_indices = which(srep_tmp$overlapNonwear != 0) + dsummary[di,fi:(fi + 4)] = c(length(sibs_indices), + length(nap_indices), + length(nonwear_indices), + length(overlapNap_indices), + length(overlapNonwear_indices)) + + ds_names[fi:(fi + 4)] = c("Nsib", "Nsrnap", "Nsrnonw", + "Noverl_sib_srnap", "Noverl_sib_srnonw") + fi = fi + 5 + if (length(sibs_indices) > 0) { + dsummary[di,fi:(fi + 1)] = c(mean(srep_tmp$duration[sibs_indices]), + sum(srep_tmp$duration[sibs_indices])) + } else { + dsummary[di,fi:(fi + 1)] = c(0, 0) + } + ds_names[fi:(fi + 1)] = c("mdur_sib", "tdur_sib") + fi = fi + 2 + + if (length(nap_indices) > 0) { + srep_tmp$duration[nap_indices] = (as.numeric(srep_tmp$end[nap_indices]) - + as.numeric(srep_tmp$start[nap_indices])) / 60 + dsummary[di,fi:(fi + 1)] = c(mean(srep_tmp$duration[nap_indices]), + sum(srep_tmp$duration[nap_indices])) + } else { + dsummary[di,fi:(fi + 1)] = c(0, 0) + } + ds_names[fi:(fi + 1)] = c("mdur_srnap", "tdur_srnap") + fi = fi + 2 + + if (length(nonwear_indices) > 0) { + dsummary[di,fi:(fi + 1)] = c(mean(srep_tmp$duration[nonwear_indices]), + sum(srep_tmp$duration[nonwear_indices])) + } else { + dsummary[di,fi:(fi + 1)] = c(0, 0) + } + ds_names[fi:(fi + 1)] = c("mdur_srnonw", "tdur_srnonw") + fi = fi + 2 + + + calcOverlapPercentage = function(overlap, duration) { + return(sum(overlap * duration) / sum(duration)) + } + if (length(overlapNap_indices) > 0) { + overlap_perc = calcOverlapPercentage(overlap = srep_tmp$overlapNap[overlapNap_indices], + duration = srep_tmp$duration[overlapNap_indices]) + dsummary[di,fi:(fi + 2)] = c(mean(srep_tmp$duration[overlapNap_indices]), + sum(srep_tmp$duration[overlapNap_indices]), + overlap_perc) + } else { + dsummary[di,fi:(fi + 2)] = c(0, 0, 0) + } + ds_names[fi:(fi + 2)] = c("mdur_overl_srnap", "tdur_overl_srnap", "perc_overl_srnap") + fi = fi + 3 + + if (length(overlapNonwear_indices) > 0) { + calcOverlapPercentage(overlap = srep_tmp$overlapNap[overlapNonwear_indices], + duration = srep_tmp$duration[overlapNonwear_indices]) + dsummary[di,fi:(fi + 2)] = c(mean(srep_tmp$duration[overlapNonwear_indices]), + sum(srep_tmp$duration[overlapNonwear_indices]), + mean(srep_tmp$overlapNonwear[overlapNonwear_indices]), + overlap_perc) + } else { + dsummary[di,fi:(fi + 2)] = c(0, 0, 0) + } + ds_names[fi:(fi + 2)] = c("mdur_overl_srnonw", "tdur_overl_srnonw", "perc_overl_srnonw") + fi = fi + 3 + rm(srep_tmp) + } else { + fi = fi + 17 + } + } + + #=============================================== # FOLDER STRUCTURE if (params_output[["storefolderstructure"]] == TRUE) { @@ -943,7 +1027,7 @@ g.part5 = function(datadir = c(), metadatadir = c(), f0=c(), f1=c(), } output = data.frame(dsummary,stringsAsFactors = FALSE) names(output) = ds_names - if (params_cleaning[["excludefirstlast.part5"]] == TRUE) { + if (params_cleaning[["excludefirstlast.part5"]] == TRUE) { output$window_number = as.numeric(output$window_number) cells2exclude = c(which(output$window_number == min(output$window_number,na.rm = TRUE)), which(output$window_number == max(output$window_number,na.rm = TRUE))) @@ -990,7 +1074,7 @@ g.part5 = function(datadir = c(), metadatadir = c(), f0=c(), f1=c(), if (length(output) > 0) { if (nrow(output) > 0) { save(output, tail_expansion_log, file = paste(metadatadir, - ms5.out, "/", fnames.ms3[i], sep = "")) + ms5.out, "/", fnames.ms3[i], sep = "")) } } } From e527e45595dff0c5f02ef7c6e25d3489822e9c04 Mon Sep 17 00:00:00 2001 From: Vincent van Hees Date: Wed, 19 Oct 2022 17:31:19 +0200 Subject: [PATCH 03/18] limit sibanalyses to waking hours only --- R/g.part5.R | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/R/g.part5.R b/R/g.part5.R index 061fbc0cb..98d7588e9 100644 --- a/R/g.part5.R +++ b/R/g.part5.R @@ -857,8 +857,8 @@ g.part5 = function(datadir = c(), metadatadir = c(), f0=c(), f1=c(), # nap/sib/nonwear overlap analysis #======================================================= if (params_output[["do.sibreport"]] == TRUE) { - srep_tmp = sibreport[which(sibreport$start > min(ts$time[sse]) & - sibreport$end < max(ts$time[sse])),] + srep_tmp = sibreport[which(sibreport$start > min(ts$time[sse[ts$diur[sse] == 0]]) & + sibreport$end < max(ts$time[sse[ts$diur[sse] == 0]])),] # account for possibility that some of these categories do not exis # identify overlapping and non-overlapping, (nap-sib, non-wear-sib, sib, nap, nonwear) # calculate for all five categories number, total duration, mean duration From 39d23a24e8b21cf23f6482b16ef5ef226623a78c Mon Sep 17 00:00:00 2001 From: Vincent van Hees Date: Fri, 21 Oct 2022 17:53:40 +0200 Subject: [PATCH 04/18] #687 updating var names, addressing timestamp issue and including extra check --- R/g.part5.R | 207 ++++++++++++++-------------- tests/testthat/test_chainof5parts.R | 2 +- 2 files changed, 104 insertions(+), 105 deletions(-) diff --git a/R/g.part5.R b/R/g.part5.R index 98d7588e9..bfefd97f4 100644 --- a/R/g.part5.R +++ b/R/g.part5.R @@ -329,7 +329,7 @@ g.part5 = function(datadir = c(), metadatadir = c(), f0=c(), f1=c(), } sibreport_fname = paste0(metadatadir,ms5.sibreport,"/sib_report_",fnames.ms3[i],"_",j,".csv") write.csv(x = sibreport, file = sibreport_fname, row.names = FALSE) - + sibreport_backup = sibreport if (length(params_sleep[["nap_model"]]) > 0) { # nap detection if (params_general[["acc.metric"]] != "ENMO" | @@ -857,120 +857,119 @@ g.part5 = function(datadir = c(), metadatadir = c(), f0=c(), f1=c(), # nap/sib/nonwear overlap analysis #======================================================= if (params_output[["do.sibreport"]] == TRUE) { - srep_tmp = sibreport[which(sibreport$start > min(ts$time[sse[ts$diur[sse] == 0]]) & - sibreport$end < max(ts$time[sse[ts$diur[sse] == 0]])),] - # account for possibility that some of these categories do not exis - # identify overlapping and non-overlapping, (nap-sib, non-wear-sib, sib, nap, nonwear) - # calculate for all five categories number, total duration, mean duration - if (nrow(srep_tmp) > 0) { - sibs = which(srep_tmp$type == "sib") - srep_tmp$overlapNonwear = 0 - srep_tmp$overlapNap = 0 - if (length(sibs) > 0) { - + longboutsi = which(sibreport$duration >= 15) + sibreport = sibreport_backup + if (length(longboutsi) > 0) { + sibreport = sibreport[longboutsi,] + srep_tmp = sibreport[which(sibreport$start > min(ts$time[sse[ts$diur[sse] == 0]]) & + sibreport$end < max(ts$time[sse[ts$diur[sse] == 0]])),] + # account for possibility that some of these categories do not exis + # identify overlapping and non-overlapping, (nap-sib, non-wear-sib, sib, nap, nonwear) + # calculate for all five categories number, total duration, mean duration + if (nrow(srep_tmp) > 0) { + sibs = which(srep_tmp$type == "sib") + srep_tmp$overlapNonwear = 0 + srep_tmp$overlapNap = 0 srep_tmp$start = as.POSIXlt(srep_tmp$start, tz = params_general[["desiredtz"]]) srep_tmp$end = as.POSIXlt(srep_tmp$end, tz = params_general[["desiredtz"]]) - - classes = unique(srep_tmp$type) - selfreport = which(srep_tmp$type == "nonwear" | srep_tmp$type == "nap") - - if (length(selfreport) > 0) { - for (si in sibs) { - for (sr in 1:length(selfreport)) { - if (srep_tmp$start[si] < srep_tmp$end[selfreport[sr]] & - srep_tmp$end[si] > srep_tmp$start[selfreport[sr]]) { - end_overlap = as.numeric(pmin(srep_tmp$end[si], srep_tmp$end[selfreport[sr]])) - start_overlap = as.numeric(pmax(srep_tmp$start[si], srep_tmp$start[selfreport[sr]])) - duration_overlap = end_overlap - start_overlap - duration_sib = as.numeric(srep_tmp$end[si]) - as.numeric(srep_tmp$start[si]) - if (srep_tmp$type[selfreport[sr]] == "nonwear") { - srep_tmp$overlapNonwear[si] = round(100 * (duration_overlap / duration_sib), digits = 1) - } else if (srep_tmp$type[selfreport[sr]] == "nap") { - srep_tmp$overlapNap[si] = round(100 * (duration_overlap / duration_sib), digits = 1) + if (length(sibs) > 0) { + classes = unique(srep_tmp$type) + selfreport = which(srep_tmp$type == "nonwear" | srep_tmp$type == "nap") + if (length(selfreport) > 0) { + for (si in sibs) { + for (sr in 1:length(selfreport)) { + if (srep_tmp$start[si] < srep_tmp$end[selfreport[sr]] & + srep_tmp$end[si] > srep_tmp$start[selfreport[sr]]) { + end_overlap = as.numeric(pmin(srep_tmp$end[si], srep_tmp$end[selfreport[sr]])) + start_overlap = as.numeric(pmax(srep_tmp$start[si], srep_tmp$start[selfreport[sr]])) + duration_overlap = end_overlap - start_overlap + duration_sib = as.numeric(srep_tmp$end[si]) - as.numeric(srep_tmp$start[si]) + if (srep_tmp$type[selfreport[sr]] == "nonwear") { + srep_tmp$overlapNonwear[si] = round(100 * (duration_overlap / duration_sib), digits = 1) + } else if (srep_tmp$type[selfreport[sr]] == "nap") { + srep_tmp$overlapNap[si] = round(100 * (duration_overlap / duration_sib), digits = 1) + } } } } } } - } - sibs_indices = which(srep_tmp$type == "sib") - nap_indices = which(srep_tmp$type == "nap") - nonwear_indices = which(srep_tmp$type == "nonwear") - overlapNap_indices = which(srep_tmp$overlapNap != 0) - overlapNonwear_indices = which(srep_tmp$overlapNonwear != 0) - dsummary[di,fi:(fi + 4)] = c(length(sibs_indices), - length(nap_indices), - length(nonwear_indices), - length(overlapNap_indices), - length(overlapNonwear_indices)) - - ds_names[fi:(fi + 4)] = c("Nsib", "Nsrnap", "Nsrnonw", - "Noverl_sib_srnap", "Noverl_sib_srnonw") - fi = fi + 5 - if (length(sibs_indices) > 0) { - dsummary[di,fi:(fi + 1)] = c(mean(srep_tmp$duration[sibs_indices]), - sum(srep_tmp$duration[sibs_indices])) - } else { - dsummary[di,fi:(fi + 1)] = c(0, 0) - } - ds_names[fi:(fi + 1)] = c("mdur_sib", "tdur_sib") - fi = fi + 2 - - if (length(nap_indices) > 0) { - srep_tmp$duration[nap_indices] = (as.numeric(srep_tmp$end[nap_indices]) - - as.numeric(srep_tmp$start[nap_indices])) / 60 - dsummary[di,fi:(fi + 1)] = c(mean(srep_tmp$duration[nap_indices]), - sum(srep_tmp$duration[nap_indices])) - } else { - dsummary[di,fi:(fi + 1)] = c(0, 0) - } - ds_names[fi:(fi + 1)] = c("mdur_srnap", "tdur_srnap") - fi = fi + 2 - - if (length(nonwear_indices) > 0) { - dsummary[di,fi:(fi + 1)] = c(mean(srep_tmp$duration[nonwear_indices]), - sum(srep_tmp$duration[nonwear_indices])) - } else { - dsummary[di,fi:(fi + 1)] = c(0, 0) - } - ds_names[fi:(fi + 1)] = c("mdur_srnonw", "tdur_srnonw") - fi = fi + 2 - - - calcOverlapPercentage = function(overlap, duration) { - return(sum(overlap * duration) / sum(duration)) - } - if (length(overlapNap_indices) > 0) { - overlap_perc = calcOverlapPercentage(overlap = srep_tmp$overlapNap[overlapNap_indices], - duration = srep_tmp$duration[overlapNap_indices]) - dsummary[di,fi:(fi + 2)] = c(mean(srep_tmp$duration[overlapNap_indices]), - sum(srep_tmp$duration[overlapNap_indices]), - overlap_perc) - } else { - dsummary[di,fi:(fi + 2)] = c(0, 0, 0) - } - ds_names[fi:(fi + 2)] = c("mdur_overl_srnap", "tdur_overl_srnap", "perc_overl_srnap") - fi = fi + 3 - - if (length(overlapNonwear_indices) > 0) { - calcOverlapPercentage(overlap = srep_tmp$overlapNap[overlapNonwear_indices], - duration = srep_tmp$duration[overlapNonwear_indices]) - dsummary[di,fi:(fi + 2)] = c(mean(srep_tmp$duration[overlapNonwear_indices]), - sum(srep_tmp$duration[overlapNonwear_indices]), - mean(srep_tmp$overlapNonwear[overlapNonwear_indices]), - overlap_perc) + + sibs_indices = which(srep_tmp$type == "sib") + nap_indices = which(srep_tmp$type == "nap") + nonwear_indices = which(srep_tmp$type == "nonwear") + overlapNap_indices = which(srep_tmp$overlapNap != 0) + overlapNonwear_indices = which(srep_tmp$overlapNonwear != 0) + dsummary[di,fi:(fi + 4)] = c(length(sibs_indices), + length(nap_indices), + length(nonwear_indices), + length(overlapNap_indices), + length(overlapNonwear_indices)) + + ds_names[fi:(fi + 4)] = c("nbouts_day_sib", "nbouts_day_srnap", "nbouts_day_srnonw", + "noverl_sib_srnap", "noverl_sib_srnonw") + fi = fi + 5 + if (length(sibs_indices) > 0) { + dsummary[di,fi:(fi + 1)] = c(mean(srep_tmp$duration[sibs_indices]), + sum(srep_tmp$duration[sibs_indices])) + } else { + dsummary[di,fi:(fi + 1)] = c(0, 0) + } + ds_names[fi:(fi + 1)] = c("frag_mean_dur_sib_day", "dur_day_sib_min") + fi = fi + 2 + if (length(nap_indices) > 0) { + srep_tmp$duration[nap_indices] = (as.numeric(srep_tmp$end[nap_indices]) - + as.numeric(srep_tmp$start[nap_indices])) / 60 + dsummary[di,fi:(fi + 1)] = c(mean(srep_tmp$duration[nap_indices]), + sum(srep_tmp$duration[nap_indices])) + } else { + dsummary[di,fi:(fi + 1)] = c(0, 0) + } + ds_names[fi:(fi + 1)] = c("frag_mean_dur_srnap_day", "dur_day_srnap_min") + fi = fi + 2 + + if (length(nonwear_indices) > 0) { + dsummary[di,fi:(fi + 1)] = c(mean(srep_tmp$duration[nonwear_indices]), + sum(srep_tmp$duration[nonwear_indices])) + } else { + dsummary[di,fi:(fi + 1)] = c(0, 0) + } + ds_names[fi:(fi + 1)] = c("frag_mean_dur_srnonw_day", "dur_day_srnonw_min") + fi = fi + 2 + + calcOverlapPercentage = function(overlap, duration) { + return(sum(overlap * duration) / sum(duration)) + } + if (length(overlapNap_indices) > 0) { + overlap_perc = calcOverlapPercentage(overlap = srep_tmp$overlapNap[overlapNap_indices], + duration = srep_tmp$duration[overlapNap_indices]) + dsummary[di,fi:(fi + 2)] = c(mean(srep_tmp$duration[overlapNap_indices]), + sum(srep_tmp$duration[overlapNap_indices]), + overlap_perc) + } else { + dsummary[di,fi:(fi + 2)] = c(0, 0, 0) + } + ds_names[fi:(fi + 2)] = c("mdur_overl_srnap", "tdur_overl_srnap", "perc_overl_srnap") + fi = fi + 3 + + if (length(overlapNonwear_indices) > 0) { + calcOverlapPercentage(overlap = srep_tmp$overlapNap[overlapNonwear_indices], + duration = srep_tmp$duration[overlapNonwear_indices]) + dsummary[di,fi:(fi + 2)] = c(mean(srep_tmp$duration[overlapNonwear_indices]), + sum(srep_tmp$duration[overlapNonwear_indices]), + mean(srep_tmp$overlapNonwear[overlapNonwear_indices]), + overlap_perc) + } else { + dsummary[di,fi:(fi + 2)] = c(0, 0, 0) + } + ds_names[fi:(fi + 2)] = c("mdur_overl_srnonw", "tdur_overl_srnonw", "perc_overl_srnonw") + fi = fi + 3 + rm(srep_tmp) } else { - dsummary[di,fi:(fi + 2)] = c(0, 0, 0) + fi = fi + 17 } - ds_names[fi:(fi + 2)] = c("mdur_overl_srnonw", "tdur_overl_srnonw", "perc_overl_srnonw") - fi = fi + 3 - rm(srep_tmp) - } else { - fi = fi + 17 } } - - #=============================================== # FOLDER STRUCTURE if (params_output[["storefolderstructure"]] == TRUE) { diff --git a/tests/testthat/test_chainof5parts.R b/tests/testthat/test_chainof5parts.R index a7ee84967..ef8df9b1c 100644 --- a/tests/testthat/test_chainof5parts.R +++ b/tests/testthat/test_chainof5parts.R @@ -135,7 +135,7 @@ test_that("chainof5parts", { expect_true(dir.exists(dirname)) expect_true(file.exists(rn[1])) expect_that(nrow(output),equals(3)) # changed because part5 now gives also first and last day - expect_that(ncol(output),equals(154)) + expect_that(ncol(output),equals(171)) expect_that(round(as.numeric(output$wakeup[2]), digits = 4), equals(36)) dirname_raw = "output_test/meta/ms5.outraw/40_100_400" rn2 = dir(dirname_raw,full.names = TRUE, recursive = T) From 409670fbb2b116b26453d84b1a0900be7917dadb Mon Sep 17 00:00:00 2001 From: Vincent van Hees Date: Wed, 26 Oct 2022 13:37:00 +0200 Subject: [PATCH 05/18] #687 sibreport now has duration for selfreported naps and nonwear and nap, more overlap variables added to output, sibreport items changed --- R/g.part5.R | 137 +++++++++++++++++++++------- R/g.sibreport.R | 15 +-- tests/testthat/test_chainof5parts.R | 2 +- tests/testthat/test_sibreport.R | 29 +++--- 4 files changed, 130 insertions(+), 53 deletions(-) diff --git a/R/g.part5.R b/R/g.part5.R index bfefd97f4..4eb0061cc 100644 --- a/R/g.part5.R +++ b/R/g.part5.R @@ -857,39 +857,71 @@ g.part5 = function(datadir = c(), metadatadir = c(), f0=c(), f1=c(), # nap/sib/nonwear overlap analysis #======================================================= if (params_output[["do.sibreport"]] == TRUE) { - longboutsi = which(sibreport$duration >= 15) sibreport = sibreport_backup + longboutsi = which(sibreport$type == "sib" | (sibreport$type != "sib" & sibreport$duration >= 15)) + # for qc purposes: + dsummary[di,fi] = length(longboutsi) + ds_names[fi] = "sibreport_n_items" + fi = fi + 1 if (length(longboutsi) > 0) { sibreport = sibreport[longboutsi,] srep_tmp = sibreport[which(sibreport$start > min(ts$time[sse[ts$diur[sse] == 0]]) & sibreport$end < max(ts$time[sse[ts$diur[sse] == 0]])),] + # account for possibility that some of these categories do not exis # identify overlapping and non-overlapping, (nap-sib, non-wear-sib, sib, nap, nonwear) # calculate for all five categories number, total duration, mean duration + # for qc purposes: + dsummary[di,fi] = nrow(srep_tmp) + ds_names[fi] = "sibreport_n_items_day" + fi = fi + 1 + if (nrow(srep_tmp) > 0) { sibs = which(srep_tmp$type == "sib") - srep_tmp$overlapNonwear = 0 - srep_tmp$overlapNap = 0 + srep_tmp$SIBoverlapNonwear = 0 + srep_tmp$SIBoverlapNap = 0 + srep_tmp$NonwearOverlapSIB = 0 + srep_tmp$NapOverlapSIB = 0 srep_tmp$start = as.POSIXlt(srep_tmp$start, tz = params_general[["desiredtz"]]) srep_tmp$end = as.POSIXlt(srep_tmp$end, tz = params_general[["desiredtz"]]) + # # for qc purposes: + # dsummary[di,fi] = nrow(srep_tmp) + # ds_names[fi] = "n_sibs_sibreport" + fi = fi + 1 if (length(sibs) > 0) { classes = unique(srep_tmp$type) selfreport = which(srep_tmp$type == "nonwear" | srep_tmp$type == "nap") if (length(selfreport) > 0) { for (si in sibs) { - for (sr in 1:length(selfreport)) { - if (srep_tmp$start[si] < srep_tmp$end[selfreport[sr]] & - srep_tmp$end[si] > srep_tmp$start[selfreport[sr]]) { - end_overlap = as.numeric(pmin(srep_tmp$end[si], srep_tmp$end[selfreport[sr]])) - start_overlap = as.numeric(pmax(srep_tmp$start[si], srep_tmp$start[selfreport[sr]])) + for (sr in selfreport) { + if (srep_tmp$start[si] <= srep_tmp$end[sr] & + srep_tmp$end[si] >= srep_tmp$start[sr]) { + end_overlap = as.numeric(pmin(srep_tmp$end[si], srep_tmp$end[sr])) + start_overlap = as.numeric(pmax(srep_tmp$start[si], srep_tmp$start[sr])) duration_overlap = end_overlap - start_overlap duration_sib = as.numeric(srep_tmp$end[si]) - as.numeric(srep_tmp$start[si]) - if (srep_tmp$type[selfreport[sr]] == "nonwear") { - srep_tmp$overlapNonwear[si] = round(100 * (duration_overlap / duration_sib), digits = 1) - } else if (srep_tmp$type[selfreport[sr]] == "nap") { - srep_tmp$overlapNap[si] = round(100 * (duration_overlap / duration_sib), digits = 1) + perc_overlap = round(100 * (duration_overlap / duration_sib), digits = 1) + if (srep_tmp$type[sr] == "nonwear") { + srep_tmp$SIBoverlapNonwear[si] = perc_overlap + } else if (srep_tmp$type[sr] == "nap") { + srep_tmp$SIBoverlapNap[si] = perc_overlap } } + + if (srep_tmp$start[sr] <= srep_tmp$end[si] & + srep_tmp$end[sr] >= srep_tmp$start[si]) { + end_overlap = as.numeric(pmin(srep_tmp$end[si], srep_tmp$end[sr])) + start_overlap = as.numeric(pmax(srep_tmp$start[si], srep_tmp$start[sr])) + duration_overlap = end_overlap - start_overlap + duration_sr = as.numeric(srep_tmp$end[sr]) - as.numeric(srep_tmp$start[sr]) + perc_overlap = round(100 * (duration_overlap / duration_sr), digits = 1) + if (srep_tmp$type[sr] == "nonwear") { + srep_tmp$NonwearOverlapSIB[sr] = perc_overlap + } else if (srep_tmp$type[sr] == "nap") { + srep_tmp$NapOverlapSIB[sr] = perc_overlap + } + } + } } } @@ -898,17 +930,22 @@ g.part5 = function(datadir = c(), metadatadir = c(), f0=c(), f1=c(), sibs_indices = which(srep_tmp$type == "sib") nap_indices = which(srep_tmp$type == "nap") nonwear_indices = which(srep_tmp$type == "nonwear") - overlapNap_indices = which(srep_tmp$overlapNap != 0) - overlapNonwear_indices = which(srep_tmp$overlapNonwear != 0) - dsummary[di,fi:(fi + 4)] = c(length(sibs_indices), + SIBoverlapNap_indices = which(srep_tmp$SIBoverlapNap != 0) + SIBoverlapNonwear_indices = which(srep_tmp$SIBoverlapNonwear != 0) + NapOverlapSIB_indices = which(srep_tmp$NapOverlapSIB != 0) + NonwearOverlapSIB_indices = which(srep_tmp$NonwearOverlapSIB != 0) + dsummary[di,fi:(fi + 6)] = c(length(sibs_indices), length(nap_indices), length(nonwear_indices), - length(overlapNap_indices), - length(overlapNonwear_indices)) + length(SIBoverlapNap_indices), + length(SIBoverlapNonwear_indices), + length(NapOverlapSIB_indices), + length(NonwearOverlapSIB_indices)) - ds_names[fi:(fi + 4)] = c("nbouts_day_sib", "nbouts_day_srnap", "nbouts_day_srnonw", - "noverl_sib_srnap", "noverl_sib_srnonw") - fi = fi + 5 + ds_names[fi:(fi + 6)] = c("nbouts_day_sib", "nbouts_day_srnap", "nbouts_day_srnonw", + "noverl_sib_srnap", "noverl_sib_srnonw", + "noverl_srnap_sib", "noverl_srnonw_sib") + fi = fi + 7 if (length(sibs_indices) > 0) { dsummary[di,fi:(fi + 1)] = c(mean(srep_tmp$duration[sibs_indices]), sum(srep_tmp$duration[sibs_indices])) @@ -925,6 +962,7 @@ g.part5 = function(datadir = c(), metadatadir = c(), f0=c(), f1=c(), } else { dsummary[di,fi:(fi + 1)] = c(0, 0) } + ds_names[fi:(fi + 1)] = c("frag_mean_dur_srnap_day", "dur_day_srnap_min") fi = fi + 2 @@ -936,37 +974,66 @@ g.part5 = function(datadir = c(), metadatadir = c(), f0=c(), f1=c(), } ds_names[fi:(fi + 1)] = c("frag_mean_dur_srnonw_day", "dur_day_srnonw_min") fi = fi + 2 - + # Overlap sib with srnap calcOverlapPercentage = function(overlap, duration) { return(sum(overlap * duration) / sum(duration)) } - if (length(overlapNap_indices) > 0) { - overlap_perc = calcOverlapPercentage(overlap = srep_tmp$overlapNap[overlapNap_indices], - duration = srep_tmp$duration[overlapNap_indices]) - dsummary[di,fi:(fi + 2)] = c(mean(srep_tmp$duration[overlapNap_indices]), - sum(srep_tmp$duration[overlapNap_indices]), + if (length(SIBoverlapNap_indices) > 0) { + overlap_perc = calcOverlapPercentage(overlap = srep_tmp$SIBoverlapNap[SIBoverlapNap_indices], + duration = srep_tmp$duration[SIBoverlapNap_indices]) + dsummary[di,fi:(fi + 2)] = c(mean(srep_tmp$duration[SIBoverlapNap_indices]), + sum(srep_tmp$duration[SIBoverlapNap_indices]), overlap_perc) } else { dsummary[di,fi:(fi + 2)] = c(0, 0, 0) } - ds_names[fi:(fi + 2)] = c("mdur_overl_srnap", "tdur_overl_srnap", "perc_overl_srnap") + ds_names[fi:(fi + 2)] = c("mdur_sib_overl_srnap", "tdur_sib_overl_srnap", "perc_sib_overl_srnap") fi = fi + 3 - if (length(overlapNonwear_indices) > 0) { - calcOverlapPercentage(overlap = srep_tmp$overlapNap[overlapNonwear_indices], - duration = srep_tmp$duration[overlapNonwear_indices]) - dsummary[di,fi:(fi + 2)] = c(mean(srep_tmp$duration[overlapNonwear_indices]), - sum(srep_tmp$duration[overlapNonwear_indices]), - mean(srep_tmp$overlapNonwear[overlapNonwear_indices]), + # Overlap srnap with sib + calcOverlapPercentage = function(overlap, duration) { + return(sum(overlap * duration) / sum(duration)) + } + if (length(SIBoverlapNap_indices) > 0) { + overlap_perc = calcOverlapPercentage(overlap = srep_tmp$NapOverlapSIB[NapOverlapSIB_indices], + duration = srep_tmp$duration[NapOverlapSIB_indices]) + dsummary[di,fi:(fi + 2)] = c(mean(srep_tmp$duration[NapOverlapSIB_indices]), + sum(srep_tmp$duration[NapOverlapSIB_indices]), + overlap_perc) + } else { + dsummary[di,fi:(fi + 2)] = c(0, 0, 0) + } + ds_names[fi:(fi + 2)] = c("mdur_srnap_overl_sib", "tdur_srnap_overl_sib", "perc_srnap_overl_sib") + fi = fi + 3 + # Overlap sib with srnonw + if (length(SIBoverlapNonwear_indices) > 0) { + calcOverlapPercentage(overlap = srep_tmp$SIBoverlapNonwear[SIBoverlapNonwear_indices], + duration = srep_tmp$duration[SIBoverlapNonwear_indices]) + dsummary[di,fi:(fi + 2)] = c(mean(srep_tmp$duration[SIBoverlapNonwear_indices]), + sum(srep_tmp$duration[SIBoverlapNonwear_indices]), + mean(srep_tmp$SIBoverlapNonwear[SIBoverlapNonwear_indices]), + overlap_perc) + } else { + dsummary[di,fi:(fi + 2)] = c(0, 0, 0) + } + ds_names[fi:(fi + 2)] = c("mdur_sib_overl_srnonw", "tdur_sib_overl_srnonw", "perc_sib_overl_srnonw") + fi = fi + 3 + # Overlap srnonw with sib + if (length(SIBoverlapNonwear_indices) > 0) { + calcOverlapPercentage(overlap = srep_tmp$NonwearOverlapSIB[NonwearOverlapSIB_indices], + duration = srep_tmp$duration[NonwearOverlapSIB_indices]) + dsummary[di,fi:(fi + 2)] = c(mean(srep_tmp$duration[NonwearOverlapSIB_indices]), + sum(srep_tmp$duration[NonwearOverlapSIB_indices]), + mean(srep_tmp$NonwearOverlapSIB[NonwearOverlapSIB_indices]), overlap_perc) } else { dsummary[di,fi:(fi + 2)] = c(0, 0, 0) } - ds_names[fi:(fi + 2)] = c("mdur_overl_srnonw", "tdur_overl_srnonw", "perc_overl_srnonw") + ds_names[fi:(fi + 2)] = c("mdur_srnonw_overl_sib", "tdur_srnonw_overl_sib", "perc_srnonw_overl_sib") fi = fi + 3 rm(srep_tmp) } else { - fi = fi + 17 + fi = fi + 27 } } } diff --git a/R/g.sibreport.R b/R/g.sibreport.R index 3733b6d51..ff86ea4d8 100644 --- a/R/g.sibreport.R +++ b/R/g.sibreport.R @@ -73,10 +73,14 @@ g.sibreport = function(ts, ID, epochlength, logs_diaries=c(), desiredtz="") { logreport_tmp = data.frame(ID = rep(ID, Nevents), type = rep(logname, Nevents), start = rep("", Nevents), - end = rep("", Nevents), stringsAsFactors = FALSE) + end = rep("", Nevents), + duration = rep(0, Nevents), stringsAsFactors = FALSE) for (bi in 1:Nevents) { - logreport_tmp$start[bi] = format(timestamps[(bi * 2) - 1]) - logreport_tmp$end[bi] = format(timestamps[(bi * 2)]) + tt1 = as.POSIXlt(timestamps[(bi * 2) - 1], tz = desiredtz) + tt2 = as.POSIXlt(timestamps[(bi * 2)], tz = desiredtz) + logreport_tmp$start[bi] = format(tt1) + logreport_tmp$end[bi] = format(tt2) + logreport_tmp$duration[bi] = abs(as.numeric(difftime(time1 = tt1, time2 = tt2, units = "mins"))) } } if (length(logreport) == 0) { @@ -95,16 +99,15 @@ g.sibreport = function(ts, ID, epochlength, logs_diaries=c(), desiredtz="") { logreport = sibreport # append all together in one output data.frame if (length(logreport) > 0 & length(naplogreport) > 0) { - logreport = merge(logreport, naplogreport, by = c("ID", "type", "start", "end"), all = TRUE) + logreport = merge(logreport, naplogreport, by = c("ID", "type", "start", "end", "duration"), all = TRUE) } else if (length(logreport) == 0 & length(naplogreport) > 0) { logreport = naplogreport } if (length(logreport) > 0 & length(nonwearlogreport) > 0) { - logreport = merge(logreport, nonwearlogreport, by = c("ID", "type", "start", "end"), all = TRUE) + logreport = merge(logreport, nonwearlogreport, by = c("ID", "type", "start", "end", "duration"), all = TRUE) } else if (length(logreport) == 0 & length(nonwearlogreport) > 0) { logreport = nonwearlogreport } - } else { logreport = sibreport } diff --git a/tests/testthat/test_chainof5parts.R b/tests/testthat/test_chainof5parts.R index ef8df9b1c..90f937036 100644 --- a/tests/testthat/test_chainof5parts.R +++ b/tests/testthat/test_chainof5parts.R @@ -135,7 +135,7 @@ test_that("chainof5parts", { expect_true(dir.exists(dirname)) expect_true(file.exists(rn[1])) expect_that(nrow(output),equals(3)) # changed because part5 now gives also first and last day - expect_that(ncol(output),equals(171)) + expect_that(ncol(output),equals(181)) expect_that(round(as.numeric(output$wakeup[2]), digits = 4), equals(36)) dirname_raw = "output_test/meta/ms5.outraw/40_100_400" rn2 = dir(dirname_raw,full.names = TRUE, recursive = T) diff --git a/tests/testthat/test_sibreport.R b/tests/testthat/test_sibreport.R index dfa0262f6..7b2680479 100644 --- a/tests/testthat/test_sibreport.R +++ b/tests/testthat/test_sibreport.R @@ -2,22 +2,29 @@ library(GGIR) context("g.sibreportc") test_that("g.sibreport creates expected object", { ts = data.frame(sibdetection = c(0, 0, 1, 1, 1, 0, 0), - diur = rep(0, 7), - time = seq(as.POSIXlt(x = "2021-3-3 15:00:00", tz="Europe/Amsterdam"), - as.POSIXlt(x = "2021-3-3 15:06:00", tz="Europe/Amsterdam"), by = 60), - ACC = 1:7, angle = 7:1) + diur = rep(0, 7), + time = seq(as.POSIXlt(x = "2021-3-3 15:00:00", tz = "Europe/Amsterdam"), + as.POSIXlt(x = "2021-3-3 15:06:00", tz = "Europe/Amsterdam"), + by = 60), + ACC = 1:7, angle = 7:1) ID = 12345 epochlength = 60 - nonwearlog = data.frame(ID=rep(ID, 5), date=seq(as.Date("2021-3-3"),as.Date("2021-3-7"), by=1), - nonwear1=rep("15:00:00", 5), nonwear2=rep("15:30:00", 5)) - naplog = data.frame(ID=rep(ID, 2), date=seq(as.Date("2021-3-3"),as.Date("2021-3-4"), by=1), - nonwear1=rep("13:00:00", 2), nonwear2=rep("13:30:00", 2)) + nonwearlog = data.frame(ID = rep(ID, 5), + date = seq(as.Date("2021-3-3"), as.Date("2021-3-7"), by = 1), + nonwear1 = rep("15:00:00", 5), + nonwear2 = rep("15:30:00", 5)) + naplog = data.frame(ID = rep(ID, 2), + date = seq(as.Date("2021-3-3"),as.Date("2021-3-4"), by = 1), + nonwear1 = rep("13:00:00", 2), + nonwear2 = rep("13:30:00", 2)) + + logs_diaries = list(sleeplog = c(), nonwearlog = nonwearlog, + naplog = naplog, dateformat = "%Y-%m-%d") - logs_diaries = list(sleeplog=c(), nonwearlog=nonwearlog, - naplog=naplog, dateformat="%Y-%m-%d") SIBREPORT = g.sibreport(ts, ID, epochlength, logs_diaries) + expect_equal(nrow(SIBREPORT), 8) - expect_equal(mean(SIBREPORT$duration, na.rm = TRUE), 3) + expect_equal(mean(SIBREPORT$duration, na.rm = TRUE), 26.625) expect_equal(mean(SIBREPORT$mean_acc_1min_before, na.rm = TRUE), 2) expect_equal(mean(SIBREPORT$mean_acc_1min_after, na.rm = TRUE), 6) expect_equal(SIBREPORT$start[8], "2021-03-03 15:02:00") From 06aeeb771c71e4e4e7b7b5993c14b34fbd06e860 Mon Sep 17 00:00:00 2001 From: Vincent van Hees Date: Wed, 26 Oct 2022 17:15:45 +0200 Subject: [PATCH 06/18] #687 move code to separate function and add unit test --- NAMESPACE | 3 +- R/g.part5.R | 190 ++---------------------- R/g.part5.analyseRest.R | 182 +++++++++++++++++++++++ man/g.part5.analyseRest.Rd | 44 ++++++ tests/testthat/test_chainof5parts.R | 2 +- tests/testthat/test_part5_analyseRest.R | 100 +++++++++++++ 6 files changed, 339 insertions(+), 182 deletions(-) create mode 100644 R/g.part5.analyseRest.R create mode 100644 man/g.part5.analyseRest.Rd create mode 100644 tests/testthat/test_part5_analyseRest.R diff --git a/NAMESPACE b/NAMESPACE index 5a7b9da8e..50ebab841 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -36,7 +36,8 @@ export(g.analyse, g.calibrate, check_params, extract_params, g.imputeTimegaps, g.part5.classifyNaps, tidyup_df, cosinorAnalyses, - ShellDoc2Vignette, parametersVignette) + ShellDoc2Vignette, parametersVignette, + g.part5.analyseRest) importFrom("grDevices", "colors", "dev.off", "pdf","rainbow","rgb") importFrom("graphics", "abline", "axis", "par", "plot", "plot.new", "rect","axis.POSIXct", "barplot", "box", "legend", diff --git a/R/g.part5.R b/R/g.part5.R index 4eb0061cc..359a0acdb 100644 --- a/R/g.part5.R +++ b/R/g.part5.R @@ -329,7 +329,7 @@ g.part5 = function(datadir = c(), metadatadir = c(), f0=c(), f1=c(), } sibreport_fname = paste0(metadatadir,ms5.sibreport,"/sib_report_",fnames.ms3[i],"_",j,".csv") write.csv(x = sibreport, file = sibreport_fname, row.names = FALSE) - sibreport_backup = sibreport + # sibreport_backup = sibreport if (length(params_sleep[["nap_model"]]) > 0) { # nap detection if (params_general[["acc.metric"]] != "ENMO" | @@ -857,185 +857,15 @@ g.part5 = function(datadir = c(), metadatadir = c(), f0=c(), f1=c(), # nap/sib/nonwear overlap analysis #======================================================= if (params_output[["do.sibreport"]] == TRUE) { - sibreport = sibreport_backup - longboutsi = which(sibreport$type == "sib" | (sibreport$type != "sib" & sibreport$duration >= 15)) - # for qc purposes: - dsummary[di,fi] = length(longboutsi) - ds_names[fi] = "sibreport_n_items" - fi = fi + 1 - if (length(longboutsi) > 0) { - sibreport = sibreport[longboutsi,] - srep_tmp = sibreport[which(sibreport$start > min(ts$time[sse[ts$diur[sse] == 0]]) & - sibreport$end < max(ts$time[sse[ts$diur[sse] == 0]])),] - - # account for possibility that some of these categories do not exis - # identify overlapping and non-overlapping, (nap-sib, non-wear-sib, sib, nap, nonwear) - # calculate for all five categories number, total duration, mean duration - # for qc purposes: - dsummary[di,fi] = nrow(srep_tmp) - ds_names[fi] = "sibreport_n_items_day" - fi = fi + 1 - - if (nrow(srep_tmp) > 0) { - sibs = which(srep_tmp$type == "sib") - srep_tmp$SIBoverlapNonwear = 0 - srep_tmp$SIBoverlapNap = 0 - srep_tmp$NonwearOverlapSIB = 0 - srep_tmp$NapOverlapSIB = 0 - srep_tmp$start = as.POSIXlt(srep_tmp$start, tz = params_general[["desiredtz"]]) - srep_tmp$end = as.POSIXlt(srep_tmp$end, tz = params_general[["desiredtz"]]) - # # for qc purposes: - # dsummary[di,fi] = nrow(srep_tmp) - # ds_names[fi] = "n_sibs_sibreport" - fi = fi + 1 - if (length(sibs) > 0) { - classes = unique(srep_tmp$type) - selfreport = which(srep_tmp$type == "nonwear" | srep_tmp$type == "nap") - if (length(selfreport) > 0) { - for (si in sibs) { - for (sr in selfreport) { - if (srep_tmp$start[si] <= srep_tmp$end[sr] & - srep_tmp$end[si] >= srep_tmp$start[sr]) { - end_overlap = as.numeric(pmin(srep_tmp$end[si], srep_tmp$end[sr])) - start_overlap = as.numeric(pmax(srep_tmp$start[si], srep_tmp$start[sr])) - duration_overlap = end_overlap - start_overlap - duration_sib = as.numeric(srep_tmp$end[si]) - as.numeric(srep_tmp$start[si]) - perc_overlap = round(100 * (duration_overlap / duration_sib), digits = 1) - if (srep_tmp$type[sr] == "nonwear") { - srep_tmp$SIBoverlapNonwear[si] = perc_overlap - } else if (srep_tmp$type[sr] == "nap") { - srep_tmp$SIBoverlapNap[si] = perc_overlap - } - } - - if (srep_tmp$start[sr] <= srep_tmp$end[si] & - srep_tmp$end[sr] >= srep_tmp$start[si]) { - end_overlap = as.numeric(pmin(srep_tmp$end[si], srep_tmp$end[sr])) - start_overlap = as.numeric(pmax(srep_tmp$start[si], srep_tmp$start[sr])) - duration_overlap = end_overlap - start_overlap - duration_sr = as.numeric(srep_tmp$end[sr]) - as.numeric(srep_tmp$start[sr]) - perc_overlap = round(100 * (duration_overlap / duration_sr), digits = 1) - if (srep_tmp$type[sr] == "nonwear") { - srep_tmp$NonwearOverlapSIB[sr] = perc_overlap - } else if (srep_tmp$type[sr] == "nap") { - srep_tmp$NapOverlapSIB[sr] = perc_overlap - } - } - - } - } - } - } - - sibs_indices = which(srep_tmp$type == "sib") - nap_indices = which(srep_tmp$type == "nap") - nonwear_indices = which(srep_tmp$type == "nonwear") - SIBoverlapNap_indices = which(srep_tmp$SIBoverlapNap != 0) - SIBoverlapNonwear_indices = which(srep_tmp$SIBoverlapNonwear != 0) - NapOverlapSIB_indices = which(srep_tmp$NapOverlapSIB != 0) - NonwearOverlapSIB_indices = which(srep_tmp$NonwearOverlapSIB != 0) - dsummary[di,fi:(fi + 6)] = c(length(sibs_indices), - length(nap_indices), - length(nonwear_indices), - length(SIBoverlapNap_indices), - length(SIBoverlapNonwear_indices), - length(NapOverlapSIB_indices), - length(NonwearOverlapSIB_indices)) - - ds_names[fi:(fi + 6)] = c("nbouts_day_sib", "nbouts_day_srnap", "nbouts_day_srnonw", - "noverl_sib_srnap", "noverl_sib_srnonw", - "noverl_srnap_sib", "noverl_srnonw_sib") - fi = fi + 7 - if (length(sibs_indices) > 0) { - dsummary[di,fi:(fi + 1)] = c(mean(srep_tmp$duration[sibs_indices]), - sum(srep_tmp$duration[sibs_indices])) - } else { - dsummary[di,fi:(fi + 1)] = c(0, 0) - } - ds_names[fi:(fi + 1)] = c("frag_mean_dur_sib_day", "dur_day_sib_min") - fi = fi + 2 - if (length(nap_indices) > 0) { - srep_tmp$duration[nap_indices] = (as.numeric(srep_tmp$end[nap_indices]) - - as.numeric(srep_tmp$start[nap_indices])) / 60 - dsummary[di,fi:(fi + 1)] = c(mean(srep_tmp$duration[nap_indices]), - sum(srep_tmp$duration[nap_indices])) - } else { - dsummary[di,fi:(fi + 1)] = c(0, 0) - } - - ds_names[fi:(fi + 1)] = c("frag_mean_dur_srnap_day", "dur_day_srnap_min") - fi = fi + 2 - - if (length(nonwear_indices) > 0) { - dsummary[di,fi:(fi + 1)] = c(mean(srep_tmp$duration[nonwear_indices]), - sum(srep_tmp$duration[nonwear_indices])) - } else { - dsummary[di,fi:(fi + 1)] = c(0, 0) - } - ds_names[fi:(fi + 1)] = c("frag_mean_dur_srnonw_day", "dur_day_srnonw_min") - fi = fi + 2 - # Overlap sib with srnap - calcOverlapPercentage = function(overlap, duration) { - return(sum(overlap * duration) / sum(duration)) - } - if (length(SIBoverlapNap_indices) > 0) { - overlap_perc = calcOverlapPercentage(overlap = srep_tmp$SIBoverlapNap[SIBoverlapNap_indices], - duration = srep_tmp$duration[SIBoverlapNap_indices]) - dsummary[di,fi:(fi + 2)] = c(mean(srep_tmp$duration[SIBoverlapNap_indices]), - sum(srep_tmp$duration[SIBoverlapNap_indices]), - overlap_perc) - } else { - dsummary[di,fi:(fi + 2)] = c(0, 0, 0) - } - ds_names[fi:(fi + 2)] = c("mdur_sib_overl_srnap", "tdur_sib_overl_srnap", "perc_sib_overl_srnap") - fi = fi + 3 - - # Overlap srnap with sib - calcOverlapPercentage = function(overlap, duration) { - return(sum(overlap * duration) / sum(duration)) - } - if (length(SIBoverlapNap_indices) > 0) { - overlap_perc = calcOverlapPercentage(overlap = srep_tmp$NapOverlapSIB[NapOverlapSIB_indices], - duration = srep_tmp$duration[NapOverlapSIB_indices]) - dsummary[di,fi:(fi + 2)] = c(mean(srep_tmp$duration[NapOverlapSIB_indices]), - sum(srep_tmp$duration[NapOverlapSIB_indices]), - overlap_perc) - } else { - dsummary[di,fi:(fi + 2)] = c(0, 0, 0) - } - ds_names[fi:(fi + 2)] = c("mdur_srnap_overl_sib", "tdur_srnap_overl_sib", "perc_srnap_overl_sib") - fi = fi + 3 - # Overlap sib with srnonw - if (length(SIBoverlapNonwear_indices) > 0) { - calcOverlapPercentage(overlap = srep_tmp$SIBoverlapNonwear[SIBoverlapNonwear_indices], - duration = srep_tmp$duration[SIBoverlapNonwear_indices]) - dsummary[di,fi:(fi + 2)] = c(mean(srep_tmp$duration[SIBoverlapNonwear_indices]), - sum(srep_tmp$duration[SIBoverlapNonwear_indices]), - mean(srep_tmp$SIBoverlapNonwear[SIBoverlapNonwear_indices]), - overlap_perc) - } else { - dsummary[di,fi:(fi + 2)] = c(0, 0, 0) - } - ds_names[fi:(fi + 2)] = c("mdur_sib_overl_srnonw", "tdur_sib_overl_srnonw", "perc_sib_overl_srnonw") - fi = fi + 3 - # Overlap srnonw with sib - if (length(SIBoverlapNonwear_indices) > 0) { - calcOverlapPercentage(overlap = srep_tmp$NonwearOverlapSIB[NonwearOverlapSIB_indices], - duration = srep_tmp$duration[NonwearOverlapSIB_indices]) - dsummary[di,fi:(fi + 2)] = c(mean(srep_tmp$duration[NonwearOverlapSIB_indices]), - sum(srep_tmp$duration[NonwearOverlapSIB_indices]), - mean(srep_tmp$NonwearOverlapSIB[NonwearOverlapSIB_indices]), - overlap_perc) - } else { - dsummary[di,fi:(fi + 2)] = c(0, 0, 0) - } - ds_names[fi:(fi + 2)] = c("mdur_srnonw_overl_sib", "tdur_srnonw_overl_sib", "perc_srnonw_overl_sib") - fi = fi + 3 - rm(srep_tmp) - } else { - fi = fi + 27 - } - } + # sibreport = sibreport_backup + restAnalyses = g.part5.analyseRest(sibreport = sibreport, dsummary = dsummary, + ds_names = ds_names, fi = fi, di = di, + time = ts$time[sse[ts$diur[sse] == 0]], + tz = params_general[["desiredtz"]]) + fi = restAnalyses$fi + di = restAnalyses$di + dsummary = restAnalyses$dsummary + ds_names = restAnalyses$ds_names } #=============================================== # FOLDER STRUCTURE diff --git a/R/g.part5.analyseRest.R b/R/g.part5.analyseRest.R new file mode 100644 index 000000000..3dd469c9d --- /dev/null +++ b/R/g.part5.analyseRest.R @@ -0,0 +1,182 @@ +g.part5.analyseRest = function(sibreport = NULL, dsummary = NULL, + ds_names = NULL, fi = NULL, di = NULL, + time = NULL, tz = NULL) { + + longboutsi = which(sibreport$type == "sib" | (sibreport$type != "sib" & sibreport$duration >= 15)) + # for qc purposes: + dsummary[di,fi] = length(longboutsi) + ds_names[fi] = "sibreport_n_items" + fi = fi + 1 + if (length(longboutsi) > 0) { + sibreport = sibreport[longboutsi,] + srep_tmp = sibreport[which(sibreport$start > min(time) & + sibreport$end < max(time)),] + + # account for possibility that some of these categories do not exis + # identify overlapping and non-overlapping, (nap-sib, non-wear-sib, sib, nap, nonwear) + # calculate for all five categories number, total duration, mean duration + # for qc purposes: + dsummary[di,fi] = nrow(srep_tmp) + ds_names[fi] = "sibreport_n_items_day" + fi = fi + 1 + + if (nrow(srep_tmp) > 0) { + sibs = which(srep_tmp$type == "sib") + srep_tmp$SIBoverlapNonwear = 0 + srep_tmp$SIBoverlapNap = 0 + srep_tmp$NonwearOverlapSIB = 0 + srep_tmp$NapOverlapSIB = 0 + srep_tmp$start = as.POSIXlt(srep_tmp$start, tz = tz) + srep_tmp$end = as.POSIXlt(srep_tmp$end, tz = tz) + # # for qc purposes: + # dsummary[di,fi] = nrow(srep_tmp) + # ds_names[fi] = "n_sibs_sibreport" + # fi = fi + 1 + if (length(sibs) > 0) { + classes = unique(srep_tmp$type) + selfreport = which(srep_tmp$type == "nonwear" | srep_tmp$type == "nap") + if (length(selfreport) > 0) { + for (si in sibs) { + for (sr in selfreport) { + if (srep_tmp$start[si] <= srep_tmp$end[sr] & + srep_tmp$end[si] >= srep_tmp$start[sr]) { + end_overlap = as.numeric(pmin(srep_tmp$end[si], srep_tmp$end[sr])) + start_overlap = as.numeric(pmax(srep_tmp$start[si], srep_tmp$start[sr])) + duration_overlap = end_overlap - start_overlap + duration_sib = as.numeric(srep_tmp$end[si]) - as.numeric(srep_tmp$start[si]) + perc_overlap = round(100 * (duration_overlap / duration_sib), digits = 1) + if (srep_tmp$type[sr] == "nonwear") { + srep_tmp$SIBoverlapNonwear[si] = perc_overlap + } else if (srep_tmp$type[sr] == "nap") { + srep_tmp$SIBoverlapNap[si] = perc_overlap + } + } + + if (srep_tmp$start[sr] <= srep_tmp$end[si] & + srep_tmp$end[sr] >= srep_tmp$start[si]) { + end_overlap = as.numeric(pmin(srep_tmp$end[si], srep_tmp$end[sr])) + start_overlap = as.numeric(pmax(srep_tmp$start[si], srep_tmp$start[sr])) + duration_overlap = end_overlap - start_overlap + duration_sr = as.numeric(srep_tmp$end[sr]) - as.numeric(srep_tmp$start[sr]) + perc_overlap = round(100 * (duration_overlap / duration_sr), digits = 1) + if (srep_tmp$type[sr] == "nonwear") { + srep_tmp$NonwearOverlapSIB[sr] = perc_overlap + } else if (srep_tmp$type[sr] == "nap") { + srep_tmp$NapOverlapSIB[sr] = perc_overlap + } + } + + } + } + } + } + + sibs_indices = which(srep_tmp$type == "sib") + nap_indices = which(srep_tmp$type == "nap") + nonwear_indices = which(srep_tmp$type == "nonwear") + SIBoverlapNap_indices = which(srep_tmp$SIBoverlapNap != 0) + SIBoverlapNonwear_indices = which(srep_tmp$SIBoverlapNonwear != 0) + NapOverlapSIB_indices = which(srep_tmp$NapOverlapSIB != 0) + NonwearOverlapSIB_indices = which(srep_tmp$NonwearOverlapSIB != 0) + dsummary[di,fi:(fi + 6)] = c(length(sibs_indices), + length(nap_indices), + length(nonwear_indices), + length(SIBoverlapNap_indices), + length(SIBoverlapNonwear_indices), + length(NapOverlapSIB_indices), + length(NonwearOverlapSIB_indices)) + + ds_names[fi:(fi + 6)] = c("nbouts_day_sib", "nbouts_day_srnap", "nbouts_day_srnonw", + "noverl_sib_srnap", "noverl_sib_srnonw", + "noverl_srnap_sib", "noverl_srnonw_sib") + fi = fi + 7 + if (length(sibs_indices) > 0) { + dsummary[di,fi:(fi + 1)] = c(mean(srep_tmp$duration[sibs_indices]), + sum(srep_tmp$duration[sibs_indices])) + } else { + dsummary[di,fi:(fi + 1)] = c(0, 0) + } + ds_names[fi:(fi + 1)] = c("frag_mean_dur_sib_day", "dur_day_sib_min") + fi = fi + 2 + if (length(nap_indices) > 0) { + srep_tmp$duration[nap_indices] = (as.numeric(srep_tmp$end[nap_indices]) - + as.numeric(srep_tmp$start[nap_indices])) / 60 + dsummary[di,fi:(fi + 1)] = c(mean(srep_tmp$duration[nap_indices]), + sum(srep_tmp$duration[nap_indices])) + } else { + dsummary[di,fi:(fi + 1)] = c(0, 0) + } + + ds_names[fi:(fi + 1)] = c("frag_mean_dur_srnap_day", "dur_day_srnap_min") + fi = fi + 2 + + if (length(nonwear_indices) > 0) { + dsummary[di,fi:(fi + 1)] = c(mean(srep_tmp$duration[nonwear_indices]), + sum(srep_tmp$duration[nonwear_indices])) + } else { + dsummary[di,fi:(fi + 1)] = c(0, 0) + } + ds_names[fi:(fi + 1)] = c("frag_mean_dur_srnonw_day", "dur_day_srnonw_min") + fi = fi + 2 + # Overlap sib with srnap + calcOverlapPercentage = function(overlap, duration) { + return(sum(overlap * duration) / sum(duration)) + } + if (length(SIBoverlapNap_indices) > 0) { + overlap_perc = calcOverlapPercentage(overlap = srep_tmp$SIBoverlapNap[SIBoverlapNap_indices], + duration = srep_tmp$duration[SIBoverlapNap_indices]) + dsummary[di,fi:(fi + 2)] = c(mean(srep_tmp$duration[SIBoverlapNap_indices]), + sum(srep_tmp$duration[SIBoverlapNap_indices]), + overlap_perc) + } else { + dsummary[di,fi:(fi + 2)] = c(0, 0, 0) + } + ds_names[fi:(fi + 2)] = c("mdur_sib_overl_srnap", "tdur_sib_overl_srnap", "perc_sib_overl_srnap") + fi = fi + 3 + + # Overlap srnap with sib + calcOverlapPercentage = function(overlap, duration) { + return(sum(overlap * duration) / sum(duration)) + } + if (length(SIBoverlapNap_indices) > 0) { + overlap_perc = calcOverlapPercentage(overlap = srep_tmp$NapOverlapSIB[NapOverlapSIB_indices], + duration = srep_tmp$duration[NapOverlapSIB_indices]) + dsummary[di,fi:(fi + 2)] = c(mean(srep_tmp$duration[NapOverlapSIB_indices]), + sum(srep_tmp$duration[NapOverlapSIB_indices]), + overlap_perc) + } else { + dsummary[di,fi:(fi + 2)] = c(0, 0, 0) + } + ds_names[fi:(fi + 2)] = c("mdur_srnap_overl_sib", "tdur_srnap_overl_sib", "perc_srnap_overl_sib") + fi = fi + 3 + # Overlap sib with srnonw + if (length(SIBoverlapNonwear_indices) > 0) { + overlap_perc = calcOverlapPercentage(overlap = srep_tmp$SIBoverlapNonwear[SIBoverlapNonwear_indices], + duration = srep_tmp$duration[SIBoverlapNonwear_indices]) + dsummary[di,fi:(fi + 2)] = c(mean(srep_tmp$duration[SIBoverlapNonwear_indices]), + sum(srep_tmp$duration[SIBoverlapNonwear_indices]), + overlap_perc) + } else { + dsummary[di,fi:(fi + 2)] = c(0, 0, 0) + } + ds_names[fi:(fi + 2)] = c("mdur_sib_overl_srnonw", "tdur_sib_overl_srnonw", "perc_sib_overl_srnonw") + fi = fi + 3 + # Overlap srnonw with sib + if (length(SIBoverlapNonwear_indices) > 0) { + overlap_perc = calcOverlapPercentage(overlap = srep_tmp$NonwearOverlapSIB[NonwearOverlapSIB_indices], + duration = srep_tmp$duration[NonwearOverlapSIB_indices]) + dsummary[di,fi:(fi + 2)] = c(mean(srep_tmp$duration[NonwearOverlapSIB_indices]), + sum(srep_tmp$duration[NonwearOverlapSIB_indices]), + overlap_perc) + } else { + dsummary[di,fi:(fi + 2)] = c(0, 0, 0) + } + ds_names[fi:(fi + 2)] = c("mdur_srnonw_overl_sib", "tdur_srnonw_overl_sib", "perc_srnonw_overl_sib") + fi = fi + 3 + rm(srep_tmp) + } else { + fi = fi + 26 + } + } + invisible(list(fi = fi, di = di, ds_names = ds_names, dsummary = dsummary)) +} diff --git a/man/g.part5.analyseRest.Rd b/man/g.part5.analyseRest.Rd new file mode 100644 index 000000000..601dcea4c --- /dev/null +++ b/man/g.part5.analyseRest.Rd @@ -0,0 +1,44 @@ +\name{g.part5.analyseRest} +\alias{g.part5.analyseRest} +\title{ + Analyse rest (internal function) +} +\description{ + Analyses overlap self-reported napping, non-wear and sib. + Internal function not intended for direct use by user. +} +\usage{ + g.part5.analyseRest(sibreport = NULL, dsummary = NULL, + ds_names = NULL, fi = NULL, + di = NULL, time = NULL, tz = NULL) +} +\arguments{ + \item{sibreport}{ + sibreport data.frame produced by \link{g.sibreport} + } + \item{dsummary}{ + matrix created internally by \link{g.part5} + } + \item{ds_names}{ + character vector with variable names corresponding to dsummary + created internally by \link{g.part5} + } + \item{fi}{ + Numeric scalar to indicate variable index, created internally by \link{g.part5} + } + \item{di}{ + Numeric scalar to indicate recording index, created internally by \link{g.part5} + } + \item{time}{ + Daytime section of time column from the ts object, created internally by \link{g.part5}, + } + \item{tz}{ + Timezone database name + } +} +\value{ + List with updated objects dsummary, ds_names, fi, and di +} +\author{ + Vincent T van Hees +} \ No newline at end of file diff --git a/tests/testthat/test_chainof5parts.R b/tests/testthat/test_chainof5parts.R index 90f937036..af00e8d88 100644 --- a/tests/testthat/test_chainof5parts.R +++ b/tests/testthat/test_chainof5parts.R @@ -134,7 +134,7 @@ test_that("chainof5parts", { expect_true(dir.exists(dirname)) expect_true(file.exists(rn[1])) - expect_that(nrow(output),equals(3)) # changed because part5 now gives also first and last day + expect_that(nrow(output), equals(3)) # changed because part5 now gives also first and last day expect_that(ncol(output),equals(181)) expect_that(round(as.numeric(output$wakeup[2]), digits = 4), equals(36)) dirname_raw = "output_test/meta/ms5.outraw/40_100_400" diff --git a/tests/testthat/test_part5_analyseRest.R b/tests/testthat/test_part5_analyseRest.R new file mode 100644 index 000000000..1115e75d8 --- /dev/null +++ b/tests/testthat/test_part5_analyseRest.R @@ -0,0 +1,100 @@ +library(GGIR) +context("g.part5.analyseRest") +tz = "Europe/Amsterdam" +ds_names = rep("", 40) +dsummary = matrix("", 1, 40) +startday = as.POSIXlt(x = "2022-06-02 08:00:00", tz = tz) +time = seq(startday, startday + (16 * 3600), by = 60) + + + +test_that("Overlap 1 nap and 1 sib", { + fi = 1 + di = 1 + sibreport = data.frame(ID = rep("test123", 2), type = c("nap", "sib"), + start = c("2022-06-02 14:00:00", "2022-06-02 14:05:00"), + end = c("2022-06-02 14:20:00", "2022-06-02 14:20:00")) + sibreport$duration = as.numeric(difftime(time1 = sibreport$end, time2 = sibreport$start, units = "mins", tz = tz)) + restAnalyses = g.part5.analyseRest(sibreport = sibreport, dsummary = dsummary, + ds_names = ds_names, fi = fi, di = di, + time = time, + tz = tz) + fi = restAnalyses$fi + di = restAnalyses$di + dsummary = restAnalyses$dsummary + ds_names = restAnalyses$ds_names + dsummary = as.numeric(dsummary[, which(ds_names != "")]) + ds_names = ds_names[ds_names != ""] + names(dsummary) = ds_names + dsummary = as.data.frame(t(dsummary)) + + expect_equal(dsummary$nbouts_day_sib, 1) + expect_equal(dsummary$nbouts_day_srnap, 1) + expect_equal(dsummary$frag_mean_dur_sib_day, 15) + expect_equal(dsummary$frag_mean_dur_srnap_day, 20) + expect_equal(dsummary$perc_sib_overl_srnap, 100) + expect_equal(dsummary$perc_srnap_overl_sib, 75) + expect_equal(sum(dsummary), 323) +}) + +test_that("Overlap 1 nonwear and 1 sib", { + fi = 1 + di = 1 + sibreport = data.frame(ID = rep("test123", 2), type = c("nonwear", "sib"), + start = c("2022-06-02 14:00:00", "2022-06-02 14:05:00"), + end = c("2022-06-02 14:20:00", "2022-06-02 14:20:00")) + sibreport$duration = as.numeric(difftime(time1 = sibreport$end, time2 = sibreport$start, units = "mins", tz = tz)) + restAnalyses = g.part5.analyseRest(sibreport = sibreport, dsummary = dsummary, + ds_names = ds_names, fi = fi, di = di, + time = time, + tz = tz) + fi = restAnalyses$fi + di = restAnalyses$di + dsummary = restAnalyses$dsummary + ds_names = restAnalyses$ds_names + dsummary = as.numeric(dsummary[, which(ds_names != "")]) + ds_names = ds_names[ds_names != ""] + names(dsummary) = ds_names + dsummary = as.data.frame(t(dsummary)) + + expect_equal(dsummary$nbouts_day_sib, 1) + expect_equal(dsummary$nbouts_day_srnonw, 1) + expect_equal(dsummary$frag_mean_dur_sib_day, 15) + expect_equal(dsummary$frag_mean_dur_srnonw_day, 20) + expect_equal(dsummary$perc_sib_overl_srnonw, 100) + expect_equal(dsummary$perc_srnonw_overl_sib, 75) + expect_equal(sum(dsummary), 323) +}) + + +test_that("No overlap 1 nonwear, 1 nap, and 1 sib", { + fi = 1 + di = 1 + sibreport = data.frame(ID = rep("test123", 3), type = c("nonwear", "nap", "sib"), + start = c("2022-06-02 12:00:00", "2022-06-02 13:00:00", "2022-06-02 15:00:00"), + end = c("2022-06-02 12:20:00", "2022-06-02 13:20:00", "2022-06-02 15:20:00")) + sibreport$duration = as.numeric(difftime(time1 = sibreport$end, time2 = sibreport$start, units = "mins", tz = tz)) + restAnalyses = g.part5.analyseRest(sibreport = sibreport, dsummary = dsummary, + ds_names = ds_names, fi = fi, di = di, + time = time, + tz = tz) + fi = restAnalyses$fi + di = restAnalyses$di + dsummary = restAnalyses$dsummary + ds_names = restAnalyses$ds_names + dsummary = as.numeric(dsummary[, which(ds_names != "")]) + ds_names = ds_names[ds_names != ""] + names(dsummary) = ds_names + dsummary = as.data.frame(t(dsummary)) + + expect_equal(dsummary$nbouts_day_sib, 1) + expect_equal(dsummary$nbouts_day_srnap, 1) + expect_equal(dsummary$nbouts_day_srnonw, 1) + expect_equal(dsummary$frag_mean_dur_sib_day, 20) + expect_equal(dsummary$frag_mean_dur_srnap_day, 20) + expect_equal(dsummary$frag_mean_dur_srnonw_day, 20) + expect_equal(dsummary$perc_sib_overl_srnap, 0) + expect_equal(dsummary$perc_sib_overl_srnonw, 0) + expect_equal(dsummary$perc_srnonw_overl_sib, 0) + expect_equal(sum(dsummary), 129) +}) From b6e61da4997c2f0b55e7c54691383e3ada288327 Mon Sep 17 00:00:00 2001 From: Vincent van Hees Date: Fri, 1 Sep 2023 08:29:20 +0200 Subject: [PATCH 07/18] Moving rest analyses to new g.part5_analyseSegment function #687 --- R/g.part5.R | 17 ++--------------- R/g.part5_analyseSegment.R | 15 ++++++++++++++- 2 files changed, 16 insertions(+), 16 deletions(-) diff --git a/R/g.part5.R b/R/g.part5.R index de4a4c752..e186635f3 100644 --- a/R/g.part5.R +++ b/R/g.part5.R @@ -503,7 +503,8 @@ g.part5 = function(datadir = c(), metadatadir = c(), f0=c(), f1=c(), fullFilename = fullfilenames[i], add_one_day_to_next_date, lightpeak_available, tail_expansion_log, - foldernamei = foldername[i]) + foldernamei = foldername[i], + sibreport = sibreport) # Extract essential object to be used as input for the next # segment indexlog = gas$indexlog @@ -523,20 +524,6 @@ g.part5 = function(datadir = c(), metadatadir = c(), f0=c(), f1=c(), ws3new = timeList$epochSize if (doNext == TRUE) next } - #======================================================= - # nap/sib/nonwear overlap analysis - #======================================================= - if (params_output[["do.sibreport"]] == TRUE) { - # sibreport = sibreport_backup - restAnalyses = g.part5.analyseRest(sibreport = sibreport, dsummary = dsummary, - ds_names = ds_names, fi = fi, di = di, - time = ts$time[sse[ts$diur[sse] == 0]], - tz = params_general[["desiredtz"]]) - fi = restAnalyses$fi - di = restAnalyses$di - dsummary = restAnalyses$dsummary - ds_names = restAnalyses$ds_names - } #=============================================== # FOLDER STRUCTURE if (params_output[["storefolderstructure"]] == TRUE) { diff --git a/R/g.part5_analyseSegment.R b/R/g.part5_analyseSegment.R index faac20ade..18493f8af 100644 --- a/R/g.part5_analyseSegment.R +++ b/R/g.part5_analyseSegment.R @@ -10,7 +10,7 @@ g.part5_analyseSegment = function(indexlog, timeList, levelList, add_one_day_to_next_date, lightpeak_available, tail_expansion_log, - foldernamei) { + foldernamei, sibreport) { # unpack list objects: # indexlog fileIndex = indexlog$fileIndex @@ -449,6 +449,19 @@ g.part5_analyseSegment = function(indexlog, timeList, levelList, fi = fi + length(luxperseg$values) } } + #======================================================= + # nap/sib/nonwear overlap analysis + #======================================================= + if (params_output[["do.sibreport"]] == TRUE) { + restAnalyses = g.part5.analyseRest(sibreport = sibreport, dsummary = dsummary, + ds_names = ds_names, fi = fi, di = si, + time = ts$time[sse[ts$diur[sse] == 0]], + tz = params_general[["desiredtz"]]) + fi = restAnalyses$fi + si = restAnalyses$di + dsummary = restAnalyses$dsummary + ds_names = restAnalyses$ds_names + } #=============================================== # FOLDER STRUCTURE if (params_output[["storefolderstructure"]] == TRUE) { From 9158f57c3966627d52ead3338425ea5bc3c8bb76 Mon Sep 17 00:00:00 2001 From: Vincent van Hees Date: Fri, 1 Sep 2023 08:30:11 +0200 Subject: [PATCH 08/18] update documentation for previous commit --- man/g.part5_analyseSegment.Rd | 3 +++ 1 file changed, 3 insertions(+) diff --git a/man/g.part5_analyseSegment.Rd b/man/g.part5_analyseSegment.Rd index cef196d71..71bef6f00 100644 --- a/man/g.part5_analyseSegment.Rd +++ b/man/g.part5_analyseSegment.Rd @@ -83,6 +83,9 @@ \item{foldernamei}{ Character with folder name in which the data file is stored. } + \item{sibreport}{ + Sibreport object as passed on from \link{g.part5} + } } \value{ List with objects: indexlog, timeList, and the matrix with the prelimenary results From b235928b9be911e23f823a444458ceae7b7fa727 Mon Sep 17 00:00:00 2001 From: Vincent van Hees Date: Fri, 1 Sep 2023 09:09:44 +0200 Subject: [PATCH 09/18] make indentation align with rest of repo --- R/g.part5.R | 14 +++++------ man/g.part5_analyseSegment.Rd | 2 +- tests/testthat/test_chainof5parts.R | 36 ++++++++++++++--------------- 3 files changed, 26 insertions(+), 26 deletions(-) diff --git a/R/g.part5.R b/R/g.part5.R index e186635f3..0e7faf6a6 100644 --- a/R/g.part5.R +++ b/R/g.part5.R @@ -63,7 +63,7 @@ g.part5 = function(datadir = c(), metadatadir = c(), f0=c(), f1=c(), #====================================================================== # compile lists of milestone data filenames fnames.ms3 = dir(paste(metadatadir, "/meta/ms3.out", sep = "")) - + fnames.ms5 = dir(paste(metadatadir, "/meta/ms5.out", sep = "")) # path to sleeplog milestonedata, if it exists: sleeplogRDA = paste(metadatadir, "/meta/sleeplog.RData", sep = "") @@ -225,9 +225,9 @@ g.part5 = function(datadir = c(), metadatadir = c(), f0=c(), f1=c(), S$sib.end.time[replaceLastWakeup] = ts$time[nrow(ts)] } } - + for (sibDef in def) { # loop through sleep definitions (defined by angle and time threshold in g.part3) - + ws3new = ws3 # reset wse3new, because if part5_agg2_60seconds is TRUE then this will have been change in the previous iteration of the loop if (params_general[["part5_agg2_60seconds"]] == TRUE) { ts = ts_backup @@ -277,7 +277,7 @@ g.part5 = function(datadir = c(), metadatadir = c(), f0=c(), f1=c(), Nepochsinhour, SPTE_end) 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 - + # only include angle if angle is present angleColName = ifelse("angle" %in% names(ts), yes = "angle", no = NULL) if (lightpeak_available == TRUE) { @@ -329,7 +329,7 @@ g.part5 = function(datadir = c(), metadatadir = c(), f0=c(), f1=c(), data.table::fwrite(x = sibreport, file = sibreport_fname, row.names = FALSE, sep = params_output[["sep_reports"]]) # nap/sib/nonwear overlap analysis - + if (length(params_sleep[["nap_model"]]) > 0) { # nap detection if (params_general[["acc.metric"]] != "ENMO" | @@ -394,7 +394,7 @@ g.part5 = function(datadir = c(), metadatadir = c(), f0=c(), f1=c(), for (TRVi in params_phyact[["threshold.vig"]]) { # derive behavioral levels (class), e.g. MVPA, inactivity bouts, etc. levelList = identify_levels(ts = ts, TRLi = TRLi, TRMi = TRMi, TRVi = TRVi, - ws3 = ws3new, params_phyact = params_phyact) + ws3 = ws3new, params_phyact = params_phyact) LEVELS = levelList$LEVELS OLEVELS = levelList$OLEVELS Lnames = levelList$Lnames @@ -583,7 +583,7 @@ g.part5 = function(datadir = c(), metadatadir = c(), f0=c(), f1=c(), } output = data.frame(dsummary,stringsAsFactors = FALSE) names(output) = ds_names - + # correct definition of sleep log availability for window = WW, because now it # also relies on sleep log from previous night whoareWW = which(output$window == "WW") # look up WW diff --git a/man/g.part5_analyseSegment.Rd b/man/g.part5_analyseSegment.Rd index 71bef6f00..d98fdee69 100644 --- a/man/g.part5_analyseSegment.Rd +++ b/man/g.part5_analyseSegment.Rd @@ -19,7 +19,7 @@ add_one_day_to_next_date, lightpeak_available, tail_expansion_log, - foldernamei) + foldernamei, sibreport) } \arguments{ \item{indexlog}{ diff --git a/tests/testthat/test_chainof5parts.R b/tests/testthat/test_chainof5parts.R index 6d267abe1..73d25d312 100644 --- a/tests/testthat/test_chainof5parts.R +++ b/tests/testthat/test_chainof5parts.R @@ -2,7 +2,7 @@ library(GGIR) context("Chainof5parts") test_that("chainof5parts", { skip_on_cran() - + Ndays = 2 create_test_acc_csv(Nmin = Ndays*1440) create_test_sleeplog_csv(advanced = FALSE) @@ -55,7 +55,7 @@ test_that("chainof5parts", { expect_equal(round(as.numeric(SUM$summary$meas_dur_def_proto_day), digits = 3), 1) expect_equal(SUM$summary$`N valid WEdays`, "1") expect_equal(SUM$summary$`N valid WKdays`, "2") - + # part 2 with strategy = 5 g.part2(datadir = fn, metadatadir = metadatadir, f0 = 1, f1 = 1, idloc = 2, desiredtz = desiredtz, ndayswindow = 1, @@ -69,7 +69,7 @@ test_that("chainof5parts", { expect_equal(round(as.numeric(SUM$summary$meas_dur_def_proto_day), digits = 3), 1) expect_equal(SUM$summary$`N valid WEdays`, "1") expect_equal(SUM$summary$`N valid WKdays`, "2") - + # part 2 with strategy = 2 and iglevels = 1 g.part2(datadir = fn, metadatadir = metadatadir, f0 = 1, f1 = 1, idloc = 2, desiredtz = desiredtz, @@ -117,10 +117,10 @@ test_that("chainof5parts", { dirname = "output_test/meta/ms3.out/" rn3 = dir(dirname,full.names = TRUE) load(rn3[1]) - + expect_true(dir.exists(dirname)) expect_that(round(sum(sib.cla.sum[,4:7]), digits = 0), equals(2957)) - + #-------------------------------------------- # part 4 g.part4(datadir = fn, metadatadir = metadatadir, f0 = 1, f1 = 1, @@ -139,7 +139,7 @@ test_that("chainof5parts", { expect_that(round(nightsummary$number_sib_wakinghours[1], digits = 4), equals(6)) expect_true(as.logical(nightsummary$acc_available[1])) expect_true(as.logical(nightsummary$sleeplog_used[1])) - + #-------------------------------------------- #part 5 g.part5(datadir = fn, metadatadir = metadatadir, f0 = 1, f1 = 1, desiredtz = desiredtz, @@ -152,13 +152,13 @@ test_that("chainof5parts", { 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"))) - + dirname = "output_test/meta/ms5.out/" rn = dir(dirname,full.names = TRUE) load(rn[1]) - + dirname = "output_test/meta/ms5.out/" - + expect_true(dir.exists(dirname)) expect_true(file.exists(rn[1])) expect_that(nrow(output),equals(3)) # changed because part5 now gives also first and last day @@ -188,7 +188,7 @@ test_that("chainof5parts", { 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")) dn = "output_test" - + #======================= # Different variations on part 4: #-------------------------------------------- @@ -199,7 +199,7 @@ test_that("chainof5parts", { nnights = 7, colid = 1, coln1 = 2, relyonguider = FALSE, desiredtz = desiredtz, storefolderstructure = TRUE, overwrite = TRUE, verbose = FALSE) - + dirname = "output_test/meta/ms4.out/" rn = dir(dirname,full.names = TRUE) load(rn[1]) @@ -207,7 +207,7 @@ test_that("chainof5parts", { expect_that(round(nightsummary$SptDuration[1], digits = 4), equals(18.075)) expect_true(as.logical(nightsummary$acc_available[1])) expect_false(as.logical(nightsummary$sleeplog_used[1])) - + #---------------------------------------------------------------------- # Part 4 - daysleeper scenario by modifying the part 3 output & criterror = 0, and relyonguider=TRUE SPTE_end = c(37, 37) # turn into midday @@ -230,7 +230,7 @@ test_that("chainof5parts", { rn = dir(dirname,full.names = TRUE) load(rn[1]) vis_sleep_file = "output_test/results/visualisation_sleep.pdf" - + expect_true(dir.exists(dirname)) expect_true(file.exists(vis_sleep_file)) expect_that(round(nightsummary$number_sib_wakinghours[1], digits = 4), equals(0)) @@ -263,7 +263,7 @@ test_that("chainof5parts", { expect_true(dir.exists(dirname)) expect_that(round(nightsummary$number_sib_wakinghours[1], digits = 4), equals(10)) expect_that(round(nightsummary$SptDuration[1], digits = 4), equals(7)) - + #---------------------------------------------------------------------- # Part 4 - DST-1 load(rn3[1]) @@ -316,7 +316,7 @@ test_that("chainof5parts", { colnames(SDF) = c("Monitor", "Day1", "Day2") selectdaysfile = "selectdaysfile.csv" write.csv(SDF, file = selectdaysfile) - + g.part1(datadir = fn, metadatadir = metadatadir, f0 = 1, f1 = 1, overwrite = TRUE, desiredtz = desiredtz, @@ -329,7 +329,7 @@ test_that("chainof5parts", { do.bfx = TRUE, do.bfy = TRUE, do.bfz = TRUE, do.hfen = TRUE, do.hfx = TRUE, do.hfy = TRUE, do.hfz = TRUE, do.lfen = TRUE, do.enmoa = TRUE, selectdaysfile = selectdaysfile, verbose = FALSE) - + rn = dir("output_test/meta/basic/", full.names = TRUE) load(rn[1]) expect_equal(ncol(M$metashort), 24) @@ -341,8 +341,8 @@ test_that("chainof5parts", { expect_equal(mean(M$metashort$roll_med_acc_z, na.rm = T), 0.007, tolerance = 3) expect_equal(mean(M$metashort$dev_roll_med_acc_x, na.rm = T), 0.007, tolerance = 3) expect_equal(mean(M$metashort$ENMOa, na.rm = T), 0.03, tolerance = 3) - - + + if (file.exists(selectdaysfile)) file.remove(selectdaysfile) if (dir.exists(dn)) unlink(dn, recursive = TRUE) if (file.exists(fn)) unlink(fn) From 129c5bca47860f95deae252e51a1e85be8db5690 Mon Sep 17 00:00:00 2001 From: Vincent van Hees Date: Tue, 12 Sep 2023 15:30:04 +0200 Subject: [PATCH 10/18] making sure rest analysis is not run when no sibreport object is available #687 and minor typo fix --- R/g.part5.R | 2 ++ R/g.part5_analyseSegment.R | 4 ++-- man/g.part5.analyseRest.Rd | 2 +- man/g.part5_analyseSegment.Rd | 2 +- 4 files changed, 6 insertions(+), 4 deletions(-) diff --git a/R/g.part5.R b/R/g.part5.R index 0e7faf6a6..7fac46c96 100644 --- a/R/g.part5.R +++ b/R/g.part5.R @@ -384,6 +384,8 @@ g.part5 = function(datadir = c(), metadatadir = c(), f0=c(), f1=c(), } } } + } else { + sibreport = NULL } ts$window = 0 # backup of nightsi outside threshold defintions to avoid diff --git a/R/g.part5_analyseSegment.R b/R/g.part5_analyseSegment.R index 18493f8af..f010d4771 100644 --- a/R/g.part5_analyseSegment.R +++ b/R/g.part5_analyseSegment.R @@ -10,7 +10,7 @@ g.part5_analyseSegment = function(indexlog, timeList, levelList, add_one_day_to_next_date, lightpeak_available, tail_expansion_log, - foldernamei, sibreport) { + foldernamei, sibreport = NULL) { # unpack list objects: # indexlog fileIndex = indexlog$fileIndex @@ -452,7 +452,7 @@ g.part5_analyseSegment = function(indexlog, timeList, levelList, #======================================================= # nap/sib/nonwear overlap analysis #======================================================= - if (params_output[["do.sibreport"]] == TRUE) { + if (params_output[["do.sibreport"]] == TRUE & !is.null(sibreport)) { restAnalyses = g.part5.analyseRest(sibreport = sibreport, dsummary = dsummary, ds_names = ds_names, fi = fi, di = si, time = ts$time[sse[ts$diur[sse] == 0]], diff --git a/man/g.part5.analyseRest.Rd b/man/g.part5.analyseRest.Rd index 601dcea4c..e6c879e77 100644 --- a/man/g.part5.analyseRest.Rd +++ b/man/g.part5.analyseRest.Rd @@ -5,7 +5,7 @@ } \description{ Analyses overlap self-reported napping, non-wear and sib. - Internal function not intended for direct use by user. + Internal function not intended for direct use by GGIR end-user. } \usage{ g.part5.analyseRest(sibreport = NULL, dsummary = NULL, diff --git a/man/g.part5_analyseSegment.Rd b/man/g.part5_analyseSegment.Rd index d98fdee69..6a24a92af 100644 --- a/man/g.part5_analyseSegment.Rd +++ b/man/g.part5_analyseSegment.Rd @@ -19,7 +19,7 @@ add_one_day_to_next_date, lightpeak_available, tail_expansion_log, - foldernamei, sibreport) + foldernamei, sibreport = NULL) } \arguments{ \item{indexlog}{ From 4ddaebec00ca637cbc1e70a8c30d588eeb905950 Mon Sep 17 00:00:00 2001 From: Vincent van Hees Date: Tue, 12 Sep 2023 16:07:40 +0200 Subject: [PATCH 11/18] Improving documentation for #687 and updating changelog --- NEWS.md | 3 +++ man/GGIR.Rd | 4 +++- vignettes/GGIR.Rmd | 13 ++++++++++++- 3 files changed, 18 insertions(+), 2 deletions(-) diff --git a/NEWS.md b/NEWS.md index eef8f677d..222ac498f 100644 --- a/NEWS.md +++ b/NEWS.md @@ -1,5 +1,8 @@ # CHANGES IN GGIR VERSION 2.10-3 +- Part 5: Performance an analysis of the overlap between self-reported naps and nonwear +with sustained inactivity bouts. When Fixes #687 + - Visualreport: Bug fixed extracting the numeric value of the days to exclude in g.plot5 #879 - Part 1: Fixed minor bug to keep calibration data in the data quality report after re-running parts 1 and 2 with overwrite = TRUE and backup.cal.coef = "retrieve" #896 diff --git a/man/GGIR.Rd b/man/GGIR.Rd index 8ee4a5f05..86698da3b 100755 --- a/man/GGIR.Rd +++ b/man/GGIR.Rd @@ -1291,7 +1291,9 @@ GGIR(mode = 1:5, \item{do.sibreport}{ Boolean (default = FALSE). In \link{g.part4}: To indicate whether to generate report for the sustained - inactivity bouts (SIB).} + inactivity bouts (SIB). If set to TRUE and when an advanced sleep diary is + available in part 4 then part 5 will use this to generate summary statistics + on the overlap between self-reported nonwear and napping with SIB} \item{do.visual}{ Boolean (default = TRUE). diff --git a/vignettes/GGIR.Rmd b/vignettes/GGIR.Rmd index 15a213d22..a9394d6d7 100644 --- a/vignettes/GGIR.Rmd +++ b/vignettes/GGIR.Rmd @@ -438,7 +438,8 @@ comes with the following changes: - You can add columns relating to self-reported napping time and nonwear time. These are not used for the sleep analysis in g.part3 and g.part4, but used in g.part5 to facilitate napping analysis, see argument - `do.sibreport`. Multiple naps and multiple nonwear periods can be + `do.sibreport` and the paragraph on naps. Multiple naps and multiple + nonwear periods can be entered per day. - Leave cells for missing values blank. - Column names are critical for the advanced sleeplog format: Date @@ -975,6 +976,16 @@ vigorous activity: | `W0.5` | Derived from power law exponent alpha, see [Chastin et al. 201 0](https://www.sciencedirect.com/science/article/abs/pii/S096663620900602X?via%3Dihub) | | nap_count | Total number of naps, only calculated when argument do.sibreport = TRUE, currently optimised for 3.5-year olds. See function documentation for function `g.part5.classifyNaps` in the [GGIR function documentation (pdf)](https:/CRAN.R-project.org/package=GGIR/GGIR.pdf). | | nap_totalduration | Total nap duration, only calculated when argument do.sibreport = TRUE, currently optimised for 3.5-year old. See function documentation for function `g.part5.classifyNaps` in the [GGIR function documentation (pdf)](https:/CRAN.R-project.org/package=GGIR/GGIR.pdf). | +| sibreport_n_items | Only created if `do.sibreport = TRUE`. Number of items in the sibreport | +| sibreport_n_items_day | Only created if `do.sibreport = TRUE`. Number of items in the sibreport for this specific day | +| nbouts_day_X | Only created if `do.sibreport = TRUE`. Number of bouts in a day of X where X can be sib (sustained inactivity bout), srnap (self-reported nap) or srnonw (self-reported nonwear) | +| noverl_X | Only created if `do.sibreport = TRUE`. Number of overlapping bouts in a day of X where X can be sib_srnap, sib_srnonw, srnap_sib, or srnonw_sib | +| frag_mean_dur_X_day | Only created if `do.sibreport = TRUE`. Mean duration of X per day, where X can be sib, srnap or srnonw | +| dur_day_X_min | Only created if `do.sibreport = TRUE`. Total duration in minutes of X per day, where X can be sib, srnap or srnonw | +| mdur_X_overl_Y | Only created if `do.sibreport = TRUE`. Mean duration of the overlap between X and Y, which are combinations of be sib, srnap or srnonw | +| tdur_X_overl_Y | Only created if `do.sibreport = TRUE`. Total duration in minutes of the overlap between X and Y, which are combinations of be sib, srnap or srnonw | +| perc_X_overl_Y | Only created if `do.sibreport = TRUE`. Percentage of overlap between X and Y, which are combinations of sib, srnap or srnonw | + **Special note if you are working on compositional data analysis:** From 6bfcc41a8b1cea516687bfba6a170cd032e46e05 Mon Sep 17 00:00:00 2001 From: Vincent van Hees Date: Tue, 12 Sep 2023 16:46:51 +0200 Subject: [PATCH 12/18] making sure code can run when nap_model is not specified --- R/g.part5.R | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) diff --git a/R/g.part5.R b/R/g.part5.R index 7fac46c96..ba06e465b 100644 --- a/R/g.part5.R +++ b/R/g.part5.R @@ -315,7 +315,7 @@ g.part5 = function(datadir = c(), metadatadir = c(), f0=c(), f1=c(), #=============================================== # Use sib.report to classify naps, non-wear and integrate these in time series # does not depend on bout detection criteria or window definitions. - if (params_output[["do.sibreport"]] == TRUE & length(params_sleep[["nap_model"]]) > 0) { + if (params_output[["do.sibreport"]] == TRUE) { IDtmp = as.character(ID) sibreport = g.sibreport(ts, ID = IDtmp, epochlength = ws3new, logs_diaries, desiredtz = params_general[["desiredtz"]]) @@ -329,7 +329,6 @@ g.part5 = function(datadir = c(), metadatadir = c(), f0=c(), f1=c(), data.table::fwrite(x = sibreport, file = sibreport_fname, row.names = FALSE, sep = params_output[["sep_reports"]]) # nap/sib/nonwear overlap analysis - if (length(params_sleep[["nap_model"]]) > 0) { # nap detection if (params_general[["acc.metric"]] != "ENMO" | From fef0bd0a01390b1608d909274b296ba7151da0df Mon Sep 17 00:00:00 2001 From: Vincent van Hees Date: Tue, 12 Sep 2023 16:47:45 +0200 Subject: [PATCH 13/18] reusing argument possible_nap_duration to specify minimum sib duration here. --- R/g.part5.analyseRest.R | 7 ++++--- R/g.part5_analyseSegment.R | 3 ++- man/g.part5_analyseSegment.Rd | 7 ++++++- 3 files changed, 12 insertions(+), 5 deletions(-) diff --git a/R/g.part5.analyseRest.R b/R/g.part5.analyseRest.R index 3dd469c9d..026d06bdf 100644 --- a/R/g.part5.analyseRest.R +++ b/R/g.part5.analyseRest.R @@ -1,8 +1,9 @@ g.part5.analyseRest = function(sibreport = NULL, dsummary = NULL, ds_names = NULL, fi = NULL, di = NULL, - time = NULL, tz = NULL) { - - longboutsi = which(sibreport$type == "sib" | (sibreport$type != "sib" & sibreport$duration >= 15)) + time = NULL, tz = NULL, possible_nap_dur = NULL) { + # Only consider sib episodes with minimum duration + longboutsi = which((sibreport$type == "sib" & sibreport$duration >= possible_nap_dur) | + (sibreport$type != "sib" & sibreport$duration >= 1)) # for qc purposes: dsummary[di,fi] = length(longboutsi) ds_names[fi] = "sibreport_n_items" diff --git a/R/g.part5_analyseSegment.R b/R/g.part5_analyseSegment.R index f010d4771..cf3cce9a0 100644 --- a/R/g.part5_analyseSegment.R +++ b/R/g.part5_analyseSegment.R @@ -456,7 +456,8 @@ g.part5_analyseSegment = function(indexlog, timeList, levelList, restAnalyses = g.part5.analyseRest(sibreport = sibreport, dsummary = dsummary, ds_names = ds_names, fi = fi, di = si, time = ts$time[sse[ts$diur[sse] == 0]], - tz = params_general[["desiredtz"]]) + tz = params_general[["desiredtz"]], + possible_nap_dur = params_sleep[["possible_nap_dur"]]) fi = restAnalyses$fi si = restAnalyses$di dsummary = restAnalyses$dsummary diff --git a/man/g.part5_analyseSegment.Rd b/man/g.part5_analyseSegment.Rd index 6a24a92af..c73d3d055 100644 --- a/man/g.part5_analyseSegment.Rd +++ b/man/g.part5_analyseSegment.Rd @@ -19,7 +19,8 @@ add_one_day_to_next_date, lightpeak_available, tail_expansion_log, - foldernamei, sibreport = NULL) + foldernamei, sibreport = NULL, + possible_nap_dur = NULL) } \arguments{ \item{indexlog}{ @@ -86,6 +87,10 @@ \item{sibreport}{ Sibreport object as passed on from \link{g.part5} } + \item{possible_nap_dur}{ + Minimum sib duration to be considered. For self-reported naps/nonwear + there is a minimum duration of 1 epoch. + } } \value{ List with objects: indexlog, timeList, and the matrix with the prelimenary results From 848ebd5bae47758b928d07960cca575a80a4c649 Mon Sep 17 00:00:00 2001 From: Vincent van Hees Date: Tue, 12 Sep 2023 17:24:34 +0200 Subject: [PATCH 14/18] Adding filter options for sibreport anlaysis #687 --- R/check_params.R | 2 +- R/g.part5.analyseRest.R | 12 ++++++++++-- R/load_params.R | 3 ++- man/GGIR.Rd | 9 +++++++-- man/g.part5.analyseRest.Rd | 11 ++++++++++- man/g.part5_analyseSegment.Rd | 7 +------ tests/testthat/test_load_check_params.R | 2 +- 7 files changed, 32 insertions(+), 14 deletions(-) diff --git a/R/check_params.R b/R/check_params.R index af67978f3..04a9f3e25 100644 --- a/R/check_params.R +++ b/R/check_params.R @@ -30,7 +30,7 @@ check_params = function(params_sleep = c(), params_metrics = c(), #----------------------------------------------------------------------------------------- if (length(params_sleep) > 0) { # Check class of sleep parameters numeric_params = c("anglethreshold", "timethreshold", "longitudinal_axis", "possible_nap_window", "possible_nap_dur", - "colid", "coln1", "def.noc.sleep", "nnights", "sleepefficiency.metric") + "colid", "coln1", "def.noc.sleep", "nnights", "sleepefficiency.metric", "possible_nap_edge_acc") boolean_params = c("ignorenonwear", "constrain2range", "HASPT.ignore.invalid", "relyonguider", "sleeplogidnum") character_params = c("HASPT.algo", "HASIB.algo", "Sadeh_axis", "nap_model", diff --git a/R/g.part5.analyseRest.R b/R/g.part5.analyseRest.R index 026d06bdf..8ae258b57 100644 --- a/R/g.part5.analyseRest.R +++ b/R/g.part5.analyseRest.R @@ -1,8 +1,16 @@ g.part5.analyseRest = function(sibreport = NULL, dsummary = NULL, ds_names = NULL, fi = NULL, di = NULL, - time = NULL, tz = NULL, possible_nap_dur = NULL) { + time = NULL, tz = NULL, possible_nap_dur = 0, + possible_nap_edge_acc = Inf) { # Only consider sib episodes with minimum duration - longboutsi = which((sibreport$type == "sib" & sibreport$duration >= possible_nap_dur) | + if (length(grep(pattern = "mean_acc_1min", x = colnames(sibreport))) > 0) { + sibreport$acc_edge = pmax(sibreport$mean_acc_1min_before, sibreport$mean_acc_1min_after) + } else { + sibreport$acc_edge = 0 + } + longboutsi = which((sibreport$type == "sib" & + sibreport$duration >= possible_nap_dur[1] & + sibreport$acc_edge <= possible_nap_edge_acc) | (sibreport$type != "sib" & sibreport$duration >= 1)) # for qc purposes: dsummary[di,fi] = length(longboutsi) diff --git a/R/load_params.R b/R/load_params.R index 5dfa437c1..c653d918f 100644 --- a/R/load_params.R +++ b/R/load_params.R @@ -26,7 +26,8 @@ load_params = function(group = c("sleep", "metrics", "rawdata", sleeplogsep = NULL, sleepwindowType = "SPT", possible_nap_window = c(9, 18), possible_nap_dur = c(15, 240), - nap_model = c(), sleepefficiency.metric = 1) + nap_model = c(), sleepefficiency.metric = 1, + possible_nap_edge_acc = Inf) } if ("metrics" %in% group) { params_metrics = list(do.anglex = FALSE, do.angley = FALSE, do.anglez = TRUE, diff --git a/man/GGIR.Rd b/man/GGIR.Rd index 86698da3b..9c7107c0d 100755 --- a/man/GGIR.Rd +++ b/man/GGIR.Rd @@ -1104,7 +1104,10 @@ GGIR(mode = 1:5, wake times). If 2, sleep efficiency is calculated as detected nocturnal sleep time divided by detected duration in sleep period time plus sleep latency (where sleep latency refers to the difference between log-derived sleep onset and - accelerometer-detected sleep onset). + accelerometer-detected sleep onset).} + + \item{possible_nap_edge_acc}{ + Minimum acceleration before or after the SIB for the nap to be considered. } } } @@ -1293,7 +1296,9 @@ GGIR(mode = 1:5, In \link{g.part4}: To indicate whether to generate report for the sustained inactivity bouts (SIB). If set to TRUE and when an advanced sleep diary is available in part 4 then part 5 will use this to generate summary statistics - on the overlap between self-reported nonwear and napping with SIB} + on the overlap between self-reported nonwear and napping with SIB. Here, + SIB cna be filter based on argument possible_nap_edge_acc and the first value + of possible_nap_dur} \item{do.visual}{ Boolean (default = TRUE). diff --git a/man/g.part5.analyseRest.Rd b/man/g.part5.analyseRest.Rd index e6c879e77..d4809f421 100644 --- a/man/g.part5.analyseRest.Rd +++ b/man/g.part5.analyseRest.Rd @@ -10,7 +10,9 @@ \usage{ g.part5.analyseRest(sibreport = NULL, dsummary = NULL, ds_names = NULL, fi = NULL, - di = NULL, time = NULL, tz = NULL) + di = NULL, time = NULL, tz = NULL, + possible_nap_dur = 0, + possible_nap_edge_acc = Inf) } \arguments{ \item{sibreport}{ @@ -35,6 +37,13 @@ \item{tz}{ Timezone database name } + \item{possible_nap_dur}{ + Minimum sib duration to be considered. For self-reported naps/nonwear + there is a minimum duration of 1 epoch. + } + \item{possible_nap_edge_acc}{ + Maximum acceleration before or after the SIB for it to be considered. + } } \value{ List with updated objects dsummary, ds_names, fi, and di diff --git a/man/g.part5_analyseSegment.Rd b/man/g.part5_analyseSegment.Rd index c73d3d055..6a24a92af 100644 --- a/man/g.part5_analyseSegment.Rd +++ b/man/g.part5_analyseSegment.Rd @@ -19,8 +19,7 @@ add_one_day_to_next_date, lightpeak_available, tail_expansion_log, - foldernamei, sibreport = NULL, - possible_nap_dur = NULL) + foldernamei, sibreport = NULL) } \arguments{ \item{indexlog}{ @@ -87,10 +86,6 @@ \item{sibreport}{ Sibreport object as passed on from \link{g.part5} } - \item{possible_nap_dur}{ - Minimum sib duration to be considered. For self-reported naps/nonwear - there is a minimum duration of 1 epoch. - } } \value{ List with objects: indexlog, timeList, and the matrix with the prelimenary results diff --git a/tests/testthat/test_load_check_params.R b/tests/testthat/test_load_check_params.R index 4ee135ad1..4f5bffac4 100644 --- a/tests/testthat/test_load_check_params.R +++ b/tests/testthat/test_load_check_params.R @@ -10,7 +10,7 @@ test_that("load_params can load parameters", { # Test length of objects expect_equal(length(params), 8) - expect_equal(length(params$params_sleep), 21) + expect_equal(length(params$params_sleep), 22) expect_equal(length(params$params_metrics), 41) expect_equal(length(params$params_rawdata), 37) expect_equal(length(params$params_247), 20) From ebbe8fa4538b65763d3eafe3e0e40fc664bd617e Mon Sep 17 00:00:00 2001 From: Vincent van Hees Date: Tue, 12 Sep 2023 17:36:18 +0200 Subject: [PATCH 15/18] Update NEWS.md --- NEWS.md | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/NEWS.md b/NEWS.md index 222ac498f..2fab960ce 100644 --- a/NEWS.md +++ b/NEWS.md @@ -1,7 +1,7 @@ # CHANGES IN GGIR VERSION 2.10-3 -- Part 5: Performance an analysis of the overlap between self-reported naps and nonwear -with sustained inactivity bouts. When Fixes #687 +- Part 5: Now able to assess overlap between self-reported naps and nonwear +with sustained inactivity bouts in order to facilitate research on nap detection. Fixes #687 - Visualreport: Bug fixed extracting the numeric value of the days to exclude in g.plot5 #879 From cf90d601493ace0bacb9c67056e50dc87d758b97 Mon Sep 17 00:00:00 2001 From: Vincent van Hees Date: Tue, 12 Sep 2023 17:58:00 +0200 Subject: [PATCH 16/18] fix typo in documentation --- vignettes/GGIR.Rmd | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/vignettes/GGIR.Rmd b/vignettes/GGIR.Rmd index a9394d6d7..c6712fc60 100644 --- a/vignettes/GGIR.Rmd +++ b/vignettes/GGIR.Rmd @@ -982,8 +982,8 @@ vigorous activity: | noverl_X | Only created if `do.sibreport = TRUE`. Number of overlapping bouts in a day of X where X can be sib_srnap, sib_srnonw, srnap_sib, or srnonw_sib | | frag_mean_dur_X_day | Only created if `do.sibreport = TRUE`. Mean duration of X per day, where X can be sib, srnap or srnonw | | dur_day_X_min | Only created if `do.sibreport = TRUE`. Total duration in minutes of X per day, where X can be sib, srnap or srnonw | -| mdur_X_overl_Y | Only created if `do.sibreport = TRUE`. Mean duration of the overlap between X and Y, which are combinations of be sib, srnap or srnonw | -| tdur_X_overl_Y | Only created if `do.sibreport = TRUE`. Total duration in minutes of the overlap between X and Y, which are combinations of be sib, srnap or srnonw | +| mdur_X_overl_Y | Only created if `do.sibreport = TRUE`. Mean duration of the overlap between X and Y, which are combinations of sib, srnap or srnonw | +| tdur_X_overl_Y | Only created if `do.sibreport = TRUE`. Total duration in minutes of the overlap between X and Y, which are combinations of sib, srnap or srnonw | | perc_X_overl_Y | Only created if `do.sibreport = TRUE`. Percentage of overlap between X and Y, which are combinations of sib, srnap or srnonw | From c80331edc8b0c5dbc462a7d308afc80a99dab9b0 Mon Sep 17 00:00:00 2001 From: Jairo H Migueles Date: Wed, 13 Sep 2023 09:57:16 +0200 Subject: [PATCH 17/18] transform time in g.part5.analyseRest if needed --- R/g.part5.analyseRest.R | 5 +++++ 1 file changed, 5 insertions(+) diff --git a/R/g.part5.analyseRest.R b/R/g.part5.analyseRest.R index 8ae258b57..ed57cedf6 100644 --- a/R/g.part5.analyseRest.R +++ b/R/g.part5.analyseRest.R @@ -2,6 +2,11 @@ g.part5.analyseRest = function(sibreport = NULL, dsummary = NULL, ds_names = NULL, fi = NULL, di = NULL, time = NULL, tz = NULL, possible_nap_dur = 0, possible_nap_edge_acc = Inf) { + # transform time to POSIX + if (is.ISO8601(time[1])) { + time = iso8601chartime2POSIX(time, tz = tz) + } + # Only consider sib episodes with minimum duration if (length(grep(pattern = "mean_acc_1min", x = colnames(sibreport))) > 0) { sibreport$acc_edge = pmax(sibreport$mean_acc_1min_before, sibreport$mean_acc_1min_after) From ee01f8037822d9283c3782f730937bb758ecc5ab Mon Sep 17 00:00:00 2001 From: Jairo H Migueles Date: Wed, 13 Sep 2023 10:14:04 +0200 Subject: [PATCH 18/18] fix issue introduced in previous commit, make sure time is character when testing if it is ISO8601 --- R/g.part5.analyseRest.R | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/R/g.part5.analyseRest.R b/R/g.part5.analyseRest.R index ed57cedf6..fefedf6c8 100644 --- a/R/g.part5.analyseRest.R +++ b/R/g.part5.analyseRest.R @@ -3,7 +3,7 @@ g.part5.analyseRest = function(sibreport = NULL, dsummary = NULL, time = NULL, tz = NULL, possible_nap_dur = 0, possible_nap_edge_acc = Inf) { # transform time to POSIX - if (is.ISO8601(time[1])) { + if (is.ISO8601(as.character(time[1]))) { time = iso8601chartime2POSIX(time, tz = tz) }