diff --git a/.DS_Store b/.DS_Store index 7afad0f..dfabc07 100644 Binary files a/.DS_Store and b/.DS_Store differ diff --git a/.gitignore b/.gitignore index 2216890..42f07ba 100644 --- a/.gitignore +++ b/.gitignore @@ -1,55 +1,40 @@ # History files .Rhistory .Rapp.history - # Session Data files .RData .RDataTmp - # User-specific files .Ruserdata - # Example code in package build process *-Ex.R - # Output files from R CMD build /*.tar.gz - # Output files from R CMD check /*.Rcheck/ - # RStudio files .Rproj.user/ - # produced vignettes vignettes/*.html vignettes/*.pdf - # OAuth2 token, see https://github.com/hadley/httr/releases/tag/v0.3 .httr-oauth - # knitr and R markdown default cache directories *_cache/ /cache/ - # Temporary files created by R markdown *.utf8.md *.knit.md - # R Environment Variables .Renviron - # pkgdown site docs/ - # translation temp files po/*~ - # RStudio Connect folder rsconnect/ .Rproj.user docs inst/doc - # DS_Store file .DS_Store diff --git a/DESCRIPTION b/DESCRIPTION index 7ce764a..fe5fac0 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -1,6 +1,6 @@ Package: LightLogR Title: Work With Data from Wearable Light Loggers and Optical Radiation Dosimeters -Version: 0.3.2 +Version: 0.3.3 Authors@R: c( person("Johannes", "Zauner", email = "johannes.zauner@tum.de", role = c("aut", "cre"), diff --git a/NAMESPACE b/NAMESPACE index 4e9fd6c..66a7158 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -7,6 +7,7 @@ export(Brown_rec) export(Datetime_breaks) export(Datetime_limits) export(aggregate_Datetime) +export(barroso_lighting_metrics) export(bright_dark_period) export(centroidLE) export(count_difftime) @@ -18,10 +19,12 @@ export(dominant_epoch) export(dst_change_handler) export(dst_change_summary) export(duration_above_threshold) +export(exponential_moving_average) export(filter_Date) export(filter_Datetime) export(filter_Datetime_multiple) export(filter_Time) +export(frequency_crossing_threshold) export(gap_finder) export(gap_handler) export(gapless_Datetimes) @@ -43,9 +46,12 @@ export(nvRC_circadianDisturbance) export(nvRC_relativeAmplitudeError) export(nvRD) export(nvRD_cumulative_response) +export(period_above_threshold) +export(pulses_above_threshold) export(sc2interval) export(sleep_int2Brown) export(symlog_trans) +export(threshold_for_duration) export(timing_above_threshold) importFrom(magrittr,"%>%") importFrom(rlang,":=") diff --git a/NEWS.md b/NEWS.md index d9e30e3..570282b 100644 --- a/NEWS.md +++ b/NEWS.md @@ -1,3 +1,7 @@ +# LightLogR 0.3.3 + +* New and updated metric functions. LightLogR now contains 16 metric families with 60 sub-metrics. + # LightLogR 0.3.2 * added import functions for `nanoLambda`and `LightWatcher` devices diff --git a/R/helper.R b/R/helper.R index 185c165..349cd21 100755 --- a/R/helper.R +++ b/R/helper.R @@ -16,9 +16,9 @@ test.Time.regex <- function(input) { if(!test){ stop(paste("input:", - rlang::as_string( - input.defused), - "needs to be in the format 'hh:mm:ss'")) + rlang::as_string( + input.defused), + "needs to be in the format 'hh:mm:ss'")) } } } @@ -45,9 +45,9 @@ compare_difftime <- function(dataset1, dataset2, Datetime.colname = Datetime, n #do a full join with every column but Quantile group_variables <- setdiff(names(Quant2), "Quant") dplyr::full_join(Quant1, Quant2, by = group_variables) %>% - dplyr::mutate( - comparison = Quant.x <= Quant.y - ) + dplyr::mutate( + comparison = Quant.x <= Quant.y + ) } #calculate whether any of the comparisons in compare_difftime is FALSE @@ -55,7 +55,7 @@ compare_difftime.any <- function(...) { comparison <- compare_difftime(...) %>% dplyr::filter(comparison == FALSE) %>% dplyr::rename(Dataset.Interval = Quant.x, - Reference.Interval = Quant.y) %>% + Reference.Interval = Quant.y) %>% dplyr::select(-comparison) if(nrow(comparison) > 0) comparison else TRUE @@ -74,10 +74,10 @@ create.Reference.label <- function(dataset, dplyr::mutate(!!Reference.label.column.str := dplyr::if_else( !is.na({{ Reference.column }}), Reference.label, NA - ) + ) ) dataset - } else dataset + } else dataset } #helper to pick the colums that are used for grouping @@ -89,9 +89,9 @@ pick.grouping.columns <- function(dataset) { # Compare with threshold compare_threshold <- function(Light.vector, - threshold, - comparison = c("above", "below"), - na.replace = FALSE){ + threshold, + comparison = c("above", "below"), + na.replace = FALSE){ comparison = match.arg(comparison) @@ -119,6 +119,105 @@ compare_threshold <- function(Light.vector, return(out) } +# Find clusters of consecutive `TRUE` values in logical vector. +# Allows to concatenate periods of consecutive values that are interrupted +# by periods of `FALSE` values. +find_clusters <- function(x, + min.length, + max.interrupt = 0, + prop.interrupt = 1, + cluster_name = "cluster") { + + stopifnot( + "`x` must be logical" = is.logical(x), + "`min.length` must be larger than 0" = min.length > 0, + "`max.interrupt` must be larger than or equal to 0" = max.interrupt >= 0, + "`prop.interrupt` must be between 0 and 1" = + prop.interrupt >= 0 & prop.interrupt <= 1 + ) + + # Replace NA with FALSE + x <- tidyr::replace_na(x, FALSE) + + # Find the start and end indices of each cluster of consecutive values + start_indices <- which(x & !dplyr::lag(x, default = FALSE)) + end_indices <- which(x & !dplyr::lead(x, default = FALSE)) + ranges <- as.numeric(matrix(rbind(start_indices, end_indices), nrow = 2)) + + # Remove ranges < min.length + intra_diff <- diff(ranges)[1:(length(ranges) - 1) %% 2 != 0] + 1 + exclude_ranges <- c( + which(intra_diff < min.length) * 2, + which(intra_diff < min.length) * 2 - 1 + ) + if (length(exclude_ranges) > 0) { + ranges <- ranges[-exclude_ranges] + } + + # Calculate cumulative ranges + intra_diff <- diff(ranges)[1:(length(ranges) - 1) %% 2 != 0] + 1 + intra_cumsum <- intra_diff[1:length(intra_diff)-1] + intra_diff[2:length(intra_diff)] + + # Inter-range differences + inter_diff <- diff(ranges)[1:(length(ranges) - 1) %% 2 == 0] - 1 + + # Proportion inter-range difference and cumulative range sums + interrupt_ratio <- inter_diff / intra_cumsum + + # Combine ranges with inter-range difference <= max.interrupt & + # interrupt ratio <= prop.interrupt + exclude_ranges <- c( + which(inter_diff <= max.interrupt & interrupt_ratio <= prop.interrupt) * 2, + which(inter_diff <= max.interrupt & interrupt_ratio <= prop.interrupt) * 2 + 1 + ) + if (length(exclude_ranges) > 0) { + ranges <- ranges[-exclude_ranges] + } + + # Make data frame with intervals + if (length(ranges) > 0) { + intervals <- + matrix(ranges, ncol = 2, byrow = TRUE) %>% + as.data.frame() %>% + magrittr::set_names(c("cluster_start", "cluster_end")) %>% + dplyr::mutate(cluster_idx = 1:length(cluster_start)) %>% + dplyr::rowwise() %>% + dplyr::mutate(row_idx = list(seq(cluster_start, cluster_end))) %>% + tidyr::unnest(cols = c(row_idx)) %>% + dplyr::ungroup() %>% + dplyr::mutate(is_cluster = TRUE) %>% + dplyr::relocate(row_idx, is_cluster, cluster_idx, cluster_start, cluster_end) + } else { + intervals <- + tibble::tibble( + row_idx = 1, is_cluster = NA, cluster_idx = NA, + cluster_start = NA, cluster_end = NA + ) + } + # Replace "cluster" with custom name + intervals <- intervals %>% + dplyr::rename_with(~gsub("cluster", cluster_name, .x)) + + return(intervals) +} + +# Convert `x` to time scale of `t` +convert_to_timescale <- function(x, t){ + if(lubridate::is.POSIXct(t)){ + x <- lubridate::as_datetime(x, tz = lubridate::tz(t)) + } + if(hms::is_hms(t)){ + x <- hms::as_hms(x) + } + if(lubridate::is.duration(t)){ + x <- lubridate::as.duration(x) + } + if(lubridate::is.difftime(t) & !hms::is_hms(t)){ + x <- lubridate::as.difftime(x, unit = units(t)) + } + return(x) +} + # Create the list of epochs, which are either dominant or not epoch_list <- function(dataset = dataset, Datetime.colname = Datetime, @@ -135,4 +234,3 @@ epoch_list <- function(dataset = dataset, epochs } - diff --git a/R/metric_EMA.R b/R/metric_EMA.R new file mode 100644 index 0000000..5889c47 --- /dev/null +++ b/R/metric_EMA.R @@ -0,0 +1,101 @@ +#' Exponential moving average filter (EMA) +#' +#' This function smoothes the data using an exponential moving average filter +#' with a specified decay half-life. +#' +#' @param Light.vector Numeric vector containing the light data. Missing values are +#' replaced by 0. +#' @param Time.vector Vector containing the time data. Can be \link[base]{POSIXct}, \link[hms]{hms}, +#' \link[lubridate]{duration}, or \link[base]{difftime}. +#' @param decay The decay half-life controlling the exponential smoothing. +#' Can be either a \link[lubridate]{duration} or a string. If it is a string, it +#' needs to be a valid \link[lubridate]{duration} string, e.g., `"1 day"` or `"10 sec"`. +#' The default is set to `"90 mins"` for a biologically relevant estimate (see +#' the reference paper). +#' @param epoch The epoch at which the data was sampled. Can be either a +#' \link[lubridate]{duration} or a string. If it is a string, it needs to be +#' either `"dominant.epoch"` (the default) for a guess based on the data, or a valid +#' \link[lubridate]{duration} string, e.g., `"1 day"` or `"10 sec"`. +#' +#' @return A numeric vector containing the smoothed light data. The output has the same +#' length as `Light.vector`. +#' +#' @export +#' +#' @family metrics +#' +#' @details The timeseries is assumed to be regular. Missing values in the +#' light data will be replaced by 0. +#' +#' @references +#' Price, L. L. A. (2014). On the Role of Exponential Smoothing in Circadian +#' Dosimetry. \emph{Photochemistry and Photobiology}, 90(5), 1184-1192. +#' \url{https://doi.org/10.1111/php.12282} +#' +#' Hartmeyer, S.L., Andersen, M. (2023). Towards a framework for light-dosimetry studies: +#' Quantification metrics. \emph{Lighting Research & Technology}. +#' \url{https://doi.org/10.1177/14771535231170500} +#' +#' @examples +#' sample.data.environment.EMA = sample.data.environment %>% +#' dplyr::filter(Id == "Participant") %>% +#' filter_Datetime(length = lubridate::days(2)) %>% +#' dplyr::mutate(MEDI.EMA = exponential_moving_average(MEDI, Datetime)) +#' +#' # Plot to compare results +#' sample.data.environment.EMA %>% +#' ggplot2::ggplot(ggplot2::aes(x = Datetime)) + +#' ggplot2::geom_line(ggplot2::aes(y = MEDI), colour = "black") + +#' ggplot2::geom_line(ggplot2::aes(y = MEDI.EMA), colour = "red") +#' +exponential_moving_average <- function(Light.vector, + Time.vector, + decay = "90 min", + epoch = "dominant.epoch") { + + # Perform argument checks + stopifnot( + "`Light.vector` must be numeric!" = is.numeric(Light.vector), + "`Time.vector` must be POSIXct, hms, duration, or difftime!" = + lubridate::is.POSIXct(Time.vector) | hms::is_hms(Time.vector) | + lubridate::is.duration(Time.vector) | lubridate::is.difftime(Time.vector), + "`Light.vector` and `Time.vector` must be same length!" = + length(Light.vector) == length(Time.vector), + "`decay` must either be a duration or a string" = + lubridate::is.duration(decay) | is.character(decay), + "`epoch` must either be a duration or a string" = + lubridate::is.duration(epoch) | is.character(epoch) + ) + + # Get the epochs based on the data + if (is.character(epoch) && epoch == "dominant.epoch") { + epoch <- count_difftime(tibble::tibble(Datetime = Time.vector))$difftime[1] + } + # If the user specified an epoch, use that instead + else { + epoch <- lubridate::as.duration(epoch) + } + + # Replace missing values + if (any(is.na(Light.vector))) { + warning("Light data contains missing values! They are replaced by 0.") + Light.vector[is.na(Light.vector)] <- 0 + } + + # Calculate smoothing factor beta + decay <- lubridate::as.duration(decay) + beta <- log(2) / (as.numeric(decay) / as.numeric(epoch)) + + # EMA filter + D <- double(length(Light.vector)) + for (idx in 1:length(Light.vector)) { + if (idx == 1) { + D[idx] <- beta * (Light.vector[idx]) + } else { + D[idx] <- D[idx - 1] + beta * (Light.vector[idx] - D[idx - 1]) + } + } + + # Return numeric vector of EMA light values + return(D) +} diff --git a/R/metric_barroso.R b/R/metric_barroso.R new file mode 100644 index 0000000..12e7743 --- /dev/null +++ b/R/metric_barroso.R @@ -0,0 +1,136 @@ +#' Circadian lighting metrics from Barroso et al. (2014) +#' +#' This function calculates the metrics proposed by Barroso et al. (2014) +#' for light-dosimetry in the context of research on the non-visual effects of light. +#' The following metrics are calculated: +#' +#' \describe{ +#' \item{`bright_threshold`}{The maximum light intensity for which at least six +#' hours of measurements are at the same or higher level.} +#' \item{`dark_threshold`}{The minimum light intensity for which at least eight +#' hours of measurements are at the same or lower level.} +#' \item{`bright_mean_level`}{The 20% trimmed mean of all light intensity measurements +#' equal or above the `bright_threshold`.} +#' \item{`dark_mean_level`}{The 20% trimmed mean of all light intensity measurements +#' equal or below the `dark_threshold`.} +#' \item{`bright_cluster`}{The longest continuous time interval above the `bright_threshold`.} +#' \item{`dark_cluster`}{The longest continuous time interval below the `dark_threshold`.} +#' \item{`circadian_variation`}{A measure of periodicity of the daily lighting +#' schedule over a given set of days. Calculated as the coefficient of variation +#' of input light data. +#' } +#' } +#' +#' @param Light.vector Numeric vector containing the light data. +#' @param Time.vector Vector containing the time data. Can be \link[base]{POSIXct}, \link[hms]{hms}, +#' \link[lubridate]{duration}, or \link[base]{difftime}. +#' @param epoch The epoch at which the data was sampled. Can be either a +#' \link[lubridate]{duration} or a string. If it is a string, it needs to be +#' either `"dominant.epoch"` (the default) for a guess based on the data, or a valid +#' \link[lubridate]{duration} string, e.g., `"1 day"` or `"10 sec"`. +#' @param loop Logical. Should the data be looped? Defaults to `FALSE`. +#' @param na.rm Logical. Should missing values (NA) be removed for the calculation? +#' Defaults to `FALSE`. If `TRUE`, for the calculation of `bright_cluster` and +#' `dark_cluster`, missing values will be replaced by 0 +#' (see \code{\link{period_above_threshold}}). +#' @param as.df Logical. Should a data frame be returned? If `TRUE`, a data +#' frame with seven columns will be returned. Defaults to `FALSE`. +#' +#' @return List or dataframe with the seven values: `bright_threshold`, `dark_threshold`, +#' `bright_mean_level`, `dark_mean_level`, `bright_cluster`, `dark_cluster`, +#' `circadian_variation`. The output type of `bright_cluster`, `dark_cluster`, +#' is a \link[lubridate]{duration} object. +#' +#' @details +#' +#' @export +#' +#' @references +#' Barroso, A., Simons, K., & Jager, P. de. (2014). Metrics of circadian +#' lighting for clinical investigations. \emph{Lighting Research & Technology}, +#' 46(6), 637–649. \url{https://doi.org/10.1177/1477153513502664} +#' +#' Hartmeyer, S.L., Andersen, M. (2023). Towards a framework for light-dosimetry studies: +#' Quantification metrics. \emph{Lighting Research & Technology}. +#' \url{https://doi.org/10.1177/14771535231170500} +#' +#' @examples +#' +#' dataset1 <- +#' tibble::tibble( +#' Id = rep("B", 60 * 24), +#' Datetime = lubridate::as_datetime(0) + lubridate::minutes(0:(60*24-1)), +#' MEDI = c(rep(sample(seq(0,1,0.1), 60*8, replace = TRUE)), +#' rep(sample(1:1000, 16, replace = TRUE), each = 60)) +#' ) +#' +#' dataset1 %>% +#' dplyr::reframe(barroso_lighting_metrics(MEDI, Datetime, as.df = TRUE)) +#' +barroso_lighting_metrics <- function(Light.vector, + Time.vector, + epoch = "dominant.epoch", + loop = FALSE, + na.rm = FALSE, + as.df = FALSE) { + + # Perform argument checks + stopifnot( + "`Light.vector` must be numeric!" = is.numeric(Light.vector), + "`Time.vector` must be POSIXct, hms, duration, or difftime!" = + lubridate::is.POSIXct(Time.vector) | hms::is_hms(Time.vector) | + lubridate::is.duration(Time.vector) | lubridate::is.difftime(Time.vector), + "`Light.vector` and `Time.vector` must be same length!" = + length(Light.vector) == length(Time.vector), + "`epoch` must either be a duration or a string" = + lubridate::is.duration(epoch) | is.character(epoch), + "`loop` must be logical!" = is.logical(loop), + "`na.rm` must be logical!" = is.logical(na.rm), + "`as.df` must be logical!" = is.logical(as.df) + ) + + # Get the epochs based on the data + if (is.character(epoch) && epoch == "dominant.epoch") { + epoch <- count_difftime(tibble::tibble(Datetime = Time.vector))$difftime[1] + } + # If the user specified an epoch, use that instead + else { + epoch <- lubridate::as.duration(epoch) + } + + # Bright/dark thresholds + bright_threshold <- + threshold_for_duration(Light.vector, Time.vector, "6 h", "above", epoch, na.rm = na.rm)[1] + dark_threshold <- + threshold_for_duration(Light.vector, Time.vector, "8 h", "below", epoch, na.rm = na.rm)[1] + + # Bright/dark mean level + bright_mean_level <- mean(Light.vector[Light.vector >= bright_threshold], trim = 0.2, na.rm = na.rm) + dark_mean_level <- mean(Light.vector[Light.vector <= dark_threshold], trim = 0.2, na.rm = na.rm) + + # Bright/dark cluster + bright_cluster <- period_above_threshold(Light.vector, Time.vector, "above", bright_threshold, + epoch = epoch, loop = loop, na.replace = na.rm)[1] + dark_cluster <- period_above_threshold(Light.vector, Time.vector, "below", dark_threshold, + epoch = epoch, loop = loop, na.replace = na.rm)[1] + + # Circadian variation + circadian_variation <- (stats::sd(Light.vector, na.rm = na.rm) / mean(Light.vector, na.rm=na.rm)) %>% + round(2) + + # Prepare output + out <- list( + "bright_threshold" = bright_threshold, + "dark_threshold" = dark_threshold, + "bright_mean_level" = bright_mean_level, + "dark_mean_level" = dark_mean_level, + "bright_cluster" = bright_cluster, + "dark_cluster" = dark_cluster, + "circadian_variation" = circadian_variation + ) + # Return data frame or list + if (as.df) { + out <- tibble::as_tibble(out) + } + return(out) +} diff --git a/R/metric_bright_dark_period.R b/R/metric_bright_dark_period.R index 5b71821..63d49f8 100755 --- a/R/metric_bright_dark_period.R +++ b/R/metric_bright_dark_period.R @@ -7,16 +7,17 @@ #' (i.e., no gaps) for correct results. #' #' @param Light.vector Numeric vector containing the light data. -#' @param Time.vector Vector containing the time data. Can be HMS or POSIXct. +#' @param Time.vector Vector containing the time data. Can be \link[base]{POSIXct}, +#' \link[hms]{hms}, \link[lubridate]{duration}, or \link[base]{difftime}. #' @param period String indicating the type of period to look for. Can be either #' `"brightest"`(the default) or `"darkest"`. #' @param timespan The timespan across which to calculate. Can be either a -#' `lubridate::duration()` or a `lubridate::duration()` string, e.g., +#' \link[lubridate]{duration} or a \link[lubridate]{duration} string, e.g., #' `"1 day"` or `"10 sec"`. #' @param epoch The epoch at which the data was sampled. Can be either a -#' `lubridate::duration()` or a string. If it is a string, it needs to be -#' either `"dominant.epoch"` (the default) for a guess based on the data, or a valid -#' `lubridate::duration()` string, e.g., `"1 day"` or `"10 sec"`. +#' \link[lubridate]{duration} or a string. If it is a string, it needs to be +#' either `"dominant.epoch"` (the default) for a guess based on the data, or a valid +#' \link[lubridate]{duration} string, e.g., `"1 day"` or `"10 sec"`. #' @param loop Logical. Should the data be looped? If `TRUE`, a full copy of the data #' will be concatenated at the end of the data. Makes only sense for 24 h data. #' Defaults to `FALSE`. @@ -29,7 +30,7 @@ #' calculated brightest or darkest period, or if `as.df == TRUE` a data frame #' with columns named `{period}_{timespan}_{metric}`. The output type corresponds #' to the type of `Time.vector`, e.g., if `Time.vector` is HMS, the timing metrics -#' will be also HMS, and vice versa for POSIXct and numeric. +#' will be also HMS, and vice versa for POSIXct. #' #' @details Assumes regular 24h light data. Otherwise, results may not be #' meaningful. Looping the data is recommended for finding the darkest period. @@ -44,7 +45,7 @@ #' @family metrics #' #' @examples -# Dataset with light > 250lx between 06:00 and 18:00 +#' # Dataset with light > 250lx between 06:00 and 18:00 #' dataset1 <- #' tibble::tibble( #' Id = rep("A", 24), @@ -58,7 +59,22 @@ #' dataset1 %>% #' dplyr::reframe(bright_dark_period(MEDI, Datetime, "darkest", "5 hours", #' loop = TRUE, as.df = TRUE)) - +#' +#' # Dataset with duration as Time.vector +#' dataset2 <- +#' tibble::tibble( +#' Id = rep("A", 24), +#' Datetime = lubridate::dhours(0:23), +#' MEDI = c(rep(1, 6), rep(250, 13), rep(1, 5)) +#' ) +#' +#' dataset2 %>% +#' dplyr::reframe(bright_dark_period(MEDI, Datetime, "brightest", "10 hours", +#' as.df = TRUE)) +#' dataset2 %>% +#' dplyr::reframe(bright_dark_period(MEDI, Datetime, "darkest", "5 hours", +#' loop = TRUE, as.df = TRUE)) +#' bright_dark_period <- function(Light.vector, Time.vector, period = c("brightest", "darkest"), @@ -73,8 +89,11 @@ bright_dark_period <- function(Light.vector, # Perform argument checks stopifnot( "`Light.vector` must be numeric!" = is.numeric(Light.vector), - "`Time.vector` must be POSIXct or HMS" = lubridate::is.POSIXct(Time.vector) | - hms::is_hms(Time.vector), + "`Time.vector` must be POSIXct, hms, duration, or difftime!" = + lubridate::is.POSIXct(Time.vector) | hms::is_hms(Time.vector) | + lubridate::is.duration(Time.vector) | lubridate::is.difftime(Time.vector), + "`Light.vector` and `Time.vector` must be same length!" = + length(Light.vector) == length(Time.vector), "`epoch` must either be a duration or a string" = lubridate::is.duration(epoch) | is.character(epoch), "`timespan` must either be a duration or a string" = @@ -89,7 +108,7 @@ bright_dark_period <- function(Light.vector, } # Get the epochs based on the data - if (epoch == "dominant.epoch") { + if (is.character(epoch) && epoch == "dominant.epoch") { epoch <- count_difftime(tibble::tibble(Datetime = Time.vector))$difftime[1] } # If the user specified an epoch, use that instead diff --git a/R/metric_centroidLE.R b/R/metric_centroidLE.R index 0c8330e..3212e39 100755 --- a/R/metric_centroidLE.R +++ b/R/metric_centroidLE.R @@ -4,16 +4,16 @@ #' time vector weighted in proportion to the corresponding binned light intensity. #' #' @param Light.vector Numeric vector containing the light data. -#' @param Time.vector Vector containing the time data. Can be numeric, HMS or POSIXct. +#' @param Time.vector Vector containing the time data. Can be \link[base]{POSIXct}, \link[hms]{hms}, +#' \link[lubridate]{duration}, or \link[base]{difftime}. #' @param bin.size Value specifying size of bins to average the light data over. -#' If `Time.vector` is of type POSIXct or HMS, `bin.size` must be either a -#' `lubridate::duration()` or a `lubridate::duration()` string, e.g., -#' `"1 day"` or `"10 sec"`. Otherwise, if `Time.vector` is numeric, `bin.size` -#' must be also numeric. If nothing is provided, no binning will be performed. +#' Must be either a \link[lubridate]{duration} or a \link[lubridate]{duration} string, e.g., +#' `"1 day"` or `"10 sec"`. If nothing is provided, no binning will be performed. #' @param na.rm Logical. Should missing values be removed for the calculation? #' Defaults to `FALSE`. -#' @param as.df Logical. Should the output be returned as a data frame? Defaults -#' to `FALSE`. +#' @param as.df Logical. Should the output be returned as a data frame? If `TRUE`, a data +#' frame with a single column named `centroidLE` will be returned. +#' Defaults to `FALSE`. #' #' @return Single column data frame or vector. #' @@ -44,7 +44,7 @@ #' "Centroid of light exposure" = centroidLE(MEDI, Datetime, "2 hours") #' ) #' -#' # Dataset with HMS time vector +#' # Dataset with hms time vector #' dataset2 <- #' tibble::tibble( #' Id = rep("A", 24), @@ -56,16 +56,16 @@ #' "Centroid of light exposure" = centroidLE(MEDI, Time, "2 hours") #' ) #' -#' # Dataset with numeric time vector +#' # Dataset with duration time vector #' dataset3 <- #' tibble::tibble( #' Id = rep("A", 24), -#' Hour = 0:23, +#' Hour = lubridate::duration(0:23, "hours"), #' MEDI = c(rep(1, 6), rep(250, 13), rep(1, 5)) #' ) #' dataset3 %>% #' dplyr::reframe( -#' "Centroid of light exposure" = centroidLE(MEDI, Hour, 2) +#' "Centroid of light exposure" = centroidLE(MEDI, Hour, "2 hours") #' ) #' centroidLE <- function(Light.vector, @@ -77,29 +77,27 @@ centroidLE <- function(Light.vector, # Perform argument checks stopifnot( "`Light.vector` must be numeric!" = is.numeric(Light.vector), - "`Time.vector` must be numeric, HMS, or POSIXct" = - is.numeric(Time.vector) | hms::is_hms(Time.vector) | lubridate::is.POSIXct(Time.vector), + "`Time.vector` must be POSIXct, hms, duration, or difftime!" = + lubridate::is.POSIXct(Time.vector) | hms::is_hms(Time.vector) | + lubridate::is.duration(Time.vector) | lubridate::is.difftime(Time.vector), + "`Light.vector` and `Time.vector` must be same length!" = + length(Light.vector) == length(Time.vector), "`na.rm` must be logical!" = is.logical(na.rm), "`as.df` must be logical!" = is.logical(as.df) ) if (!is.null(bin.size)) { - if (lubridate::is.POSIXct(Time.vector) | hms::is_hms(Time.vector)) { - stopifnot("`bin.size` must be a either a `lubridate::duration` object or `lubridate::duration` string, because `Time.vector` is HMS or POSIXct" = - !is.numeric(bin.size)) - bin.size <- lubridate::duration(bin.size) - stopifnot("`bin.size` must be a either a `lubridate::duration` object or `lubridate::duration` string, because `Time.vector` is HMS or POSIXct" = - !is.na(bin.size) & lubridate::is.duration(bin.size)) - bin.size <- lubridate::as.period(bin.size) - } - else { - stopifnot("`bin.size` must be numeric because `Time.vector` is numeric" = - is.numeric(bin.size)) - } + stopifnot("`bin.size` must be a either a duration or a string" = + lubridate::is.duration(bin.size) | is.character(bin.size)) + bin.size <- lubridate::as.period(bin.size) } # Make tibble df <- tibble::tibble(Light = Light.vector, Time = Time.vector) + if(na.rm){ + df <- df %>% dplyr::filter(!is.na(Light)) + } + # Average into bins if(!is.null(bin.size)) { if (lubridate::is.POSIXct(Time.vector)){ @@ -115,25 +113,19 @@ centroidLE <- function(Light.vector, ) %>% dplyr::summarise(Light = mean(Light, na.rm = na.rm)) } - if (is.numeric(Time.vector)) { + if (lubridate::is.duration(Time.vector) | lubridate::is.difftime(Time.vector)) { df <- df %>% - dplyr::group_by(Time = (Time - Time %% bin.size)) %>% + dplyr::group_by(Time = (as.numeric(Time) - as.numeric(Time) %% as.numeric(bin.size))) %>% dplyr::summarise(Light = mean(Light, na.rm = na.rm)) } } # Calculate weighted mean weights <- (df$Light / sum(df$Light, na.rm = na.rm)) - centroidLE <- sum(as.numeric(df$Time) * weights, na.rm = na.rm) + centroidLE <- sum(as.numeric(df$Time) * weights, na.rm = na.rm) %>% round() - # Convert to right time class - if(hms::is_hms(Time.vector)) { - centroidLE <- centroidLE %>% round() %>% hms::as_hms() - } - if(lubridate::is.POSIXct(Time.vector)){ - centroidLE <- centroidLE %>% round() %>% - lubridate::as_datetime(tz = lubridate::tz(Time.vector)) - } + # Convert to corresponding time scale + centroidLE <- centroidLE %>% convert_to_timescale(Time.vector) # Return data frame or numeric vector if (as.df) { diff --git a/R/metric_disparity_index.R b/R/metric_disparity_index.R index 586689d..79d06f3 100644 --- a/R/metric_disparity_index.R +++ b/R/metric_disparity_index.R @@ -6,8 +6,9 @@ #' #' @param Light.vector Numeric vector containing the light data. #' @param na.rm Logical. Should missing values be removed? Defaults to FALSE -#' @param as.df Logical. Should the output be returned as a data frame? Defaults -#' to FALSE +#' @param as.df Logical. Should the output be returned as a data frame? If `TRUE`, a data +#' frame with a single column named `disparity_index` will be returned. +#' Defaults to `FALSE`. #' #' @return Single column data frame or vector. #' @@ -41,20 +42,31 @@ disparity_index <- function(Light.vector, na.rm = FALSE, as.df = FALSE) { + # Perform argument checks + stopifnot( + "`Light.vector` must be numeric!" = is.numeric(Light.vector), + "`na.rm` must be logical!" = is.logical(na.rm), + "`as.df` must be logical!" = is.logical(as.df) + ) + # Remove NAs if (na.rm) { Light.vector <- Light.vector[!is.na(Light.vector)] } - if (length(Light.vector) == 1) { - di <- 0 - } else { - # Calculate disparity index - fractions <- (Light.vector[2:length(Light.vector)] + 1) / - (Light.vector[1:length(Light.vector) - 1] + 1) - di <- 1 / (length(Light.vector) - 1) * sum(abs(log(fractions))) + if (any(is.na(Light.vector))){ + di <- NA + } + else{ + if (length(Light.vector) == 1) { + di <- 0 + } else { + # Calculate disparity index + fractions <- (Light.vector[2:length(Light.vector)] + 1) / + (Light.vector[1:length(Light.vector) - 1] + 1) + di <- 1 / (length(Light.vector) - 1) * sum(abs(log(fractions))) + } } - # Return as data frame or numeric vector if (as.df) { return(tibble::tibble("disparity_index" = di)) diff --git a/R/metric_duration_above_threshold.R b/R/metric_duration_above_threshold.R index 0f2998e..63fba2f 100755 --- a/R/metric_duration_above_threshold.R +++ b/R/metric_duration_above_threshold.R @@ -1,10 +1,11 @@ -#' Time above/below threshold or within threshold range +#' Duration above/below threshold or within threshold range #' -#' This function calculates the time spent above/below a specified threshold +#' This function calculates the duration spent above/below a specified threshold #' light level or within a specified range of light levels. #' #' @param Light.vector Numeric vector containing the light data. -#' @param Time.vector Vector containing the time data. Can be numeric, HMS or POSIXct. +#' @param Time.vector Vector containing the time data. Can be \link[base]{POSIXct}, \link[hms]{hms}, +#' \link[lubridate]{duration}, or \link[base]{difftime}. #' @param comparison String specifying whether the time above or below threshold #' should be calculated. Can be either `"above"` (the default) or `"below"`. If #' two values are provided for `threshold`, this argument will be ignored. @@ -12,17 +13,16 @@ #' threshold light level(s) to compare with. If a vector with two values is provided, #' the time within the two thresholds will be calculated. #' @param epoch The epoch at which the data was sampled. Can be either a -#' `lubridate::duration()` or a string. If it is a string, it needs to be +#' \link[lubridate]{duration} or a string. If it is a string, it needs to be #' either `"dominant.epoch"` (the default) for a guess based on the data, or a valid -#' `lubridate::duration()` string, e.g., `"1 day"` or `"10 sec"`. +#' \link[lubridate]{duration} string, e.g., `"1 day"` or `"10 sec"`. #' @param na.rm Logical. Should missing values (NA) be removed for the calculation? #' Defaults to `FALSE`. #' @param as.df Logical. Should a data frame with be returned? If `TRUE`, a data -#' frame with a single column named `TAT_{threshold}` will be returned. +#' frame with a single column named `duration_{comparison}_{threshold}` will be returned. #' Defaults to `FALSE`. #' -#' @return A duration object (see \code{\link[lubridate]{duration}}) as single value, -#' or single column data frame. +#' @return A \link[lubridate]{duration} object as single value, or single column data frame. #' #' @references #' Hartmeyer, S.L., Andersen, M. (2023). Towards a framework for light-dosimetry studies: @@ -34,7 +34,7 @@ #' @family metrics #' #' @examples -#' N <- 50 +#' N <- 60 #' # Dataset with epoch = 1min #' dataset1 <- #' tibble::tibble( @@ -75,8 +75,11 @@ duration_above_threshold <- function(Light.vector, # Perform argument checks stopifnot( "`Light.vector` must be numeric!" = is.numeric(Light.vector), - "`Time.vector` must be numeric, HMS, or POSIXct" = - is.numeric(Time.vector) | hms::is_hms(Time.vector) | lubridate::is.POSIXct(Time.vector), + "`Time.vector` must be POSIXct, hms, duration, or difftime!" = + lubridate::is.POSIXct(Time.vector) | hms::is_hms(Time.vector) | + lubridate::is.duration(Time.vector) | lubridate::is.difftime(Time.vector), + "`Light.vector` and `Time.vector` must be same length!" = + length(Light.vector) == length(Time.vector), "`threshold` must be numeric!" = is.numeric(threshold), "`threshold` must be either one or two values!" = length(threshold) %in% c(1, 2), "`epoch` must either be a duration or a string" = @@ -86,7 +89,7 @@ duration_above_threshold <- function(Light.vector, ) # Get the epochs based on the data - if (epoch == "dominant.epoch") { + if (is.character(epoch) && epoch == "dominant.epoch") { epoch <- count_difftime(tibble::tibble(Datetime = Time.vector))$difftime[1] } # If the user specified an epoch, use that instead @@ -102,6 +105,9 @@ duration_above_threshold <- function(Light.vector, # Return data frame or numeric value if (as.df) { + if(length(threshold) == 2){ + comparison <- "within" + } threshold <- stringr::str_flatten(sort(threshold), collapse = "-") return(tibble::tibble("duration_{comparison}_{threshold}" := tat)) } else { diff --git a/R/metric_frequency_crossing_threshold.R b/R/metric_frequency_crossing_threshold.R new file mode 100644 index 0000000..14f1d03 --- /dev/null +++ b/R/metric_frequency_crossing_threshold.R @@ -0,0 +1,79 @@ +#' Frequency of crossing light threshold +#' +#' This functions calculates the number of times a given threshold +#' light level is crossed. +#' +#' @param Light.vector Numeric vector containing the light data. +#' @param threshold Single numeric value specifying the threshold light level to compare with. +#' @param na.rm Logical. Should missing light values be removed? Defaults to `FALSE`. +#' @param as.df Logical. Should the output be returned as a data frame? If `TRUE`, a data +#' frame with a single column named `frequency_crossing_{threshold}` will be returned. +#' Defaults to `FALSE`. +#' +#' @return Data frame or matrix with pairs of threshold and calculated values. +#' +#' @export +#' +#' @family metrics +#' +#' @references +#' Alvarez, A. A., & Wildsoet, C. F. (2013). Quantifying light +#' exposure patterns in young adult students. \emph{Journal of Modern Optics}, +#' 60(14), 1200–1208. \url{https://doi.org/10.1080/09500340.2013.845700} +#' +#' Hartmeyer, S.L., Andersen, M. (2023). Towards a framework for light-dosimetry studies: +#' Quantification metrics. \emph{Lighting Research & Technology}. +#' \url{https://doi.org/10.1177/14771535231170500} +#' +#' @examples +# +#' N = 60 +#' dataset1 <- +#' tibble::tibble( +#' Id = rep("A", N), +#' Datetime = lubridate::as_datetime(0) + lubridate::minutes(1:N), +#' MEDI = sample(c(sample(1:249, N / 2), sample(250:1000, N / 2))), +#' ) +#' +#' dataset1 %>% +#' dplyr::reframe("Frequency crossing 250lx" = frequency_crossing_threshold(MEDI, threshold = 250)) +#' +#' dataset1 %>% +#' dplyr::reframe(frequency_crossing_threshold(MEDI, threshold = 250, as.df = TRUE)) +#' +frequency_crossing_threshold <- function(Light.vector, + threshold, + na.rm = FALSE, + as.df = FALSE) { + + # Perform argument checks + stopifnot( + "`Light.vector` must be numeric!" = is.numeric(Light.vector), + "`threshold` must be numeric!" = is.numeric(threshold), + "`threshold` must be one value!" = length(threshold) == 1, + "`na.rm` must be logical!" = is.logical(na.rm), + "`as.df` must be logical!" = is.logical(as.df) + ) + + # Remove NAs + if (na.rm) { + Light.vector <- Light.vector[!is.na(Light.vector)] + } + + if (any(is.na(Light.vector))){ + fic <- NA + } + else{ + # Calculate FIC + fic <- sum(abs(diff(compare_threshold(Light.vector, threshold)))) + } + + # Return data frame or numeric value + if (as.df) { + threshold <- as.character(threshold) + return(tibble::tibble("frequency_crossing_{threshold}" := fic)) + } else { + return(fic) + } +} + diff --git a/R/metric_interdaily_stability.R b/R/metric_interdaily_stability.R index b9965b8..8076256 100755 --- a/R/metric_interdaily_stability.R +++ b/R/metric_interdaily_stability.R @@ -10,8 +10,9 @@ #' @param Light.vector Numeric vector containing the light data. #' @param Datetime.vector Vector containing the time data. Must be POSIXct. #' @param na.rm Logical. Should missing values be removed? Defaults to `FALSE`. -#' @param as.df Logical. Should the output be returned as a data frame? Defaults -#' to `FALSE`. +#' @param as.df Logical. Should the output be returned as a data frame? If `TRUE`, a data +#' frame with a single column named `interdaily_stability` will be returned. +#' Defaults to `FALSE`. #' #' @return Numeric value or dataframe with column 'IS'. #' @@ -53,6 +54,8 @@ interdaily_stability <- function(Light.vector, stopifnot( "`Light.vector` must be numeric!" = is.numeric(Light.vector), "`Datetime.vector` must be POSIXct!" = lubridate::is.POSIXct(Datetime.vector), + "`Light.vector` and `Datetime.vector` must be same length!" = + length(Light.vector) == length(Datetime.vector), "`na.rm` must be logical!" = is.logical(na.rm), "`as.df` must be logical!" = is.logical(as.df) ) diff --git a/R/metric_intradaily_variability.R b/R/metric_intradaily_variability.R index ef7dd79..7aab154 100755 --- a/R/metric_intradaily_variability.R +++ b/R/metric_intradaily_variability.R @@ -8,8 +8,9 @@ #' @param Light.vector Numeric vector containing the light data. #' @param Datetime.vector Vector containing the time data. Must be POSIXct. #' @param na.rm Logical. Should missing values be removed? Defaults to `FALSE`. -#' @param as.df Logical. Should the output be returned as a data frame? Defaults -#' to `FALSE`. +#' @param as.df Logical. Should the output be returned as a data frame? If `TRUE`, a data +#' frame with a single column named `intradaily_variability` will be returned. +#' Defaults to `FALSE`. #' #' @return Numeric value or dataframe with column 'IV'. #' @@ -51,6 +52,8 @@ intradaily_variability <- function(Light.vector, stopifnot( "`Light.vector` must be numeric!" = is.numeric(Light.vector), "`Datetime.vector` must be POSIXct!" = lubridate::is.POSIXct(Datetime.vector), + "`Light.vector` and `Datetime.vector` must be same length!" = + length(Light.vector) == length(Datetime.vector), "`na.rm` must be logical!" = is.logical(na.rm), "`as.df` must be logical!" = is.logical(as.df) ) diff --git a/R/metric_midpointCE.R b/R/metric_midpointCE.R index 6496f69..03a1904 100644 --- a/R/metric_midpointCE.R +++ b/R/metric_midpointCE.R @@ -3,13 +3,14 @@ #' This function calculates the timing corresponding to half of the cumulative #' light exposure within the given time series. #' -#' @param Light.vector Numeric vector containing the light data. Missing values are -#' replaced with 0. -#' @param Time.vector Vector containing the time data. Can be numeric, HMS or POSIXct. -#' @param na.rm Logical. Should missing values be removed for the calculation? +#' @param Light.vector Numeric vector containing the light data. +#' @param Time.vector Vector containing the time data. Can be \link[base]{POSIXct}, \link[hms]{hms}, +#' \link[lubridate]{duration}, or \link[base]{difftime}. +#' @param na.rm Logical. Should missing values be removed for the calculation? If `TRUE`, +#' missing values will be replaced by zero. Defaults to `FALSE`. +#' @param as.df Logical. Should the output be returned as a data frame? If `TRUE`, a data +#' frame with a single column named `midpointCE` will be returned. #' Defaults to `FALSE`. -#' @param as.df Logical. Should the output be returned as a data frame? Defaults -#' to `FALSE`. #' #' @return Single column data frame or vector. #' @@ -28,7 +29,7 @@ #' \url{https://doi.org/10.1177/14771535231170500} #' #' @examples -#' # Dataset with POSIXct time vector +# Dataset with POSIXct time vector #' dataset1 <- #' tibble::tibble( #' Id = rep("A", 24), @@ -52,11 +53,11 @@ #' "Midpoint of cmulative exposure" = midpointCE(MEDI, Time) #' ) #' -#' # Dataset with numeric time vector +#' # Dataset with duration time vector #' dataset3 <- #' tibble::tibble( #' Id = rep("A", 24), -#' Hour = 0:23, +#' Hour = lubridate::duration(0:23, "hours"), #' MEDI = c(rep(1, 6), rep(250, 13), rep(1, 5)) #' ) #' dataset3 %>% @@ -72,8 +73,11 @@ midpointCE <- function(Light.vector, # Perform argument checks stopifnot( "`Light.vector` must be numeric!" = is.numeric(Light.vector), - "`Time.vector` must be numeric, HMS, or POSIXct" = - is.numeric(Time.vector) | hms::is_hms(Time.vector) | lubridate::is.POSIXct(Time.vector), + "`Time.vector` must be POSIXct, hms, duration, or difftime!" = + lubridate::is.POSIXct(Time.vector) | hms::is_hms(Time.vector) | + lubridate::is.duration(Time.vector) | lubridate::is.difftime(Time.vector), + "`Light.vector` and `Time.vector` must be same length!" = + length(Light.vector) == length(Time.vector), "`na.rm` must be logical!" = is.logical(na.rm), "`as.df` must be logical!" = is.logical(as.df) ) @@ -83,15 +87,21 @@ midpointCE <- function(Light.vector, Light.vector[is.na(Light.vector)] <- 0 } - # Find midpoint of CE - cumsum <- cumsum(Light.vector) - halfSum <- cumsum[length(cumsum)] / 2 - midpointCE <- which.min(abs(cumsum - halfSum)) + # If any value is NA, return NA + if(any(is.na(Light.vector))){ + midpointCE = convert_to_timescale(NA, Time.vector) + } + else{ + # Find midpoint of CE + cumsum <- cumsum(Light.vector) + halfSum <- cumsum[length(cumsum)] / 2 + midpointCE = Time.vector[which.min(abs(cumsum - halfSum))] + } # Return as data frame or numeric vector if (as.df) { - return(tibble::tibble("midpointCE" = Time.vector[midpointCE])) + return(tibble::tibble("midpointCE" = midpointCE)) } else { - return(Time.vector[midpointCE]) + return(midpointCE) } } diff --git a/R/metric_nvR.R b/R/metric_nvR.R index 588d47a..0796fb8 100644 --- a/R/metric_nvR.R +++ b/R/metric_nvR.R @@ -10,11 +10,12 @@ #' #' @param MEDI.vector Numeric vector containing the melanopic EDI data. #' @param Illuminance.vector Numeric vector containing the Illuminance data. -#' @param Time.vector Vector containing the time data. Can be numeric, HMS or POSIXct. +#' @param Time.vector Vector containing the time data. Can be \link[base]{POSIXct}, \link[hms]{hms}, +#' \link[lubridate]{duration}, or \link[base]{difftime}. #' @param epoch The epoch at which the data was sampled. Can be either a -#' `lubridate::duration()` or a string. If it is a string, it needs to be +#' \link[lubridate]{duration} or a string. If it is a string, it needs to be #' either `"dominant.epoch"` (the default) for a guess based on the data, or a valid -#' `lubridate::duration()` string, e.g., `"1 day"` or `"10 sec"`. +#' \link[lubridate]{duration} string, e.g., `"1 day"` or `"10 sec"`. #' #' @return A numeric vector containing the nvRD data. The output has the same #' length as `Time.vector`. @@ -119,8 +120,11 @@ nvRD <- function(MEDI.vector, stopifnot( "`MEDI.vector` must be numeric!" = is.numeric(MEDI.vector), "`Illuminance.vector` must be numeric!" = is.numeric(Illuminance.vector), - "`Time.vector` must be numeric, HMS, or POSIXct" = - is.numeric(Time.vector) | hms::is_hms(Time.vector) | lubridate::is.POSIXct(Time.vector), + "`Time.vector` must be POSIXct, hms, duration, or difftime!" = + lubridate::is.POSIXct(Time.vector) | hms::is_hms(Time.vector) | + lubridate::is.duration(Time.vector) | lubridate::is.difftime(Time.vector), + "`MEDI.vector` and `Time.vector` must be same length!" = + length(MEDI.vector) == length(Time.vector), "`epoch` must either be a duration or a string" = lubridate::is.duration(epoch) | is.character(epoch) ) @@ -140,7 +144,7 @@ nvRD <- function(MEDI.vector, slope <- 3.55 # Get the epochs based on the data - if (epoch == "dominant.epoch") { + if (is.character(epoch) && epoch == "dominant.epoch") { epoch <- count_difftime(tibble::tibble(Datetime = Time.vector))$difftime[1] } # If the user specified an epoch, use that instead @@ -161,7 +165,7 @@ nvRD <- function(MEDI.vector, dFL1 <- ifelse(delta < 0.3, round(0.3 / delta), 1) # Filter LH - dFLH <- ifelse(delta < 0.7, round(0.7 / delta), 1) + dFLH <- ifelse(delta < 1.7, round(1.7 / delta), 1) # MODEL u <- nvR_filterSMA(dFL1, Ieff) @@ -179,13 +183,14 @@ nvRD <- function(MEDI.vector, #' #' @param nvRD Numeric vector containing the non-visual direct response. #' See \code{\link{nvRD}}. -#' @param Time.vector Vector containing the time data. Can be numeric, HMS or POSIXct. +#' @param Time.vector Vector containing the time data. Can be \link[base]{POSIXct}, \link[hms]{hms}, +#' \link[lubridate]{duration}, or \link[base]{difftime}. #' @param epoch The epoch at which the data was sampled. Can be either a -#' `lubridate::duration()` or a string. If it is a string, it needs to be +#' \link[lubridate]{duration} or a string. If it is a string, it needs to be #' either `"dominant.epoch"` (the default) for a guess based on the data, or a valid -#' `lubridate::duration()` string, e.g., `"1 day"` or `"10 sec"`. +#' \link[lubridate]{duration} string, e.g., `"1 day"` or `"10 sec"`. #' @param as.df Logical. Should a data frame with be returned? If `TRUE`, a data -#' frame with a single column named `TAT_{threshold}` will be returned. +#' frame with a single column named `nvRD_cumulative` will be returned. #' Defaults to `FALSE`. #' #' @return A numeric value or single column data frame. @@ -222,15 +227,18 @@ nvRD_cumulative_response <- function(nvRD, # Perform argument checks stopifnot( "`nvRD` must be numeric!" = is.numeric(nvRD), - "`Time.vector` must be numeric, HMS, or POSIXct" = - is.numeric(Time.vector) | hms::is_hms(Time.vector) | lubridate::is.POSIXct(Time.vector), + "`Time.vector` must be POSIXct, hms, duration, or difftime!" = + lubridate::is.POSIXct(Time.vector) | hms::is_hms(Time.vector) | + lubridate::is.duration(Time.vector) | lubridate::is.difftime(Time.vector), + "`nvRD` and `Time.vector` must be same length!" = + length(nvRD) == length(Time.vector), "`epoch` must either be a duration or a string" = lubridate::is.duration(epoch) | is.character(epoch), "`as.df` must be logical!" = is.logical(as.df) ) # Get the epochs based on the data - if (epoch == "dominant.epoch") { + if (is.character(epoch) && epoch == "dominant.epoch") { epoch <- count_difftime(tibble::tibble(Datetime = Time.vector))$difftime[1] } # If the user specified an epoch, use that instead @@ -260,11 +268,12 @@ nvRD_cumulative_response <- function(nvRD, #' #' @param MEDI.vector Numeric vector containing the melanopic EDI data. #' @param Illuminance.vector Numeric vector containing the Illuminance data. -#' @param Time.vector Vector containing the time data. Can be numeric, HMS or POSIXct. +#' @param Time.vector Vector containing the time data. Can be \link[base]{POSIXct}, +#' \link[hms]{hms}, \link[lubridate]{duration}, or \link[base]{difftime}. #' @param epoch The epoch at which the data was sampled. Can be either a -#' `lubridate::duration()` or a string. If it is a string, it needs to be +#' \link[lubridate]{duration} or a string. If it is a string, it needs to be #' either `"dominant.epoch"` (the default) for a guess based on the data, or a valid -#' `lubridate::duration()` string, e.g., `"1 day"` or `"10 sec"`. +#' \link[lubridate]{duration} string, e.g., `"1 day"` or `"10 sec"`. #' @param sleep.onset The time of habitual sleep onset. Can be HMS, numeric, or NULL. #' If NULL (the default), then the data is assumed to start at habitual sleep onset. #' If `Time.vector` is HMS or POSIXct, `sleep.onset` must be HMS. Likewise, if @@ -353,6 +362,7 @@ nvRD_cumulative_response <- function(nvRD, #' \url{http://dx.doi.org/10.5075/epfl-thesis-7146} #' #' @examples +#' #' dataset1 <- #' tibble::tibble( #' Id = rep("B", 60 * 48), @@ -361,12 +371,20 @@ nvRD_cumulative_response <- function(nvRD, #' rep(0, 60*8), rep(sample(1:1000, 16, replace = TRUE), each = 60)), #' MEDI = Illuminance * rep(sample(0.5:1.5, 48, replace = TRUE), each = 60) #' ) -#' +#' # Time.vector as POSIXct #' dataset1.nvRC <- dataset1 %>% #' dplyr::mutate( #' nvRC = nvRC(MEDI, Illuminance, Datetime, sleep.onset = hms::as_hms("22:00:00")) #' ) #' +#' # Time.vector as difftime +#' dataset2 <- dataset1 %>% +#' dplyr::mutate(Datetime = Datetime - lubridate::as_datetime(lubridate::dhours(22))) +#' dataset2.nvRC <- dataset2 %>% +#' dplyr::mutate( +#' nvRC = nvRC(MEDI, Illuminance, Datetime, sleep.onset = lubridate::dhours(0)) +#' ) +#' nvRC <- function(MEDI.vector, Illuminance.vector, Time.vector, @@ -377,11 +395,16 @@ nvRC <- function(MEDI.vector, stopifnot( "`MEDI.vector` must be numeric!" = is.numeric(MEDI.vector), "`Illuminance.vector` must be numeric!" = is.numeric(Illuminance.vector), - "`Time.vector` must be numeric, HMS, or POSIXct" = - is.numeric(Time.vector) | hms::is_hms(Time.vector) | lubridate::is.POSIXct(Time.vector), + "`Time.vector` must be POSIXct, hms, duration, or difftime!" = + lubridate::is.POSIXct(Time.vector) | hms::is_hms(Time.vector) | + lubridate::is.duration(Time.vector) | lubridate::is.difftime(Time.vector), + "`MEDI.vector` and `Time.vector` must be same length!" = + length(MEDI.vector) == length(Time.vector), "`epoch` must either be a duration or a string" = lubridate::is.duration(epoch) | is.character(epoch), - "`sleep.onset` must be numeric, HMS, or NULL" = is.null(sleep.onset) | hms::is_hms(sleep.onset) + "`sleep.onset` must be hms, duration, difftime, or NULL" = + is.null(sleep.onset) | hms::is_hms(sleep.onset) | + lubridate::is.duration(sleep.onset) | lubridate::is.difftime(sleep.onset) ) # Replace missing values @@ -399,7 +422,7 @@ nvRC <- function(MEDI.vector, slope <- 1.42 # Get the epochs based on the data - if (epoch == "dominant.epoch") { + if (is.character(epoch) && epoch == "dominant.epoch") { epoch <- count_difftime(tibble::tibble(Datetime = Time.vector))$difftime[1] } # If the user specified an epoch, use that instead @@ -413,36 +436,35 @@ nvRC <- function(MEDI.vector, # Convert to effective irradiance Ieff <- nvR_getIeff(MEDI.vector, Illuminance.vector, delta) - if(is.numeric(Time.vector) & is.numeric(sleep.onset)) { - bt <- sleep.onset - dt <- Time.vector + if(!is.null(sleep.onset)){ + bt = as.numeric(sleep.onset) - t <- (dt - bt) - } - else{ - # Align circadian sensitivity curve with sleep onset - if (hms::is_hms(sleep.onset)) { - # Check that Time.vector is HMS or POSIXct + if(hms::is.hms(sleep.onset)){ stopifnot("`sleep.onset` is HMS but `Time.vector` is not HMS or POSIXct" = hms::is_hms(Time.vector) | lubridate::is.POSIXct(Time.vector)) - - bt <- as.numeric(sleep.onset) + # Normalise datetime vector dt <- as.numeric(Time.vector) - as.numeric(Time.vector)[1] + as.numeric(hms::as_hms(Time.vector[1])) - - # Seconds to hours - t = (dt - bt) / 3600 - - } else { - # Check whether first three hours are darkness --> timeseries should start at - # habitual sleep onset. - if (mean(Illuminance.vector[1:(3 / delta)]) > 10) { - warning("Average Illuminance.vector across the first three hours is higher than 10lx. Does the timeseries really start at habitual sleep onset?") - } - t <- seq(0, (length(Ieff) * delta) - delta, delta) } + else{ + stopifnot("`sleep.onset` is duration or difftime but `Time.vector` is not duration or difftime" = + lubridate::is.duration(Time.vector) | lubridate::is.difftime(Time.vector)) + dt <- as.numeric(Time.vector) + } + + # Seconds to hours + t = (dt - bt) / 3600 + } + else { + # Check whether first three hours are darkness --> timeseries should start at + # habitual sleep onset. + if (mean(Illuminance.vector[1:(3 / delta)]) > 10) { + warning("Average Illuminance.vector across the first three hours is higher than 10lx. Does the timeseries really start at habitual sleep onset?") + } + t <- seq(0, (length(Ieff) * delta) - delta, delta) } + # Get circadian modulation Rmax <- nvR_circadianModulator(t) diff --git a/R/metric_period_above_threshold.R b/R/metric_period_above_threshold.R new file mode 100644 index 0000000..75b56ee --- /dev/null +++ b/R/metric_period_above_threshold.R @@ -0,0 +1,141 @@ +#' Length of longest continuous period above/below threshold +#' +#' This function finds the length of the longest continous period above/below +#' a specified threshold light level or within a specified range of light levels. +#' +#' @param Light.vector Numeric vector containing the light data. +#' @param Time.vector Vector containing the time data. Can be \link[base]{POSIXct}, +#' \link[hms]{hms}, \link[lubridate]{duration}, or \link[base]{difftime}. +#' @param comparison String specifying whether the period of light levels above or +#' below threshold should be calculated. Can be either `"above"` (the default) +#' or `"below"`. If two values are provided for `threshold`, this argument will be ignored. +#' @param threshold Single numeric value or two numeric values specifying the +#' threshold light level(s) to compare with. If a vector with two values is provided, +#' the period of light levels within the two thresholds will be calculated. +#' @param epoch The epoch at which the data was sampled. Can be either a +#' \link[lubridate]{duration} or a string. If it is a string, it needs to be +#' either `"dominant.epoch"` (the default) for a guess based on the data, or a valid +#' \link[lubridate]{duration} string, e.g., `"1 day"` or `"10 sec"`. +#' @param loop Logical. Should the data be looped? Defaults to `FALSE`. +#' @param na.replace Logical. Should missing values (NA) be replaced +#' for the calculation? If `TRUE` missing values will not be removed but will +#' result in `FALSE` when comparing `Light.vector` with `threshold`. +#' Defaults to `FALSE`. +#' @param na.rm Logical. Should missing values (NA) be removed for the calculation? +#' If `TRUE`, this argument will override `na.replace`. Defaults to `FALSE`. +#' @param as.df Logical. Should a data frame be returned? If `TRUE`, a data +#' frame with a single column named `period_{comparison}_{threshold}` will be returned. +#' Defaults to `FALSE`. +#' +#' @return A duration object (see \code{\link[lubridate]{duration}}) as single value, +#' or single column data frame. +#' +#' @export +#' +#' @family metrics +#' +#' @examples +#' +#' N <- 60 +#' # Dataset with continous period of >250lx for 35min +#' dataset1 <- +#' tibble::tibble( +#' Id = rep("A", N), +#' Datetime = lubridate::as_datetime(0) + lubridate::minutes(1:N), +#' MEDI = c(sample(1:249, N-35, replace = TRUE), +#' sample(250:1000, 35, replace = TRUE)) +#' ) +#' +#' dataset1 %>% +#' dplyr::reframe("Period >250lx" = period_above_threshold(MEDI, Datetime, threshold = 250)) +#' +#' dataset1 %>% +#' dplyr::reframe("Period <250lx" = period_above_threshold(MEDI, Datetime, "below", threshold = 250)) +#' +#' # Dataset with continous period of 100-250lx for 20min +#' dataset2 <- +#' tibble::tibble( +#' Id = rep("B", N), +#' Datetime = lubridate::as_datetime(0) + lubridate::minutes(1:N), +#' MEDI = c(sample(c(1:99, 251-1000), N-20, replace = TRUE), +#' sample(100:250, 20, replace = TRUE)), +#' ) +#' dataset2 %>% +#' dplyr::reframe("Period 250lx" = period_above_threshold(MEDI, Datetime, threshold = c(100,250))) +#' +#' # Return data frame +#' dataset1 %>% +#' dplyr::reframe(period_above_threshold(MEDI, Datetime, threshold = 250, as.df = TRUE)) +#' +period_above_threshold <- function(Light.vector, + Time.vector, + comparison = c("above", "below"), + threshold, + epoch = "dominant.epoch", + loop = FALSE, + na.replace = FALSE, + na.rm = FALSE, + as.df = FALSE) { + + # Match input arguments + comparison <- match.arg(comparison) + + # Perform argument checks + stopifnot( + "`Light.vector` must be numeric!" = is.numeric(Light.vector), + "`Time.vector` must be POSIXct, hms, duration, or difftime!" = + lubridate::is.POSIXct(Time.vector) | hms::is_hms(Time.vector) | + lubridate::is.duration(Time.vector) | lubridate::is.difftime(Time.vector), + "`Light.vector` and `Time.vector` must be same length!" = + length(Light.vector) == length(Time.vector), + "`threshold` must be numeric!" = is.numeric(threshold), + "`threshold` must be either one or two values!" = length(threshold) %in% c(1, 2), + "`epoch` must either be a duration or a string" = + lubridate::is.duration(epoch) | is.character(epoch), + "`loop` must be logical!" = is.logical(loop), + "`na.rm` must be logical!" = is.logical(na.rm), + "`as.df` must be logical!" = is.logical(as.df) + ) + + # Get the epochs based on the data + if (is.character(epoch) && epoch == "dominant.epoch") { + epoch <- count_difftime(tibble::tibble(Datetime = Time.vector))$difftime[1] + } + # If the user specified an epoch, use that instead + else { + epoch <- lubridate::as.duration(epoch) + } + + # Loop data + if (loop) { + Light.vector <- c(Light.vector, Light.vector) + } + + # Remove NAs + if(na.rm){ + Light.vector <- Light.vector[!is.na(Light.vector)] + } + + # Find longest period + longest_period <- function(x){ + z <- c(x, 0) + z <- (cumsum(z) * c(diff(z) < 0, 0)) + max(diff(c(0, 0, z[z != 0]))) + } + out <- compare_threshold(Light.vector, threshold, comparison, na.replace) %>% + longest_period() + + # As duration object + out <- lubridate::as.duration(out * epoch) + + # Return data frame or numeric value + if (as.df) { + if(length(threshold) == 2){ + comparison <- "within" + } + threshold <- stringr::str_flatten(sort(threshold), collapse = "-") + return(tibble::tibble("period_{comparison}_{threshold}" := out)) + } else { + return(out) + } +} diff --git a/R/metric_pulses_above_threshold.R b/R/metric_pulses_above_threshold.R new file mode 100644 index 0000000..6cb8dca --- /dev/null +++ b/R/metric_pulses_above_threshold.R @@ -0,0 +1,210 @@ +#' Pulses above threshold +#' +#' This function clusters the light data into continuous clusters (pulses) of +#' light above/below a given threshold. Clustering may be fine-tuned by setting +#' the minimum length of the clusters and by allowing brief interruptions to be +#' included in a single cluster, with a specified maximum length of interruption +#' episodes and proportion of total amount of interruptions to light above +#' threshold. +#' +#' @param Light.vector Numeric vector containing the light data. Missing values will +#' be considered as `FALSE` when comparing light levels against the threshold. +#' @param Time.vector Vector containing the time data. Can be \link[base]{POSIXct}, +#' \link[hms]{hms}, \link[lubridate]{duration}, or \link[base]{difftime}. +#' @param comparison String specifying whether the time above or below threshold +#' should be calculated. Can be either `"above"` (the default) or `"below"`. If +#' two values are provided for `threshold`, this argument will be ignored. +#' @param threshold Single numeric value or two numeric values specifying the +#' threshold light level(s) to compare with. If a vector with two values is provided, +#' the timing corresponding to light levels between the two thresholds will be +#' calculated. +#' @param min.length The minimum length of a pulse. Can be either a +#' \link[lubridate]{duration} or a string. If it is a string, it needs to be a valid +#' \link[lubridate]{duration} string, e.g., `"1 day"` or `"10 sec"`. Defaults to +#' `"8 mins"` as in Wilson et al. (2018). +#' @param max.interrupt Maximum length of each episode of interruptions. Can be either a +#' \link[lubridate]{duration} or a string. If it is a string, it needs to be a valid +#' \link[lubridate]{duration} string, e.g., `"1 day"` or `"10 sec"`. Defaults to +#' `"2 mins"` as in Wilson et al. (2018). +#' @param prop.interrupt Numeric value between `0` and `1` specifying the +#' maximum proportion of the total number of interruptions. Defaults to `0.25` +#' as in Wilson et al. (2018). +#' @param epoch The epoch at which the data was sampled. Can be either a +#' \link[lubridate]{duration} or a string. If it is a string, it needs to be +#' either `"dominant.epoch"` (the default) for a guess based on the data, or a valid +#' \link[lubridate]{duration} string, e.g., `"1 day"` or `"10 sec"`. +#' @param return.indices Logical. Should the cluster indices be returned? Only works if +#' `as.df` is `FALSE`. Defaults to `FALSE`. +#' @param na.rm Logical. Should missing values be removed for the calculation of +#' pulse metrics? Defaults to `FALSE`. +#' @param as.df Logical. Should a data frame be returned? If `TRUE`, a data +#' frame with seven columns ("n", "mean_level", "mean_duration", "total_duration", +#' "mean_onset", "mean_midpoint", "mean_offset") and the threshold (e.g., `_{threshold}`) +#' will be returned. Defaults to `FALSE`. +#' +#' @return List or data frame with calculated values. +#' +#' @details The timeseries is assumed to be regular. Missing values in the +#' light data will be replaced by 0. +#' +#' @export +#' +#' @family metrics +#' +#' @references Wilson, J., Reid, K. J., Braun, R. I., Abbott, S. M., & Zee, P. C. +#' (2018). Habitual light exposure relative to circadian timing in delayed +#' sleep-wake phase disorder. \emph{Sleep}, 41(11). +#' \url{https://doi.org/10.1093/sleep/zsy166} +#' +#' @examples +#' # Sample data +#' data = sample.data.environment %>% +#' dplyr::filter(Id == "Participant") %>% +#' filter_Datetime(length = lubridate::days(1)) %>% +#' dplyr::mutate( +#' Time = hms::as_hms(Datetime), +#' ) +#' +#' # Time vector as datetime +#' data %>% +#' dplyr::reframe(pulses_above_threshold(MEDI, Datetime, threshold = 250, as.df = TRUE)) +#' +#' # Time vector as hms time +#' data %>% +#' dplyr::reframe(pulses_above_threshold(MEDI, Time, threshold = 250, as.df = TRUE)) +#' +#' # Pulses below threshold +#' data %>% +#' dplyr::reframe(pulses_above_threshold(MEDI, Datetime, "below", threshold = 250, as.df = TRUE)) +#' +#' # Pulses within threshold range +#' data %>% +#' dplyr::reframe(pulses_above_threshold(MEDI, Datetime, threshold = c(250,1000), as.df = TRUE)) +#' +pulses_above_threshold <- function(Light.vector, + Time.vector, + comparison = c("above", "below"), + threshold, + min.length = "8 mins", + max.interrupt = "2 mins", + prop.interrupt = 0.25, + epoch = "dominant.epoch", + return.indices = FALSE, + na.rm = FALSE, + as.df = FALSE) { + + # Match input arguments + comparison <- match.arg(comparison) + + # Perform argument checks + stopifnot( + "`Light.vector` must be numeric!" = is.numeric(Light.vector), + "`Time.vector` must be POSIXct, hms, duration, or difftime!" = + lubridate::is.POSIXct(Time.vector) | hms::is_hms(Time.vector) | + lubridate::is.duration(Time.vector) | lubridate::is.difftime(Time.vector), + "`Light.vector` and `Time.vector` must be same length!" = + length(Light.vector) == length(Time.vector), + "`threshold` must be numeric!" = is.numeric(threshold), + "`threshold` must be either one or two values!" = length(threshold) %in% c(1, 2), + "`min.length` must either be a duration or a valid duration string" = + lubridate::is.duration(min.length) | is.character(min.length), + "`max.interrupt` must either be a duration or a valid duration string" = + lubridate::is.duration(max.interrupt) | is.character(max.interrupt), + "`prop.interrupt` must be numeric" = is.numeric(prop.interrupt), + "`prop.interrupt` must be between 0 and 1" = prop.interrupt >= 0 & prop.interrupt <= 1, + "`epoch` must either be a duration or a string" = + lubridate::is.duration(epoch) | is.character(epoch), + "`na.rm` must be logical!" = is.logical(na.rm), + "`as.df` must be logical!" = is.logical(as.df) + ) + + # Get the epochs based on the data + if (is.character(epoch) && epoch == "dominant.epoch") { + epoch <- count_difftime(tibble::tibble(Datetime = Time.vector))$difftime[1] + } + # If the user specified an epoch, use that instead + else { + epoch <- lubridate::as.duration(epoch) + } + + # Parse duration inputs + min.length <- lubridate::as.duration(min.length) + max.interrupt <- lubridate::as.duration(max.interrupt) + + # Additional checks + stopifnot( + "Time parameters must be equal to or longer than the epoch." = + any(c(min.length, max.interrupt) >= epoch) + ) + + # Convert to sample counts + min.length <- round(min.length / epoch) + max.interrupt <- round(max.interrupt / epoch) + + # Find pulses + pulses <- find_clusters( + compare_threshold(Light.vector, threshold, comparison), + min.length, max.interrupt, prop.interrupt, "pulse" + ) + + # Summarise pulse metrics + options(dplyr.summarise.inform = FALSE) + data.pulses <- + tibble::tibble( + row_idx = 1:length(Light.vector), + light = Light.vector, + time = as.numeric(Time.vector) + ) %>% + dplyr::left_join(pulses, by = "row_idx") %>% + dplyr::filter(is_pulse) %>% + + # Summarise per pulse + dplyr::group_by(pulse_idx) %>% + dplyr::summarise( + level = mean(light, na.rm = na.rm), + duration = dplyr::n()*epoch, + onset = dplyr::first(time), + offset = dplyr::last(time), + midpoint = mean(time, na.rm = na.rm) + ) %>% + dplyr::ungroup() %>% + + # Summarise across pulses + dplyr::summarise( + n = dplyr::n(), + mean_level = mean(level), + mean_duration = mean(duration) %>% round() %>% lubridate::as.duration(), + total_duration = mean_duration*n, + mean_onset = mean(onset) %>% round(), + mean_midpoint = mean(midpoint) %>% round(), + mean_offset = mean(offset) %>% round(), + ) %>% + + # Convert timing metrics to corresponding time scale + dplyr::mutate( + dplyr::across( + c(mean_onset, mean_midpoint, mean_offset), + ~convert_to_timescale(.x, Time.vector) + ) + ) + + # Return data frame or list + if (as.df) { + if(length(threshold) == 2){ + comparison <- "within" + } + threshold <- stringr::str_flatten(sort(threshold), collapse = "-") + data.pulses <- data.pulses %>% + dplyr::rename_with(~paste(.x, "pulses", comparison, threshold, sep = "_")) + if(return.indices){ + warning("`return.indices`is `TRUE` but indices cannot be returned if `as.df` is `TRUE`. Please use `as.df = FALSE`") + } + } + else{ + data.pulses <- data.pulses %>% as.list() + if(return.indices){ + data.pulses$pulse_indices <- pulses$row_idx + } + } + return(data.pulses) +} diff --git a/R/metric_threshold_for_duration.R b/R/metric_threshold_for_duration.R new file mode 100644 index 0000000..0dc677a --- /dev/null +++ b/R/metric_threshold_for_duration.R @@ -0,0 +1,116 @@ +#' Find threshold for given duration +#' +#' This function finds the threshold for which light levels are above/below for +#' a given duration. This function can be considered as the inverse of +#' \code{\link{duration_above_threshold}}. +#' +#' @param Light.vector Numeric vector containing the light data. +#' @param Time.vector Vector containing the time data. Can be \link[base]{POSIXct}, +#' \link[hms]{hms}, \link[lubridate]{duration}, or \link[base]{difftime}. +#' @param duration The duration for which the threshold should be found. Can be either a +#' \link[lubridate]{duration} or a string. If it is a string, it needs to be a valid +#' \link[lubridate]{duration} string, e.g., `"1 day"` or `"10 sec"`. +#' @param comparison String specifying whether light levels above or below the threshold +#' should be considered. Can be either `"above"` (the default) or `"below"`. +#' @param epoch The epoch at which the data was sampled. Can be either a +#' \link[lubridate]{duration} or a string. If it is a string, it needs to be +#' either `"dominant.epoch"` (the default) for a guess based on the data, or a valid +#' \link[lubridate]{duration} string, e.g., `"1 day"` or `"10 sec"`. +#' @param na.rm Logical. Should missing values (NA) be removed for the calculation? +#' Defaults to `FALSE`. +#' @param as.df Logical. Should a data frame with be returned? If `TRUE`, a data +#' frame with a single column named `threshold_{comparison}_for_{duration}` will be returned. +#' Defaults to `FALSE`. +#' +#' @return Single numeric value or single column data frame. +#' +#' @export +#' +#' @family metrics +#' +#' @examples +#' N <- 60 +#' # Dataset with 30 min < 250lx and 30min > 250lx +#' dataset1 <- +#' tibble::tibble( +#' Id = rep("A", N), +#' Datetime = lubridate::as_datetime(0) + lubridate::minutes(1:N), +#' MEDI = sample(c(sample(1:249, N / 2, replace = TRUE), +#' sample(250:1000, N / 2, replace = TRUE))), +#' ) +#' +#' dataset1 %>% +#' dplyr::reframe("Threshold above which for 30 mins" = +#' threshold_for_duration(MEDI, Datetime, duration = "30 mins")) +#' +#' dataset1 %>% +#' dplyr::reframe("Threshold below which for 30 mins" = +#' threshold_for_duration(MEDI, Datetime, duration = "30 mins", +#' comparison = "below")) +#' +#' dataset1 %>% +#' dplyr::reframe(threshold_for_duration(MEDI, Datetime, duration = "30 mins", +#' as.df = TRUE)) +#' +threshold_for_duration <- function(Light.vector, + Time.vector, + duration, + comparison = c("above", "below"), + epoch = "dominant.epoch", + na.rm = FALSE, + as.df = FALSE) { + + + # Match input arguments + comparison <- match.arg(comparison) + + # Perform argument checks + stopifnot( + "`Light.vector` must be numeric!" = is.numeric(Light.vector), + "`Time.vector` must be POSIXct, hms, duration, or difftime!" = + lubridate::is.POSIXct(Time.vector) | hms::is_hms(Time.vector) | + lubridate::is.duration(Time.vector) | lubridate::is.difftime(Time.vector), + "`Light.vector` and `Time.vector` must be same length!" = + length(Light.vector) == length(Time.vector), + "`duration` must either be duration or a string!" = + lubridate::is.duration(duration) | is.character(duration), + "`epoch` must either be a duration or a string" = + lubridate::is.duration(epoch) | is.character(epoch), + "`as.df` must be logical!" = is.logical(as.df) + ) + + # Get the epochs based on the data + if (is.character(epoch) && epoch == "dominant.epoch") { + epoch <- count_difftime(tibble::tibble(Datetime = Time.vector))$difftime[1] + } + # If the user specified an epoch, use that instead + else { + epoch <- lubridate::as.duration(epoch) + } + + # Duration parameter as duration object + duration <- lubridate::as.duration(duration) + + if(na.rm){ + Light.vector = Light.vector[!is.na(Light.vector)] + } + + # Find threshold for given duration + idx = floor(duration / epoch) + sorted <- sort(Light.vector, decreasing = (comparison == "above")) + threshold <- sorted[idx] + + # Return NA if missing values present + if(any(is.na(Light.vector))){ + threshold <- as.double(NA) + } + + # Return data frame or numeric value + if (as.df) { + duration <- stringr::str_extract(as.character(duration), "~.*") %>% + stringr::str_remove("~") %>% stringr::str_remove("\\)") %>% stringr::str_replace_all(" ", "_") + return(tibble::tibble("threshold_{comparison}_for_{duration}" := threshold)) + } else { + return(threshold) + } +} diff --git a/R/metric_timing_above_threshold.R b/R/metric_timing_above_threshold.R index e2d5040..9d1b72e 100755 --- a/R/metric_timing_above_threshold.R +++ b/R/metric_timing_above_threshold.R @@ -5,7 +5,8 @@ #' time interval. #' #' @param Light.vector Numeric vector containing the light data. -#' @param Time.vector Vector containing the time data. Can be numeric, HMS or POSIXct. +#' @param Time.vector Vector containing the time data. Can be \link[base]{POSIXct}, +#' \link[hms]{hms}, \link[lubridate]{duration}, or \link[base]{difftime}. #' @param comparison String specifying whether the time above or below threshold #' should be calculated. Can be either `"above"` (the default) or `"below"`. If #' two values are provided for `threshold`, this argument will be ignored. @@ -15,7 +16,7 @@ #' calculated. #' @param na.rm Logical. Should missing values be removed for the calculation? #' Defaults to `FALSE`. -#' @param as.df Logical. Should a data frame with be returned? If `TRUE`, a data +#' @param as.df Logical. Should a data frame be returned? If `TRUE`, a data #' frame with three columns (MLiT, FLiT, LLiT) and the threshold (e.g., `MLiT_{threshold}`) #' will be returned. Defaults to `FALSE`. #' @@ -71,8 +72,11 @@ timing_above_threshold <- function(Light.vector, # Perform argument checks stopifnot( "`Light.vector` must be numeric!" = is.numeric(Light.vector), - "`Time.vector` must be numeric, HMS, or POSIXct" = - is.numeric(Time.vector) | hms::is_hms(Time.vector) | lubridate::is.POSIXct(Time.vector), + "`Time.vector` must be POSIXct, hms, duration, or difftime!" = + lubridate::is.POSIXct(Time.vector) | hms::is_hms(Time.vector) | + lubridate::is.duration(Time.vector) | lubridate::is.difftime(Time.vector), + "`Light.vector` and `Time.vector` must be same length!" = + length(Light.vector) == length(Time.vector), "`threshold` must be numeric!" = is.numeric(threshold), "`threshold` must be either one or two values!" = length(threshold) %in% c(1, 2), "`na.rm` must be logical!" = is.logical(na.rm), @@ -85,14 +89,8 @@ timing_above_threshold <- function(Light.vector, flit = t %>% dplyr::first() llit = t %>% dplyr::last() - # Convert to HMS - if(hms::is_hms(Time.vector)) { - mlit <- mlit %>% hms::as_hms() - } - if(lubridate::is.POSIXct(Time.vector)){ - mlit <- mlit %>% - lubridate::as_datetime(tz = lubridate::tz(Time.vector)) - } + # Convert to corresponding time scale + mlit <- mlit %>% convert_to_timescale(Time.vector) # Prepare output out <- list( @@ -103,6 +101,9 @@ timing_above_threshold <- function(Light.vector, # Return data frame or list if (as.df) { + if(length(threshold) == 2){ + comparison <- "within" + } threshold <- stringr::str_flatten(sort(threshold), collapse = "-") out <- tibble::as_tibble(out) %>% dplyr::rename_with(~paste(.x, "timing", comparison, threshold, sep = "_")) diff --git a/_pkgdown.yml b/_pkgdown.yml index e7f1f60..28610f4 100644 --- a/_pkgdown.yml +++ b/_pkgdown.yml @@ -71,10 +71,13 @@ reference: desc: > This section includes functions to calculate light exposure metrics. contents: + - barroso_lighting_metrics - bright_dark_period - centroidLE - disparity_index - duration_above_threshold + - exponential_moving_average + - frequency_crossing_threshold - intradaily_variability - interdaily_stability - midpointCE @@ -84,8 +87,12 @@ reference: - nvRC_relativeAmplitudeError - nvRD - nvRD_cumulative_response + - period_above_threshold + - pulses_above_threshold + - threshold_for_duration - timing_above_threshold + - title: Helpers desc: > This section includes helper functions that are used in the other sections. diff --git a/man/barroso_lighting_metrics.Rd b/man/barroso_lighting_metrics.Rd new file mode 100644 index 0000000..e5beb58 --- /dev/null +++ b/man/barroso_lighting_metrics.Rd @@ -0,0 +1,88 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/metric_barroso.R +\name{barroso_lighting_metrics} +\alias{barroso_lighting_metrics} +\title{Circadian lighting metrics from Barroso et al. (2014)} +\usage{ +barroso_lighting_metrics( + Light.vector, + Time.vector, + epoch = "dominant.epoch", + loop = FALSE, + na.rm = FALSE, + as.df = FALSE +) +} +\arguments{ +\item{Light.vector}{Numeric vector containing the light data.} + +\item{Time.vector}{Vector containing the time data. Can be \link[base]{POSIXct}, \link[hms]{hms}, +\link[lubridate]{duration}, or \link[base]{difftime}.} + +\item{epoch}{The epoch at which the data was sampled. Can be either a +\link[lubridate]{duration} or a string. If it is a string, it needs to be +either \code{"dominant.epoch"} (the default) for a guess based on the data, or a valid +\link[lubridate]{duration} string, e.g., \code{"1 day"} or \code{"10 sec"}.} + +\item{loop}{Logical. Should the data be looped? Defaults to \code{FALSE}.} + +\item{na.rm}{Logical. Should missing values (NA) be removed for the calculation? +Defaults to \code{FALSE}. If \code{TRUE}, for the calculation of \code{bright_cluster} and +\code{dark_cluster}, missing values will be replaced by 0 +(see \code{\link{period_above_threshold}}).} + +\item{as.df}{Logical. Should a data frame be returned? If \code{TRUE}, a data +frame with seven columns will be returned. Defaults to \code{FALSE}.} +} +\value{ +List or dataframe with the seven values: \code{bright_threshold}, \code{dark_threshold}, +\code{bright_mean_level}, \code{dark_mean_level}, \code{bright_cluster}, \code{dark_cluster}, +\code{circadian_variation}. The output type of \code{bright_cluster}, \code{dark_cluster}, +is a \link[lubridate]{duration} object. +} +\description{ +This function calculates the metrics proposed by Barroso et al. (2014) +for light-dosimetry in the context of research on the non-visual effects of light. +The following metrics are calculated: +} +\details{ +\describe{ +\item{\code{bright_threshold}}{The maximum light intensity for which at least six +hours of measurements are at the same or higher level.} +\item{\code{dark_threshold}}{The minimum light intensity for which at least eight +hours of measurements are at the same or lower level.} +\item{\code{bright_mean_level}}{The 20\% trimmed mean of all light intensity measurements +equal or above the \code{bright_threshold}.} +\item{\code{dark_mean_level}}{The 20\% trimmed mean of all light intensity measurements +equal or below the \code{dark_threshold}.} +\item{\code{bright_cluster}}{The longest continuous time interval above the \code{bright_threshold}.} +\item{\code{dark_cluster}}{The longest continuous time interval below the \code{dark_threshold}.} +\item{\code{circadian_variation}}{A measure of periodicity of the daily lighting +schedule over a given set of days. Calculated as the coefficient of variation +of input light data. +} +} +} +\examples{ + +dataset1 <- + tibble::tibble( + Id = rep("B", 60 * 24), + Datetime = lubridate::as_datetime(0) + lubridate::minutes(0:(60*24-1)), + MEDI = c(rep(sample(seq(0,1,0.1), 60*8, replace = TRUE)), + rep(sample(1:1000, 16, replace = TRUE), each = 60)) + ) + +dataset1 \%>\% + dplyr::reframe(barroso_lighting_metrics(MEDI, Datetime, as.df = TRUE)) + +} +\references{ +Barroso, A., Simons, K., & Jager, P. de. (2014). Metrics of circadian +lighting for clinical investigations. \emph{Lighting Research & Technology}, +46(6), 637–649. \url{https://doi.org/10.1177/1477153513502664} + +Hartmeyer, S.L., Andersen, M. (2023). Towards a framework for light-dosimetry studies: +Quantification metrics. \emph{Lighting Research & Technology}. +\url{https://doi.org/10.1177/14771535231170500} +} diff --git a/man/bright_dark_period.Rd b/man/bright_dark_period.Rd index 740039e..4c84673 100755 --- a/man/bright_dark_period.Rd +++ b/man/bright_dark_period.Rd @@ -18,19 +18,20 @@ bright_dark_period( \arguments{ \item{Light.vector}{Numeric vector containing the light data.} -\item{Time.vector}{Vector containing the time data. Can be HMS or POSIXct.} +\item{Time.vector}{Vector containing the time data. Can be \link[base]{POSIXct}, +\link[hms]{hms}, \link[lubridate]{duration}, or \link[base]{difftime}.} \item{period}{String indicating the type of period to look for. Can be either \code{"brightest"}(the default) or \code{"darkest"}.} \item{timespan}{The timespan across which to calculate. Can be either a -\code{lubridate::duration()} or a \code{lubridate::duration()} string, e.g., +\link[lubridate]{duration} or a \link[lubridate]{duration} string, e.g., \code{"1 day"} or \code{"10 sec"}.} \item{epoch}{The epoch at which the data was sampled. Can be either a -\code{lubridate::duration()} or a string. If it is a string, it needs to be +\link[lubridate]{duration} or a string. If it is a string, it needs to be either \code{"dominant.epoch"} (the default) for a guess based on the data, or a valid -\code{lubridate::duration()} string, e.g., \code{"1 day"} or \code{"10 sec"}.} +\link[lubridate]{duration} string, e.g., \code{"1 day"} or \code{"10 sec"}.} \item{loop}{Logical. Should the data be looped? If \code{TRUE}, a full copy of the data will be concatenated at the end of the data. Makes only sense for 24 h data. @@ -47,7 +48,7 @@ A named list with the \code{mean}, \code{onset}, \code{midpoint}, and \code{offs calculated brightest or darkest period, or if \code{as.df == TRUE} a data frame with columns named \verb{\{period\}_\{timespan\}_\{metric\}}. The output type corresponds to the type of \code{Time.vector}, e.g., if \code{Time.vector} is HMS, the timing metrics -will be also HMS, and vice versa for POSIXct and numeric. +will be also HMS, and vice versa for POSIXct. } \description{ This function finds the brightest or darkest continuous period of a given @@ -61,6 +62,7 @@ Assumes regular 24h light data. Otherwise, results may not be meaningful. Looping the data is recommended for finding the darkest period. } \examples{ +# Dataset with light > 250lx between 06:00 and 18:00 dataset1 <- tibble::tibble( Id = rep("A", 24), @@ -74,6 +76,22 @@ dataset1 \%>\% dataset1 \%>\% dplyr::reframe(bright_dark_period(MEDI, Datetime, "darkest", "5 hours", loop = TRUE, as.df = TRUE)) + +# Dataset with duration as Time.vector +dataset2 <- + tibble::tibble( + Id = rep("A", 24), + Datetime = lubridate::dhours(0:23), + MEDI = c(rep(1, 6), rep(250, 13), rep(1, 5)) + ) + +dataset2 \%>\% + dplyr::reframe(bright_dark_period(MEDI, Datetime, "brightest", "10 hours", + as.df = TRUE)) +dataset2 \%>\% + dplyr::reframe(bright_dark_period(MEDI, Datetime, "darkest", "5 hours", + loop = TRUE, as.df = TRUE)) + } \references{ Hartmeyer, S.L., Andersen, M. (2023). Towards a framework for light-dosimetry studies: @@ -85,12 +103,17 @@ Other metrics: \code{\link{centroidLE}()}, \code{\link{disparity_index}()}, \code{\link{duration_above_threshold}()}, +\code{\link{exponential_moving_average}()}, +\code{\link{frequency_crossing_threshold}()}, \code{\link{interdaily_stability}()}, \code{\link{intradaily_variability}()}, \code{\link{midpointCE}()}, \code{\link{nvRC}()}, \code{\link{nvRD}()}, \code{\link{nvRD_cumulative_response}()}, +\code{\link{period_above_threshold}()}, +\code{\link{pulses_above_threshold}()}, +\code{\link{threshold_for_duration}()}, \code{\link{timing_above_threshold}()} } \concept{metrics} diff --git a/man/centroidLE.Rd b/man/centroidLE.Rd index 8444589..d4e15e2 100644 --- a/man/centroidLE.Rd +++ b/man/centroidLE.Rd @@ -15,19 +15,19 @@ centroidLE( \arguments{ \item{Light.vector}{Numeric vector containing the light data.} -\item{Time.vector}{Vector containing the time data. Can be numeric, HMS or POSIXct.} +\item{Time.vector}{Vector containing the time data. Can be \link[base]{POSIXct}, \link[hms]{hms}, +\link[lubridate]{duration}, or \link[base]{difftime}.} \item{bin.size}{Value specifying size of bins to average the light data over. -If \code{Time.vector} is of type POSIXct or HMS, \code{bin.size} must be either a -\code{lubridate::duration()} or a \code{lubridate::duration()} string, e.g., -\code{"1 day"} or \code{"10 sec"}. Otherwise, if \code{Time.vector} is numeric, \code{bin.size} -must be also numeric. If nothing is provided, no binning will be performed.} +Must be either a \link[lubridate]{duration} or a \link[lubridate]{duration} string, e.g., +\code{"1 day"} or \code{"10 sec"}. If nothing is provided, no binning will be performed.} \item{na.rm}{Logical. Should missing values be removed for the calculation? Defaults to \code{FALSE}.} -\item{as.df}{Logical. Should the output be returned as a data frame? Defaults -to \code{FALSE}.} +\item{as.df}{Logical. Should the output be returned as a data frame? If \code{TRUE}, a data +frame with a single column named \code{centroidLE} will be returned. +Defaults to \code{FALSE}.} } \value{ Single column data frame or vector. @@ -49,7 +49,7 @@ dataset1 \%>\% "Centroid of light exposure" = centroidLE(MEDI, Datetime, "2 hours") ) -# Dataset with HMS time vector +# Dataset with hms time vector dataset2 <- tibble::tibble( Id = rep("A", 24), @@ -61,16 +61,16 @@ dataset2 \%>\% "Centroid of light exposure" = centroidLE(MEDI, Time, "2 hours") ) -# Dataset with numeric time vector +# Dataset with duration time vector dataset3 <- tibble::tibble( Id = rep("A", 24), - Hour = 0:23, + Hour = lubridate::duration(0:23, "hours"), MEDI = c(rep(1, 6), rep(250, 13), rep(1, 5)) ) dataset3 \%>\% dplyr::reframe( - "Centroid of light exposure" = centroidLE(MEDI, Hour, 2) + "Centroid of light exposure" = centroidLE(MEDI, Hour, "2 hours") ) } @@ -90,12 +90,17 @@ Other metrics: \code{\link{bright_dark_period}()}, \code{\link{disparity_index}()}, \code{\link{duration_above_threshold}()}, +\code{\link{exponential_moving_average}()}, +\code{\link{frequency_crossing_threshold}()}, \code{\link{interdaily_stability}()}, \code{\link{intradaily_variability}()}, \code{\link{midpointCE}()}, \code{\link{nvRC}()}, \code{\link{nvRD}()}, \code{\link{nvRD_cumulative_response}()}, +\code{\link{period_above_threshold}()}, +\code{\link{pulses_above_threshold}()}, +\code{\link{threshold_for_duration}()}, \code{\link{timing_above_threshold}()} } \concept{metrics} diff --git a/man/disparity_index.Rd b/man/disparity_index.Rd index 2da1011..70297e5 100644 --- a/man/disparity_index.Rd +++ b/man/disparity_index.Rd @@ -11,8 +11,9 @@ disparity_index(Light.vector, na.rm = FALSE, as.df = FALSE) \item{na.rm}{Logical. Should missing values be removed? Defaults to FALSE} -\item{as.df}{Logical. Should the output be returned as a data frame? Defaults -to FALSE} +\item{as.df}{Logical. Should the output be returned as a data frame? If \code{TRUE}, a data +frame with a single column named \code{disparity_index} will be returned. +Defaults to \code{FALSE}.} } \value{ Single column data frame or vector. @@ -50,12 +51,17 @@ Other metrics: \code{\link{bright_dark_period}()}, \code{\link{centroidLE}()}, \code{\link{duration_above_threshold}()}, +\code{\link{exponential_moving_average}()}, +\code{\link{frequency_crossing_threshold}()}, \code{\link{interdaily_stability}()}, \code{\link{intradaily_variability}()}, \code{\link{midpointCE}()}, \code{\link{nvRC}()}, \code{\link{nvRD}()}, \code{\link{nvRD_cumulative_response}()}, +\code{\link{period_above_threshold}()}, +\code{\link{pulses_above_threshold}()}, +\code{\link{threshold_for_duration}()}, \code{\link{timing_above_threshold}()} } \concept{metrics} diff --git a/man/duration_above_threshold.Rd b/man/duration_above_threshold.Rd index 06f9592..66a9315 100644 --- a/man/duration_above_threshold.Rd +++ b/man/duration_above_threshold.Rd @@ -2,7 +2,7 @@ % Please edit documentation in R/metric_duration_above_threshold.R \name{duration_above_threshold} \alias{duration_above_threshold} -\title{Time above/below threshold or within threshold range} +\title{Duration above/below threshold or within threshold range} \usage{ duration_above_threshold( Light.vector, @@ -17,7 +17,8 @@ duration_above_threshold( \arguments{ \item{Light.vector}{Numeric vector containing the light data.} -\item{Time.vector}{Vector containing the time data. Can be numeric, HMS or POSIXct.} +\item{Time.vector}{Vector containing the time data. Can be \link[base]{POSIXct}, \link[hms]{hms}, +\link[lubridate]{duration}, or \link[base]{difftime}.} \item{comparison}{String specifying whether the time above or below threshold should be calculated. Can be either \code{"above"} (the default) or \code{"below"}. If @@ -28,27 +29,26 @@ threshold light level(s) to compare with. If a vector with two values is provide the time within the two thresholds will be calculated.} \item{epoch}{The epoch at which the data was sampled. Can be either a -\code{lubridate::duration()} or a string. If it is a string, it needs to be +\link[lubridate]{duration} or a string. If it is a string, it needs to be either \code{"dominant.epoch"} (the default) for a guess based on the data, or a valid -\code{lubridate::duration()} string, e.g., \code{"1 day"} or \code{"10 sec"}.} +\link[lubridate]{duration} string, e.g., \code{"1 day"} or \code{"10 sec"}.} \item{na.rm}{Logical. Should missing values (NA) be removed for the calculation? Defaults to \code{FALSE}.} \item{as.df}{Logical. Should a data frame with be returned? If \code{TRUE}, a data -frame with a single column named \verb{TAT_\{threshold\}} will be returned. +frame with a single column named \verb{duration_\{comparison\}_\{threshold\}} will be returned. Defaults to \code{FALSE}.} } \value{ -A duration object (see \code{\link[lubridate]{duration}}) as single value, -or single column data frame. +A \link[lubridate]{duration} object as single value, or single column data frame. } \description{ -This function calculates the time spent above/below a specified threshold +This function calculates the duration spent above/below a specified threshold light level or within a specified range of light levels. } \examples{ -N <- 50 +N <- 60 # Dataset with epoch = 1min dataset1 <- tibble::tibble( @@ -87,12 +87,17 @@ Other metrics: \code{\link{bright_dark_period}()}, \code{\link{centroidLE}()}, \code{\link{disparity_index}()}, +\code{\link{exponential_moving_average}()}, +\code{\link{frequency_crossing_threshold}()}, \code{\link{interdaily_stability}()}, \code{\link{intradaily_variability}()}, \code{\link{midpointCE}()}, \code{\link{nvRC}()}, \code{\link{nvRD}()}, \code{\link{nvRD_cumulative_response}()}, +\code{\link{period_above_threshold}()}, +\code{\link{pulses_above_threshold}()}, +\code{\link{threshold_for_duration}()}, \code{\link{timing_above_threshold}()} } \concept{metrics} diff --git a/man/exponential_moving_average.Rd b/man/exponential_moving_average.Rd new file mode 100644 index 0000000..f15e305 --- /dev/null +++ b/man/exponential_moving_average.Rd @@ -0,0 +1,84 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/metric_EMA.R +\name{exponential_moving_average} +\alias{exponential_moving_average} +\title{Exponential moving average filter (EMA)} +\usage{ +exponential_moving_average( + Light.vector, + Time.vector, + decay = "90 min", + epoch = "dominant.epoch" +) +} +\arguments{ +\item{Light.vector}{Numeric vector containing the light data. Missing values are +replaced by 0.} + +\item{Time.vector}{Vector containing the time data. Can be \link[base]{POSIXct}, \link[hms]{hms}, +\link[lubridate]{duration}, or \link[base]{difftime}.} + +\item{decay}{The decay half-life controlling the exponential smoothing. +Can be either a \link[lubridate]{duration} or a string. If it is a string, it +needs to be a valid \link[lubridate]{duration} string, e.g., \code{"1 day"} or \code{"10 sec"}. +The default is set to \code{"90 mins"} for a biologically relevant estimate (see +the reference paper).} + +\item{epoch}{The epoch at which the data was sampled. Can be either a +\link[lubridate]{duration} or a string. If it is a string, it needs to be +either \code{"dominant.epoch"} (the default) for a guess based on the data, or a valid +\link[lubridate]{duration} string, e.g., \code{"1 day"} or \code{"10 sec"}.} +} +\value{ +A numeric vector containing the smoothed light data. The output has the same +length as \code{Light.vector}. +} +\description{ +This function smoothes the data using an exponential moving average filter +with a specified decay half-life. +} +\details{ +The timeseries is assumed to be regular. Missing values in the +light data will be replaced by 0. +} +\examples{ +sample.data.environment.EMA = sample.data.environment \%>\% + dplyr::filter(Id == "Participant") \%>\% + filter_Datetime(length = lubridate::days(2)) \%>\% + dplyr::mutate(MEDI.EMA = exponential_moving_average(MEDI, Datetime)) + +# Plot to compare results +sample.data.environment.EMA \%>\% + ggplot2::ggplot(ggplot2::aes(x = Datetime)) + + ggplot2::geom_line(ggplot2::aes(y = MEDI), colour = "black") + + ggplot2::geom_line(ggplot2::aes(y = MEDI.EMA), colour = "red") + +} +\references{ +Price, L. L. A. (2014). On the Role of Exponential Smoothing in Circadian +Dosimetry. \emph{Photochemistry and Photobiology}, 90(5), 1184-1192. +\url{https://doi.org/10.1111/php.12282} + +Hartmeyer, S.L., Andersen, M. (2023). Towards a framework for light-dosimetry studies: +Quantification metrics. \emph{Lighting Research & Technology}. +\url{https://doi.org/10.1177/14771535231170500} +} +\seealso{ +Other metrics: +\code{\link{bright_dark_period}()}, +\code{\link{centroidLE}()}, +\code{\link{disparity_index}()}, +\code{\link{duration_above_threshold}()}, +\code{\link{frequency_crossing_threshold}()}, +\code{\link{interdaily_stability}()}, +\code{\link{intradaily_variability}()}, +\code{\link{midpointCE}()}, +\code{\link{nvRC}()}, +\code{\link{nvRD}()}, +\code{\link{nvRD_cumulative_response}()}, +\code{\link{period_above_threshold}()}, +\code{\link{pulses_above_threshold}()}, +\code{\link{threshold_for_duration}()}, +\code{\link{timing_above_threshold}()} +} +\concept{metrics} diff --git a/man/frequency_crossing_threshold.Rd b/man/frequency_crossing_threshold.Rd new file mode 100644 index 0000000..44c927f --- /dev/null +++ b/man/frequency_crossing_threshold.Rd @@ -0,0 +1,75 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/metric_frequency_crossing_threshold.R +\name{frequency_crossing_threshold} +\alias{frequency_crossing_threshold} +\title{Frequency of crossing light threshold} +\usage{ +frequency_crossing_threshold( + Light.vector, + threshold, + na.rm = FALSE, + as.df = FALSE +) +} +\arguments{ +\item{Light.vector}{Numeric vector containing the light data.} + +\item{threshold}{Single numeric value specifying the threshold light level to compare with.} + +\item{na.rm}{Logical. Should missing light values be removed? Defaults to \code{FALSE}.} + +\item{as.df}{Logical. Should the output be returned as a data frame? If \code{TRUE}, a data +frame with a single column named \verb{frequency_crossing_\{threshold\}} will be returned. +Defaults to \code{FALSE}.} +} +\value{ +Data frame or matrix with pairs of threshold and calculated values. +} +\description{ +This functions calculates the number of times a given threshold +light level is crossed. +} +\examples{ +N = 60 +dataset1 <- + tibble::tibble( + Id = rep("A", N), + Datetime = lubridate::as_datetime(0) + lubridate::minutes(1:N), + MEDI = sample(c(sample(1:249, N / 2), sample(250:1000, N / 2))), + ) + +dataset1 \%>\% + dplyr::reframe("Frequency crossing 250lx" = frequency_crossing_threshold(MEDI, threshold = 250)) + +dataset1 \%>\% + dplyr::reframe(frequency_crossing_threshold(MEDI, threshold = 250, as.df = TRUE)) + +} +\references{ +Alvarez, A. A., & Wildsoet, C. F. (2013). Quantifying light +exposure patterns in young adult students. \emph{Journal of Modern Optics}, +60(14), 1200–1208. \url{https://doi.org/10.1080/09500340.2013.845700} + +Hartmeyer, S.L., Andersen, M. (2023). Towards a framework for light-dosimetry studies: +Quantification metrics. \emph{Lighting Research & Technology}. +\url{https://doi.org/10.1177/14771535231170500} +} +\seealso{ +Other metrics: +\code{\link{bright_dark_period}()}, +\code{\link{centroidLE}()}, +\code{\link{disparity_index}()}, +\code{\link{duration_above_threshold}()}, +\code{\link{exponential_moving_average}()}, +\code{\link{interdaily_stability}()}, +\code{\link{intradaily_variability}()}, +\code{\link{midpointCE}()}, +\code{\link{nvRC}()}, +\code{\link{nvRD}()}, +\code{\link{nvRD_cumulative_response}()}, +\code{\link{period_above_threshold}()}, +\code{\link{pulses_above_threshold}()}, +\code{\link{threshold_for_duration}()}, +\code{\link{timing_above_threshold}()} +} +\concept{metrics} diff --git a/man/import_Dataset.Rd b/man/import_Dataset.Rd index 204e75d..698b499 100644 --- a/man/import_Dataset.Rd +++ b/man/import_Dataset.Rd @@ -142,7 +142,7 @@ dataset <- import_Dataset("LYS", filepath, auto.plot = FALSE) #> #> Successfully read in 11'422 observations across 1 Ids from 1 LYS-file(s). #> Timezone set is UTC. -#> The system timezone is Europe/Berlin. Please correct if necessary! +#> The system timezone is Europe/Zurich. Please correct if necessary! #> #> First Observation: 2023-06-21 00:00:12 #> Last Observation: 2023-06-22 23:59:48 @@ -163,7 +163,7 @@ dataset <- import$ActLumus(filepath, auto.plot = FALSE) #> #> Successfully read in 61'016 observations across 1 Ids from 1 ActLumus-file(s). #> Timezone set is UTC. -#> The system timezone is Europe/Berlin. Please correct if necessary! +#> The system timezone is Europe/Zurich. Please correct if necessary! #> #> First Observation: 2023-08-28 08:47:54 #> Last Observation: 2023-09-04 10:17:04 diff --git a/man/interdaily_stability.Rd b/man/interdaily_stability.Rd index 504376f..8c50fd3 100755 --- a/man/interdaily_stability.Rd +++ b/man/interdaily_stability.Rd @@ -18,8 +18,9 @@ interdaily_stability( \item{na.rm}{Logical. Should missing values be removed? Defaults to \code{FALSE}.} -\item{as.df}{Logical. Should the output be returned as a data frame? Defaults -to \code{FALSE}.} +\item{as.df}{Logical. Should the output be returned as a data frame? If \code{TRUE}, a data +frame with a single column named \code{interdaily_stability} will be returned. +Defaults to \code{FALSE}.} } \value{ Numeric value or dataframe with column 'IS'. @@ -67,11 +68,16 @@ Other metrics: \code{\link{centroidLE}()}, \code{\link{disparity_index}()}, \code{\link{duration_above_threshold}()}, +\code{\link{exponential_moving_average}()}, +\code{\link{frequency_crossing_threshold}()}, \code{\link{intradaily_variability}()}, \code{\link{midpointCE}()}, \code{\link{nvRC}()}, \code{\link{nvRD}()}, \code{\link{nvRD_cumulative_response}()}, +\code{\link{period_above_threshold}()}, +\code{\link{pulses_above_threshold}()}, +\code{\link{threshold_for_duration}()}, \code{\link{timing_above_threshold}()} } \concept{metrics} diff --git a/man/intradaily_variability.Rd b/man/intradaily_variability.Rd index be25fdc..87dbb43 100755 --- a/man/intradaily_variability.Rd +++ b/man/intradaily_variability.Rd @@ -18,8 +18,9 @@ intradaily_variability( \item{na.rm}{Logical. Should missing values be removed? Defaults to \code{FALSE}.} -\item{as.df}{Logical. Should the output be returned as a data frame? Defaults -to \code{FALSE}.} +\item{as.df}{Logical. Should the output be returned as a data frame? If \code{TRUE}, a data +frame with a single column named \code{intradaily_variability} will be returned. +Defaults to \code{FALSE}.} } \value{ Numeric value or dataframe with column 'IV'. @@ -64,11 +65,16 @@ Other metrics: \code{\link{centroidLE}()}, \code{\link{disparity_index}()}, \code{\link{duration_above_threshold}()}, +\code{\link{exponential_moving_average}()}, +\code{\link{frequency_crossing_threshold}()}, \code{\link{interdaily_stability}()}, \code{\link{midpointCE}()}, \code{\link{nvRC}()}, \code{\link{nvRD}()}, \code{\link{nvRD_cumulative_response}()}, +\code{\link{period_above_threshold}()}, +\code{\link{pulses_above_threshold}()}, +\code{\link{threshold_for_duration}()}, \code{\link{timing_above_threshold}()} } \concept{metrics} diff --git a/man/midpointCE.Rd b/man/midpointCE.Rd index 7b73154..7344af1 100644 --- a/man/midpointCE.Rd +++ b/man/midpointCE.Rd @@ -7,16 +7,17 @@ midpointCE(Light.vector, Time.vector, na.rm = FALSE, as.df = FALSE) } \arguments{ -\item{Light.vector}{Numeric vector containing the light data. Missing values are -replaced with 0.} +\item{Light.vector}{Numeric vector containing the light data.} -\item{Time.vector}{Vector containing the time data. Can be numeric, HMS or POSIXct.} +\item{Time.vector}{Vector containing the time data. Can be \link[base]{POSIXct}, \link[hms]{hms}, +\link[lubridate]{duration}, or \link[base]{difftime}.} -\item{na.rm}{Logical. Should missing values be removed for the calculation? -Defaults to \code{FALSE}.} +\item{na.rm}{Logical. Should missing values be removed for the calculation? If \code{TRUE}, +missing values will be replaced by zero. Defaults to \code{FALSE}.} -\item{as.df}{Logical. Should the output be returned as a data frame? Defaults -to \code{FALSE}.} +\item{as.df}{Logical. Should the output be returned as a data frame? If \code{TRUE}, a data +frame with a single column named \code{midpointCE} will be returned. +Defaults to \code{FALSE}.} } \value{ Single column data frame or vector. @@ -26,7 +27,6 @@ This function calculates the timing corresponding to half of the cumulative light exposure within the given time series. } \examples{ -# Dataset with POSIXct time vector dataset1 <- tibble::tibble( Id = rep("A", 24), @@ -50,11 +50,11 @@ dataset2 \%>\% "Midpoint of cmulative exposure" = midpointCE(MEDI, Time) ) -# Dataset with numeric time vector +# Dataset with duration time vector dataset3 <- tibble::tibble( Id = rep("A", 24), - Hour = 0:23, + Hour = lubridate::duration(0:23, "hours"), MEDI = c(rep(1, 6), rep(250, 13), rep(1, 5)) ) dataset3 \%>\% @@ -79,11 +79,16 @@ Other metrics: \code{\link{centroidLE}()}, \code{\link{disparity_index}()}, \code{\link{duration_above_threshold}()}, +\code{\link{exponential_moving_average}()}, +\code{\link{frequency_crossing_threshold}()}, \code{\link{interdaily_stability}()}, \code{\link{intradaily_variability}()}, \code{\link{nvRC}()}, \code{\link{nvRD}()}, \code{\link{nvRD_cumulative_response}()}, +\code{\link{period_above_threshold}()}, +\code{\link{pulses_above_threshold}()}, +\code{\link{threshold_for_duration}()}, \code{\link{timing_above_threshold}()} } \concept{metrics} diff --git a/man/nvRC.Rd b/man/nvRC.Rd index 8bc09e9..d1db526 100644 --- a/man/nvRC.Rd +++ b/man/nvRC.Rd @@ -17,12 +17,13 @@ nvRC( \item{Illuminance.vector}{Numeric vector containing the Illuminance data.} -\item{Time.vector}{Vector containing the time data. Can be numeric, HMS or POSIXct.} +\item{Time.vector}{Vector containing the time data. Can be \link[base]{POSIXct}, +\link[hms]{hms}, \link[lubridate]{duration}, or \link[base]{difftime}.} \item{epoch}{The epoch at which the data was sampled. Can be either a -\code{lubridate::duration()} or a string. If it is a string, it needs to be +\link[lubridate]{duration} or a string. If it is a string, it needs to be either \code{"dominant.epoch"} (the default) for a guess based on the data, or a valid -\code{lubridate::duration()} string, e.g., \code{"1 day"} or \code{"10 sec"}.} +\link[lubridate]{duration} string, e.g., \code{"1 day"} or \code{"10 sec"}.} \item{sleep.onset}{The time of habitual sleep onset. Can be HMS, numeric, or NULL. If NULL (the default), then the data is assumed to start at habitual sleep onset. @@ -112,6 +113,7 @@ exposure, to determine the final output \eqn{r_{C}(t)}. responses. } \examples{ + dataset1 <- tibble::tibble( Id = rep("B", 60 * 48), @@ -120,12 +122,20 @@ dataset1 <- rep(0, 60*8), rep(sample(1:1000, 16, replace = TRUE), each = 60)), MEDI = Illuminance * rep(sample(0.5:1.5, 48, replace = TRUE), each = 60) ) - +# Time.vector as POSIXct dataset1.nvRC <- dataset1 \%>\% dplyr::mutate( nvRC = nvRC(MEDI, Illuminance, Datetime, sleep.onset = hms::as_hms("22:00:00")) ) +# Time.vector as difftime +dataset2 <- dataset1 \%>\% + dplyr::mutate(Datetime = Datetime - lubridate::as_datetime(lubridate::dhours(22))) +dataset2.nvRC <- dataset2 \%>\% + dplyr::mutate( + nvRC = nvRC(MEDI, Illuminance, Datetime, sleep.onset = lubridate::dhours(0)) + ) + } \references{ Amundadottir, M.L. (2016). Light-driven model for identifying @@ -139,11 +149,16 @@ Other metrics: \code{\link{centroidLE}()}, \code{\link{disparity_index}()}, \code{\link{duration_above_threshold}()}, +\code{\link{exponential_moving_average}()}, +\code{\link{frequency_crossing_threshold}()}, \code{\link{interdaily_stability}()}, \code{\link{intradaily_variability}()}, \code{\link{midpointCE}()}, \code{\link{nvRD}()}, \code{\link{nvRD_cumulative_response}()}, +\code{\link{period_above_threshold}()}, +\code{\link{pulses_above_threshold}()}, +\code{\link{threshold_for_duration}()}, \code{\link{timing_above_threshold}()} } \concept{metrics} diff --git a/man/nvRD.Rd b/man/nvRD.Rd index d7825d1..af27031 100644 --- a/man/nvRD.Rd +++ b/man/nvRD.Rd @@ -11,12 +11,13 @@ nvRD(MEDI.vector, Illuminance.vector, Time.vector, epoch = "dominant.epoch") \item{Illuminance.vector}{Numeric vector containing the Illuminance data.} -\item{Time.vector}{Vector containing the time data. Can be numeric, HMS or POSIXct.} +\item{Time.vector}{Vector containing the time data. Can be \link[base]{POSIXct}, \link[hms]{hms}, +\link[lubridate]{duration}, or \link[base]{difftime}.} \item{epoch}{The epoch at which the data was sampled. Can be either a -\code{lubridate::duration()} or a string. If it is a string, it needs to be +\link[lubridate]{duration} or a string. If it is a string, it needs to be either \code{"dominant.epoch"} (the default) for a guess based on the data, or a valid -\code{lubridate::duration()} string, e.g., \code{"1 day"} or \code{"10 sec"}.} +\link[lubridate]{duration} string, e.g., \code{"1 day"} or \code{"10 sec"}.} } \value{ A numeric vector containing the nvRD data. The output has the same @@ -125,11 +126,16 @@ Other metrics: \code{\link{centroidLE}()}, \code{\link{disparity_index}()}, \code{\link{duration_above_threshold}()}, +\code{\link{exponential_moving_average}()}, +\code{\link{frequency_crossing_threshold}()}, \code{\link{interdaily_stability}()}, \code{\link{intradaily_variability}()}, \code{\link{midpointCE}()}, \code{\link{nvRC}()}, \code{\link{nvRD_cumulative_response}()}, +\code{\link{period_above_threshold}()}, +\code{\link{pulses_above_threshold}()}, +\code{\link{threshold_for_duration}()}, \code{\link{timing_above_threshold}()} } \concept{metrics} diff --git a/man/nvRD_cumulative_response.Rd b/man/nvRD_cumulative_response.Rd index f1a082a..37a9451 100644 --- a/man/nvRD_cumulative_response.Rd +++ b/man/nvRD_cumulative_response.Rd @@ -15,15 +15,16 @@ nvRD_cumulative_response( \item{nvRD}{Numeric vector containing the non-visual direct response. See \code{\link{nvRD}}.} -\item{Time.vector}{Vector containing the time data. Can be numeric, HMS or POSIXct.} +\item{Time.vector}{Vector containing the time data. Can be \link[base]{POSIXct}, \link[hms]{hms}, +\link[lubridate]{duration}, or \link[base]{difftime}.} \item{epoch}{The epoch at which the data was sampled. Can be either a -\code{lubridate::duration()} or a string. If it is a string, it needs to be +\link[lubridate]{duration} or a string. If it is a string, it needs to be either \code{"dominant.epoch"} (the default) for a guess based on the data, or a valid -\code{lubridate::duration()} string, e.g., \code{"1 day"} or \code{"10 sec"}.} +\link[lubridate]{duration} string, e.g., \code{"1 day"} or \code{"10 sec"}.} \item{as.df}{Logical. Should a data frame with be returned? If \code{TRUE}, a data -frame with a single column named \verb{TAT_\{threshold\}} will be returned. +frame with a single column named \code{nvRD_cumulative} will be returned. Defaults to \code{FALSE}.} } \value{ @@ -63,11 +64,16 @@ Other metrics: \code{\link{centroidLE}()}, \code{\link{disparity_index}()}, \code{\link{duration_above_threshold}()}, +\code{\link{exponential_moving_average}()}, +\code{\link{frequency_crossing_threshold}()}, \code{\link{interdaily_stability}()}, \code{\link{intradaily_variability}()}, \code{\link{midpointCE}()}, \code{\link{nvRC}()}, \code{\link{nvRD}()}, +\code{\link{period_above_threshold}()}, +\code{\link{pulses_above_threshold}()}, +\code{\link{threshold_for_duration}()}, \code{\link{timing_above_threshold}()} } \concept{metrics} diff --git a/man/period_above_threshold.Rd b/man/period_above_threshold.Rd new file mode 100644 index 0000000..617150d --- /dev/null +++ b/man/period_above_threshold.Rd @@ -0,0 +1,112 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/metric_period_above_threshold.R +\name{period_above_threshold} +\alias{period_above_threshold} +\title{Length of longest continuous period above/below threshold} +\usage{ +period_above_threshold( + Light.vector, + Time.vector, + comparison = c("above", "below"), + threshold, + epoch = "dominant.epoch", + loop = FALSE, + na.replace = FALSE, + na.rm = FALSE, + as.df = FALSE +) +} +\arguments{ +\item{Light.vector}{Numeric vector containing the light data.} + +\item{Time.vector}{Vector containing the time data. Can be \link[base]{POSIXct}, +\link[hms]{hms}, \link[lubridate]{duration}, or \link[base]{difftime}.} + +\item{comparison}{String specifying whether the period of light levels above or +below threshold should be calculated. Can be either \code{"above"} (the default) +or \code{"below"}. If two values are provided for \code{threshold}, this argument will be ignored.} + +\item{threshold}{Single numeric value or two numeric values specifying the +threshold light level(s) to compare with. If a vector with two values is provided, +the period of light levels within the two thresholds will be calculated.} + +\item{epoch}{The epoch at which the data was sampled. Can be either a +\link[lubridate]{duration} or a string. If it is a string, it needs to be +either \code{"dominant.epoch"} (the default) for a guess based on the data, or a valid +\link[lubridate]{duration} string, e.g., \code{"1 day"} or \code{"10 sec"}.} + +\item{loop}{Logical. Should the data be looped? Defaults to \code{FALSE}.} + +\item{na.replace}{Logical. Should missing values (NA) be replaced +for the calculation? If \code{TRUE} missing values will not be removed but will +result in \code{FALSE} when comparing \code{Light.vector} with \code{threshold}. +Defaults to \code{FALSE}.} + +\item{na.rm}{Logical. Should missing values (NA) be removed for the calculation? +If \code{TRUE}, this argument will override \code{na.replace}. Defaults to \code{FALSE}.} + +\item{as.df}{Logical. Should a data frame be returned? If \code{TRUE}, a data +frame with a single column named \verb{period_\{comparison\}_\{threshold\}} will be returned. +Defaults to \code{FALSE}.} +} +\value{ +A duration object (see \code{\link[lubridate]{duration}}) as single value, +or single column data frame. +} +\description{ +This function finds the length of the longest continous period above/below +a specified threshold light level or within a specified range of light levels. +} +\examples{ + +N <- 60 +# Dataset with continous period of >250lx for 35min +dataset1 <- + tibble::tibble( + Id = rep("A", N), + Datetime = lubridate::as_datetime(0) + lubridate::minutes(1:N), + MEDI = c(sample(1:249, N-35, replace = TRUE), + sample(250:1000, 35, replace = TRUE)) + ) + +dataset1 \%>\% + dplyr::reframe("Period >250lx" = period_above_threshold(MEDI, Datetime, threshold = 250)) + +dataset1 \%>\% + dplyr::reframe("Period <250lx" = period_above_threshold(MEDI, Datetime, "below", threshold = 250)) + +# Dataset with continous period of 100-250lx for 20min +dataset2 <- + tibble::tibble( + Id = rep("B", N), + Datetime = lubridate::as_datetime(0) + lubridate::minutes(1:N), + MEDI = c(sample(c(1:99, 251-1000), N-20, replace = TRUE), + sample(100:250, 20, replace = TRUE)), + ) +dataset2 \%>\% + dplyr::reframe("Period 250lx" = period_above_threshold(MEDI, Datetime, threshold = c(100,250))) + +# Return data frame +dataset1 \%>\% + dplyr::reframe(period_above_threshold(MEDI, Datetime, threshold = 250, as.df = TRUE)) + +} +\seealso{ +Other metrics: +\code{\link{bright_dark_period}()}, +\code{\link{centroidLE}()}, +\code{\link{disparity_index}()}, +\code{\link{duration_above_threshold}()}, +\code{\link{exponential_moving_average}()}, +\code{\link{frequency_crossing_threshold}()}, +\code{\link{interdaily_stability}()}, +\code{\link{intradaily_variability}()}, +\code{\link{midpointCE}()}, +\code{\link{nvRC}()}, +\code{\link{nvRD}()}, +\code{\link{nvRD_cumulative_response}()}, +\code{\link{pulses_above_threshold}()}, +\code{\link{threshold_for_duration}()}, +\code{\link{timing_above_threshold}()} +} +\concept{metrics} diff --git a/man/pulses_above_threshold.Rd b/man/pulses_above_threshold.Rd new file mode 100644 index 0000000..2991193 --- /dev/null +++ b/man/pulses_above_threshold.Rd @@ -0,0 +1,132 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/metric_pulses_above_threshold.R +\name{pulses_above_threshold} +\alias{pulses_above_threshold} +\title{Pulses above threshold} +\usage{ +pulses_above_threshold( + Light.vector, + Time.vector, + comparison = c("above", "below"), + threshold, + min.length = "8 mins", + max.interrupt = "2 mins", + prop.interrupt = 0.25, + epoch = "dominant.epoch", + return.indices = FALSE, + na.rm = FALSE, + as.df = FALSE +) +} +\arguments{ +\item{Light.vector}{Numeric vector containing the light data. Missing values will +be considered as \code{FALSE} when comparing light levels against the threshold.} + +\item{Time.vector}{Vector containing the time data. Can be \link[base]{POSIXct}, +\link[hms]{hms}, \link[lubridate]{duration}, or \link[base]{difftime}.} + +\item{comparison}{String specifying whether the time above or below threshold +should be calculated. Can be either \code{"above"} (the default) or \code{"below"}. If +two values are provided for \code{threshold}, this argument will be ignored.} + +\item{threshold}{Single numeric value or two numeric values specifying the +threshold light level(s) to compare with. If a vector with two values is provided, +the timing corresponding to light levels between the two thresholds will be +calculated.} + +\item{min.length}{The minimum length of a pulse. Can be either a +\link[lubridate]{duration} or a string. If it is a string, it needs to be a valid +\link[lubridate]{duration} string, e.g., \code{"1 day"} or \code{"10 sec"}. Defaults to +\code{"8 mins"} as in Wilson et al. (2018).} + +\item{max.interrupt}{Maximum length of each episode of interruptions. Can be either a +\link[lubridate]{duration} or a string. If it is a string, it needs to be a valid +\link[lubridate]{duration} string, e.g., \code{"1 day"} or \code{"10 sec"}. Defaults to +\code{"2 mins"} as in Wilson et al. (2018).} + +\item{prop.interrupt}{Numeric value between \code{0} and \code{1} specifying the +maximum proportion of the total number of interruptions. Defaults to \code{0.25} +as in Wilson et al. (2018).} + +\item{epoch}{The epoch at which the data was sampled. Can be either a +\link[lubridate]{duration} or a string. If it is a string, it needs to be +either \code{"dominant.epoch"} (the default) for a guess based on the data, or a valid +\link[lubridate]{duration} string, e.g., \code{"1 day"} or \code{"10 sec"}.} + +\item{return.indices}{Logical. Should the cluster indices be returned? Only works if +\code{as.df} is \code{FALSE}. Defaults to \code{FALSE}.} + +\item{na.rm}{Logical. Should missing values be removed for the calculation of +pulse metrics? Defaults to \code{FALSE}.} + +\item{as.df}{Logical. Should a data frame be returned? If \code{TRUE}, a data +frame with seven columns ("n", "mean_level", "mean_duration", "total_duration", +"mean_onset", "mean_midpoint", "mean_offset") and the threshold (e.g., \verb{_\{threshold\}}) +will be returned. Defaults to \code{FALSE}.} +} +\value{ +List or data frame with calculated values. +} +\description{ +This function clusters the light data into continuous clusters (pulses) of +light above/below a given threshold. Clustering may be fine-tuned by setting +the minimum length of the clusters and by allowing brief interruptions to be +included in a single cluster, with a specified maximum length of interruption +episodes and proportion of total amount of interruptions to light above +threshold. +} +\details{ +The timeseries is assumed to be regular. Missing values in the +light data will be replaced by 0. +} +\examples{ +# Sample data +data = sample.data.environment \%>\% + dplyr::filter(Id == "Participant") \%>\% + filter_Datetime(length = lubridate::days(1)) \%>\% + dplyr::mutate( + Time = hms::as_hms(Datetime), + ) + +# Time vector as datetime +data \%>\% + dplyr::reframe(pulses_above_threshold(MEDI, Datetime, threshold = 250, as.df = TRUE)) + +# Time vector as hms time +data \%>\% + dplyr::reframe(pulses_above_threshold(MEDI, Time, threshold = 250, as.df = TRUE)) + +# Pulses below threshold +data \%>\% + dplyr::reframe(pulses_above_threshold(MEDI, Datetime, "below", threshold = 250, as.df = TRUE)) + +# Pulses within threshold range +data \%>\% + dplyr::reframe(pulses_above_threshold(MEDI, Datetime, threshold = c(250,1000), as.df = TRUE)) + +} +\references{ +Wilson, J., Reid, K. J., Braun, R. I., Abbott, S. M., & Zee, P. C. +(2018). Habitual light exposure relative to circadian timing in delayed +sleep-wake phase disorder. \emph{Sleep}, 41(11). +\url{https://doi.org/10.1093/sleep/zsy166} +} +\seealso{ +Other metrics: +\code{\link{bright_dark_period}()}, +\code{\link{centroidLE}()}, +\code{\link{disparity_index}()}, +\code{\link{duration_above_threshold}()}, +\code{\link{exponential_moving_average}()}, +\code{\link{frequency_crossing_threshold}()}, +\code{\link{interdaily_stability}()}, +\code{\link{intradaily_variability}()}, +\code{\link{midpointCE}()}, +\code{\link{nvRC}()}, +\code{\link{nvRD}()}, +\code{\link{nvRD_cumulative_response}()}, +\code{\link{period_above_threshold}()}, +\code{\link{threshold_for_duration}()}, +\code{\link{timing_above_threshold}()} +} +\concept{metrics} diff --git a/man/threshold_for_duration.Rd b/man/threshold_for_duration.Rd new file mode 100644 index 0000000..2fe4676 --- /dev/null +++ b/man/threshold_for_duration.Rd @@ -0,0 +1,93 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/metric_threshold_for_duration.R +\name{threshold_for_duration} +\alias{threshold_for_duration} +\title{Find threshold for given duration} +\usage{ +threshold_for_duration( + Light.vector, + Time.vector, + duration, + comparison = c("above", "below"), + epoch = "dominant.epoch", + na.rm = FALSE, + as.df = FALSE +) +} +\arguments{ +\item{Light.vector}{Numeric vector containing the light data.} + +\item{Time.vector}{Vector containing the time data. Can be \link[base]{POSIXct}, +\link[hms]{hms}, \link[lubridate]{duration}, or \link[base]{difftime}.} + +\item{duration}{The duration for which the threshold should be found. Can be either a +\link[lubridate]{duration} or a string. If it is a string, it needs to be a valid +\link[lubridate]{duration} string, e.g., \code{"1 day"} or \code{"10 sec"}.} + +\item{comparison}{String specifying whether light levels above or below the threshold +should be considered. Can be either \code{"above"} (the default) or \code{"below"}.} + +\item{epoch}{The epoch at which the data was sampled. Can be either a +\link[lubridate]{duration} or a string. If it is a string, it needs to be +either \code{"dominant.epoch"} (the default) for a guess based on the data, or a valid +\link[lubridate]{duration} string, e.g., \code{"1 day"} or \code{"10 sec"}.} + +\item{na.rm}{Logical. Should missing values (NA) be removed for the calculation? +Defaults to \code{FALSE}.} + +\item{as.df}{Logical. Should a data frame with be returned? If \code{TRUE}, a data +frame with a single column named \verb{threshold_\{comparison\}_for_\{duration\}} will be returned. +Defaults to \code{FALSE}.} +} +\value{ +Single numeric value or single column data frame. +} +\description{ +This function finds the threshold for which light levels are above/below for +a given duration. This function can be considered as the inverse of +\code{\link{duration_above_threshold}}. +} +\examples{ +N <- 60 +# Dataset with 30 min < 250lx and 30min > 250lx +dataset1 <- + tibble::tibble( + Id = rep("A", N), + Datetime = lubridate::as_datetime(0) + lubridate::minutes(1:N), + MEDI = sample(c(sample(1:249, N / 2, replace = TRUE), + sample(250:1000, N / 2, replace = TRUE))), + ) + +dataset1 \%>\% + dplyr::reframe("Threshold above which for 30 mins" = + threshold_for_duration(MEDI, Datetime, duration = "30 mins")) + +dataset1 \%>\% + dplyr::reframe("Threshold below which for 30 mins" = + threshold_for_duration(MEDI, Datetime, duration = "30 mins", + comparison = "below")) + +dataset1 \%>\% + dplyr::reframe(threshold_for_duration(MEDI, Datetime, duration = "30 mins", + as.df = TRUE)) + +} +\seealso{ +Other metrics: +\code{\link{bright_dark_period}()}, +\code{\link{centroidLE}()}, +\code{\link{disparity_index}()}, +\code{\link{duration_above_threshold}()}, +\code{\link{exponential_moving_average}()}, +\code{\link{frequency_crossing_threshold}()}, +\code{\link{interdaily_stability}()}, +\code{\link{intradaily_variability}()}, +\code{\link{midpointCE}()}, +\code{\link{nvRC}()}, +\code{\link{nvRD}()}, +\code{\link{nvRD_cumulative_response}()}, +\code{\link{period_above_threshold}()}, +\code{\link{pulses_above_threshold}()}, +\code{\link{timing_above_threshold}()} +} +\concept{metrics} diff --git a/man/timing_above_threshold.Rd b/man/timing_above_threshold.Rd index f54f864..1ea388d 100755 --- a/man/timing_above_threshold.Rd +++ b/man/timing_above_threshold.Rd @@ -16,7 +16,8 @@ timing_above_threshold( \arguments{ \item{Light.vector}{Numeric vector containing the light data.} -\item{Time.vector}{Vector containing the time data. Can be numeric, HMS or POSIXct.} +\item{Time.vector}{Vector containing the time data. Can be \link[base]{POSIXct}, +\link[hms]{hms}, \link[lubridate]{duration}, or \link[base]{difftime}.} \item{comparison}{String specifying whether the time above or below threshold should be calculated. Can be either \code{"above"} (the default) or \code{"below"}. If @@ -30,7 +31,7 @@ calculated.} \item{na.rm}{Logical. Should missing values be removed for the calculation? Defaults to \code{FALSE}.} -\item{as.df}{Logical. Should a data frame with be returned? If \code{TRUE}, a data +\item{as.df}{Logical. Should a data frame be returned? If \code{TRUE}, a data frame with three columns (MLiT, FLiT, LLiT) and the threshold (e.g., \verb{MLiT_\{threshold\}}) will be returned. Defaults to \code{FALSE}.} } @@ -83,11 +84,16 @@ Other metrics: \code{\link{centroidLE}()}, \code{\link{disparity_index}()}, \code{\link{duration_above_threshold}()}, +\code{\link{exponential_moving_average}()}, +\code{\link{frequency_crossing_threshold}()}, \code{\link{interdaily_stability}()}, \code{\link{intradaily_variability}()}, \code{\link{midpointCE}()}, \code{\link{nvRC}()}, \code{\link{nvRD}()}, -\code{\link{nvRD_cumulative_response}()} +\code{\link{nvRD_cumulative_response}()}, +\code{\link{period_above_threshold}()}, +\code{\link{pulses_above_threshold}()}, +\code{\link{threshold_for_duration}()} } \concept{metrics} diff --git a/tests/testthat/test-helper.R b/tests/testthat/test-helper.R index f0f55d4..5335545 100644 --- a/tests/testthat/test-helper.R +++ b/tests/testthat/test-helper.R @@ -7,3 +7,25 @@ test_that("count_difftime works", { dataset, c("15s", "16s", "17s", "18s") %>% lubridate::as.duration() ) }) + +test_that("find_clusters works", { + data <- as.logical(c(1,1,1,0,0,1,0,1,1,0,1,1)) + reference.1 <- tibble::tibble( + row_idx = c(1,2,3,8,9,11,12), + is_cluster = c(TRUE,TRUE,TRUE,TRUE,TRUE,TRUE,TRUE), + cluster_idx = c(1,1,1,2,2,3,3), + cluster_start = c(1,1,1,8,8,11,11), + cluster_end = c(3,3,3,9,9,12,12) + ) + reference.2 <- tibble::tibble( + row_idx = c(1,2,3,8,9,10,11,12), + is_cluster = c(TRUE,TRUE,TRUE,TRUE,TRUE,TRUE,TRUE,TRUE), + cluster_idx = c(1,1,1,2,2,2,2,2), + cluster_start = c(1,1,1,8,8,8,8,8), + cluster_end = c(3,3,3,12,12,12,12,12) + ) + expect_equal(find_clusters(data, min.length = 2), reference.1) + expect_equal(find_clusters(data, min.length = 2, max.interrupt = 1), reference.2) + expect_equal(find_clusters(data, min.length = 2, max.interrupt = 1, prop.interrupt = 0.1), reference.1) + expect_equal(find_clusters(data, min.length = 2, max.interrupt = 1, prop.interrupt = 0.25), reference.2) +}) diff --git a/tests/testthat/test-metric_EMA.R b/tests/testthat/test-metric_EMA.R new file mode 100644 index 0000000..6514c66 --- /dev/null +++ b/tests/testthat/test-metric_EMA.R @@ -0,0 +1,35 @@ +test_that("Half life works", { + light <- c(rep(10,90), rep(0,90)) + datetime <- lubridate::as_datetime(0) + lubridate::minutes(1:180) + ema <- exponential_moving_average(light, datetime, decay = "90 mins") + expect_equal(round(ema[90]/2, 1), round(ema[180], 1)) +}) + +test_that("Exponential works", { + light <- c(rep(10,90), rep(0,90)) + datetime <- lubridate::as_datetime(0) + lubridate::minutes(1:180) + ema <- exponential_moving_average(light, datetime, decay = "90 mins") + expect_equal(cor(log(ema[91:180]), 1:90), -1) +}) + +test_that("Handling missing values works", { + light <- c(rep(10,90), rep(NA, 45), rep(0,45)) + datetime <- lubridate::as_datetime(0) + lubridate::minutes(1:180) + expect_warning(ema <- exponential_moving_average(light, datetime, decay = "90 mins"), + "Light data contains missing values! They are replaced by 0.") + expect_equal(cor(log(ema[91:180]), 1:90), -1) + expect_equal(any(is.na(ema)), FALSE) +}) + +test_that("Input checks", { + expect_error(exponential_moving_average("1",lubridate::as_datetime(0)), + "`Light.vector` must be numeric!") + expect_error(exponential_moving_average(1, "20/01/2022 10:00:00"), + "`Time.vector` must be POSIXct, hms, duration, or difftime!") + expect_error(exponential_moving_average(1, 12), + "`Time.vector` must be POSIXct, hms, duration, or difftime!") + expect_error(exponential_moving_average(1, lubridate::as_datetime(0), decay = 1), + "`decay` must either be a duration or a string") + expect_error(exponential_moving_average(1, lubridate::as_datetime(0), epoch = 1), + "`epoch` must either be a duration or a string") +}) \ No newline at end of file diff --git a/tests/testthat/test-metric_barroso.R b/tests/testthat/test-metric_barroso.R new file mode 100644 index 0000000..e888b9b --- /dev/null +++ b/tests/testthat/test-metric_barroso.R @@ -0,0 +1,74 @@ +test_that("Calculation works", { + MEDI = c(rep(0,6), rep(250,4), rep(500,8), rep(0,6)) + Datetime = lubridate::as_datetime(lubridate::dhours(0:23)) + expect_equal( + barroso_lighting_metrics(MEDI, Datetime), + list( + "bright_threshold" = 500 , + "dark_threshold" = 0, + "bright_mean_level" = 500, + "dark_mean_level" = 0, + "bright_cluster" = lubridate::dhours(8), + "dark_cluster" = lubridate::dhours(6), + "circadian_variation" = 1.10 + ) + ) + expect_equal( + barroso_lighting_metrics(MEDI, Datetime, loop = TRUE), + list( + "bright_threshold" = 500 , + "dark_threshold" = 0, + "bright_mean_level" = 500, + "dark_mean_level" = 0, + "bright_cluster" = lubridate::dhours(8), + "dark_cluster" = lubridate::dhours(12), + "circadian_variation" = 1.10 + ) + ) +}) + +test_that("Handling of missing values works", { + MEDI = c(rep(0,6), rep(250,4), rep(500,4), NA, rep(500,3), rep(0,6)) + Datetime = lubridate::as_datetime(lubridate::dhours(0:23)) + expect_equal( + barroso_lighting_metrics(MEDI, Datetime), + list( + "bright_threshold" = as.double(NA), + "dark_threshold" = as.double(NA), + "bright_mean_level" = as.double(NA), + "dark_mean_level" = as.double(NA), + "bright_cluster" = lubridate::dhours(NA), + "dark_cluster" = lubridate::dhours(NA), + "circadian_variation" = as.double(NA) + ) + ) + expect_equal( + barroso_lighting_metrics(MEDI, Datetime, na.rm = TRUE), + list( + "bright_threshold" = 500, + "dark_threshold" = 0, + "bright_mean_level" = 500, + "dark_mean_level" = 0, + "bright_cluster" = lubridate::dhours(4), + "dark_cluster" = lubridate::dhours(6), + "circadian_variation" = 1.15 + ) + ) +}) + +test_that("Return data frame works", { + MEDI = c(rep(0,6), rep(250,4), rep(500,8), rep(0,6)) + Datetime = lubridate::as_datetime(lubridate::dhours(0:23)) + expect_equal( + barroso_lighting_metrics(MEDI, Datetime, as.df = TRUE), + tibble::tibble( + "bright_threshold" = 500 , + "dark_threshold" = 0, + "bright_mean_level" = 500, + "dark_mean_level" = 0, + "bright_cluster" = lubridate::dhours(8), + "dark_cluster" = lubridate::dhours(6), + "circadian_variation" = 1.10 + ) + ) +}) diff --git a/tests/testthat/test-metric_centroidLE.R b/tests/testthat/test-metric_centroidLE.R new file mode 100644 index 0000000..00ddcb5 --- /dev/null +++ b/tests/testthat/test-metric_centroidLE.R @@ -0,0 +1,35 @@ +test_that("Calculation works", { + MEDI = c(rep(0, 6), rep(250, 12), rep(0, 6)) + Datetime = lubridate::as_datetime(lubridate::dhours(0:23), tz = "UTC") + expect_equal(centroidLE(MEDI, Datetime), lubridate::as_datetime(lubridate::dhours(11.5), tz = "UTC")) + expect_equal(centroidLE(MEDI, Datetime, bin.size = "2 hours"), lubridate::as_datetime(lubridate::dhours(11), tz = "UTC")) +}) + +test_that("Missing values", { + MEDI = c(rep(0, 6), rep(250, 12), NA, NA, rep(0, 4)) + Datetime = lubridate::as_datetime(lubridate::dhours(0:23), tz = "UTC") + expect_equal(centroidLE(MEDI, Datetime), lubridate::as_datetime(NA)) + expect_equal(centroidLE(MEDI, Datetime, na.rm = TRUE), lubridate::as_datetime(lubridate::dhours(11.5), tz = "UTC")) +}) + +test_that("Works with different date time units", { + MEDI = c(rep(0, 6), rep(250, 12), rep(0, 6)) + Datetime = lubridate::as_datetime(lubridate::dhours(0:23), tz = "UTC") + HMS = hms::as_hms(Datetime) + Duration = lubridate::as.duration(HMS) + Difftime = Datetime-Datetime[1] + expect_equal(centroidLE(MEDI, HMS), hms::hms(0,30,11)) + expect_equal(centroidLE(MEDI, Duration), lubridate::dhours(11.5)) + expect_equal(centroidLE(MEDI, Difftime), lubridate::as.difftime(lubridate::dhours(11.5))) +}) + +test_that("Return data frame works", { + MEDI = c(rep(0, 6), rep(250, 12), rep(0, 6)) + Datetime = lubridate::as_datetime(lubridate::dhours(0:23), tz = "UTC") + expect_equal( + centroidLE(MEDI, Datetime, as.df = TRUE), + tibble::tibble( + "centroidLE" = lubridate::as_datetime(lubridate::dhours(11.5), tz = "UTC"), + ) + ) +}) \ No newline at end of file diff --git a/tests/testthat/test-metric_disparity_index.R b/tests/testthat/test-metric_disparity_index.R new file mode 100644 index 0000000..50d466a --- /dev/null +++ b/tests/testthat/test-metric_disparity_index.R @@ -0,0 +1,25 @@ + +test_that("Calculation works", { + expect_equal(disparity_index(c(1)), 0) + expect_equal(disparity_index(c(1,1,1)), 0) + expect_equal(round(disparity_index(c(100,100,100,100,50,50,50,50)),2), 0.10) + expect_equal(round(disparity_index(c(100,50,100,50,100,50,100,50)),2), 0.68) +}) + +test_that("Handling of missing values", { + expect_equal(disparity_index(c(1,NA,1)), NA) + expect_equal(round(disparity_index(c(100,100,100,100,NA,50,50,50,50), na.rm = TRUE), 2), 0.10) +}) + +test_that("Return data frame", { + expect_equal(round(disparity_index(c(100,100,100,100,50,50,50,50), as.df = TRUE), 2), + tibble::tibble("disparity_index" = 0.10)) +}) + +test_that("Input checks", { + expect_error(disparity_index(c("0","2","1","0")), "`Light.vector` must be numeric!") + expect_error(disparity_index(c(0,2,1,0,0,2,1,0), na.rm="FALSE"), + "`na.rm` must be logical!") + expect_error(disparity_index(c(0,2,1,0,0,2,1,0), as.df="TRUE"), + "`as.df` must be logical!") +}) diff --git a/tests/testthat/test-metric_duration_above_threshold.R b/tests/testthat/test-metric_duration_above_threshold.R new file mode 100644 index 0000000..c180524 --- /dev/null +++ b/tests/testthat/test-metric_duration_above_threshold.R @@ -0,0 +1,73 @@ +test_that("Calculation works", { + MEDI = c(0,100,99,0,0,101,100,0) + datetime = lubridate::as_datetime((1:8)*60) + + expect_equal( + duration_above_threshold(MEDI, datetime, threshold = 100), lubridate::dminutes(3) + ) + expect_equal( + duration_above_threshold(MEDI, datetime, threshold = 102), lubridate::dminutes(0) + ) + expect_equal( + duration_above_threshold(MEDI, datetime, "below", 100), lubridate::dminutes(7) + ) + expect_equal( + duration_above_threshold(MEDI, datetime, threshold = c(0,99)), lubridate::dminutes(5) + ) +}) + +test_that("Handling of missing values", { + MEDI = c(0,100,99,0,NA,101,100,0) + datetime = lubridate::as_datetime((1:8)*60) + expect_equal( + duration_above_threshold(MEDI, datetime, threshold = 100), lubridate::dminutes(NA) + ) + expect_equal( + duration_above_threshold(MEDI, datetime, threshold = 100, na.rm = TRUE), + lubridate::dminutes(3) + ) +}) + +test_that("Return data frame", { + MEDI = c(0,100,99,0,0,101,100,0) + datetime = lubridate::as_datetime((1:8)*60) + expect_equal( + duration_above_threshold(MEDI, datetime, threshold = 100, as.df = TRUE), + tibble::tibble("duration_above_100" = lubridate::dminutes(3)) + ) + expect_equal( + duration_above_threshold(MEDI, datetime, comparison = "below", threshold = 100, as.df = TRUE), + tibble::tibble("duration_below_100" = lubridate::dminutes(7)) + ) +}) + +test_that("Input checks", { + expect_error( + duration_above_threshold("100", lubridate::as_datetime(1), threshold = 100), + "`Light.vector` must be numeric!" + ) + expect_error( + duration_above_threshold(100, 1, threshold = 100), + "`Time.vector` must be POSIXct, hms, duration, or difftime!" + ) + expect_error( + duration_above_threshold(100, lubridate::as_datetime(1), threshold = "100"), + "`threshold` must be numeric!" + ) + expect_error( + duration_above_threshold(100, lubridate::as_datetime(1), threshold = c(1,2,3)), + "`threshold` must be either one or two values!" + ) + expect_error( + duration_above_threshold(100, lubridate::as_datetime(1), threshold = 100, epoch = 60), + "`epoch` must either be a duration or a string" + ) + expect_error( + duration_above_threshold(100, lubridate::as_datetime(1), threshold = 100, na.rm="FALSE"), + "`na.rm` must be logical!" + ) + expect_error( + duration_above_threshold(100, lubridate::as_datetime(1), threshold = 100, as.df="TRUE"), + "`as.df` must be logical!" + ) +}) diff --git a/tests/testthat/test-metric_frequency_crossing_threshold.R b/tests/testthat/test-metric_frequency_crossing_threshold.R new file mode 100644 index 0000000..3dc43ac --- /dev/null +++ b/tests/testthat/test-metric_frequency_crossing_threshold.R @@ -0,0 +1,28 @@ + +test_that("Calculation works", { + expect_equal(frequency_crossing_threshold(c(0,2,1,0,0,2,1,0), 1), 4) + expect_equal(frequency_crossing_threshold(c(3,2,1,3,3,2,3,1), 1), 0) +}) + +test_that("Handling of missing values", { + expect_equal(frequency_crossing_threshold(c(0,2,1,NA,NA,2,1,0), 1), NA) + expect_equal(frequency_crossing_threshold(c(0,2,1,NA,NA,2,1,0), 1, na.rm = TRUE), 2) +}) + +test_that("Return data frame", { + expect_equal(frequency_crossing_threshold(c(0,2,1,0,0,2,1,0), 1, as.df = TRUE), + tibble::tibble("frequency_crossing_1" = 4)) +}) + +test_that("Input checks", { + expect_error(frequency_crossing_threshold(c("0","2","1","0"), 1), + "`Light.vector` must be numeric!") + expect_error(frequency_crossing_threshold(c(0,2,1,0,0,2,1,0), "1"), + "`threshold` must be numeric!") + expect_error(frequency_crossing_threshold(c(0,2,1,0,0,2,1,0), c(1,2)), + "`threshold` must be one value!") + expect_error(frequency_crossing_threshold(c(0,2,1,0,0,2,1,0), 1, na.rm="FALSE"), + "`na.rm` must be logical!") + expect_error(frequency_crossing_threshold(c(0,2,1,0,0,2,1,0), 1, as.df="TRUE"), + "`as.df` must be logical!") +}) diff --git a/tests/testthat/test-metric_midpointCE.R b/tests/testthat/test-metric_midpointCE.R new file mode 100644 index 0000000..1d651eb --- /dev/null +++ b/tests/testthat/test-metric_midpointCE.R @@ -0,0 +1,34 @@ +test_that("Calculation works", { + MEDI = c(rep(0, 6), rep(250, 12), rep(0, 6)) + Datetime = lubridate::as_datetime(lubridate::dhours(0:23), tz = "UTC") + expect_equal(midpointCE(MEDI, Datetime), lubridate::as_datetime(lubridate::dhours(11), tz = "UTC")) +}) + +test_that("Missing values", { + MEDI = c(rep(0, 6), rep(250, 10), NA, NA, rep(0, 6)) + Datetime = lubridate::as_datetime(lubridate::dhours(0:23), tz = "UTC") + expect_equal(midpointCE(MEDI, Datetime), lubridate::as_datetime(NA)) + expect_equal(midpointCE(MEDI, Datetime, na.rm = TRUE), lubridate::as_datetime(lubridate::dhours(10), tz = "UTC")) +}) + +test_that("Works with different date time units", { + MEDI = c(rep(0, 6), rep(250, 12), rep(0, 6)) + Datetime = lubridate::as_datetime(lubridate::dhours(0:23), tz = "UTC") + HMS = hms::as_hms(Datetime) + Duration = lubridate::as.duration(HMS) + Difftime = Datetime-Datetime[1] + expect_equal(midpointCE(MEDI, HMS), hms::hms(0,0,11)) + expect_equal(midpointCE(MEDI, Duration), lubridate::dhours(11)) + expect_equal(midpointCE(MEDI, Difftime), lubridate::as.difftime(lubridate::dhours(11))) +}) + +test_that("Return data frame works", { + MEDI = c(rep(0, 6), rep(250, 12), rep(0, 6)) + Datetime = lubridate::as_datetime(lubridate::dhours(0:23), tz = "UTC") + expect_equal( + midpointCE(MEDI, Datetime, as.df = TRUE), + tibble::tibble( + "midpointCE" = lubridate::as_datetime(lubridate::dhours(11), tz = "UTC"), + ) + ) +}) \ No newline at end of file diff --git a/tests/testthat/test-metric_period_above_threshold.R b/tests/testthat/test-metric_period_above_threshold.R new file mode 100644 index 0000000..d45d014 --- /dev/null +++ b/tests/testthat/test-metric_period_above_threshold.R @@ -0,0 +1,46 @@ +test_that("Calculation works", { + MEDI = c(10,10,10,0,0,11,12,13,14) + datetime = lubridate::as_datetime(lubridate::dminutes(1:9)) + + expect_equal( + period_above_threshold(MEDI, datetime, threshold = 10), lubridate::dminutes(4) + ) + expect_equal( + period_above_threshold(MEDI, datetime, "below", 10), lubridate::dminutes(5) + ) + expect_equal( + period_above_threshold(MEDI, datetime, threshold = c(12,13)), lubridate::dminutes(2) + ) + expect_equal( + period_above_threshold(MEDI, datetime, threshold = 10, loop = TRUE), lubridate::dminutes(7) + ) +}) + +test_that("Handling of missing values", { + MEDI = c(0,10,10,0,0,NA,12,13,14) + datetime = lubridate::as_datetime(lubridate::dminutes(1:9)) + expect_equal( + period_above_threshold(MEDI, datetime, threshold = 10), lubridate::dminutes(NA) + ) + expect_equal( + period_above_threshold(MEDI, datetime, threshold = 10, na.rm = TRUE), + lubridate::dminutes(3) + ) + expect_equal( + period_above_threshold(MEDI, datetime, "below", 10, na.replace = TRUE), + lubridate::dminutes(5) + ) +}) + +test_that("Return data frame", { + MEDI = c(0,10,10,0,0,11,12,13,14) + datetime = lubridate::as_datetime(lubridate::dminutes(1:9)) + expect_equal( + period_above_threshold(MEDI, datetime, threshold = 10, as.df = TRUE), + tibble::tibble("period_above_10" = lubridate::dminutes(4)) + ) + expect_equal( + period_above_threshold(MEDI, datetime, "below", threshold = 10, as.df = TRUE), + tibble::tibble("period_below_10" = lubridate::dminutes(5)) + ) +}) \ No newline at end of file diff --git a/tests/testthat/test-metric_pulses_above_threshold.R b/tests/testthat/test-metric_pulses_above_threshold.R new file mode 100644 index 0000000..096f166 --- /dev/null +++ b/tests/testthat/test-metric_pulses_above_threshold.R @@ -0,0 +1,179 @@ +test_that("calculation works", { + MEDI = c(10,20,30,0,0,50,60) + Time = lubridate::dminutes(1:7) + expect_equal( + pulses_above_threshold(MEDI,Time,"above",10,"1 mins", "2 mins", 0, as.df = TRUE), + tibble::tibble( + "n_pulses_above_10" = 2, + "mean_level_pulses_above_10" = 37.5, + "mean_duration_pulses_above_10" = lubridate::dminutes(2.5), + "total_duration_pulses_above_10" = lubridate::dminutes(5), + "mean_onset_pulses_above_10" = lubridate::dminutes(3.5), + "mean_midpoint_pulses_above_10" = lubridate::dminutes(4.25), + "mean_offset_pulses_above_10" = lubridate::dminutes(5), + ) + ) + expect_equal( + pulses_above_threshold(MEDI,Time,"above",10,"3 mins", "1 mins", 0, as.df = TRUE), + tibble::tibble( + "n_pulses_above_10" = 1, + "mean_level_pulses_above_10" = 20, + "mean_duration_pulses_above_10" = lubridate::dminutes(3), + "total_duration_pulses_above_10" = lubridate::dminutes(3), + "mean_onset_pulses_above_10" = lubridate::dminutes(1), + "mean_midpoint_pulses_above_10" = lubridate::dminutes(2), + "mean_offset_pulses_above_10" = lubridate::dminutes(3), + ) + ) + expect_equal( + pulses_above_threshold(MEDI,Time,"above",10,"2 mins", "2 mins", 0.4, as.df = TRUE), + tibble::tibble( + "n_pulses_above_10" = 1, + "mean_level_pulses_above_10" = mean(MEDI), + "mean_duration_pulses_above_10" = lubridate::dminutes(7), + "total_duration_pulses_above_10" = lubridate::dminutes(7), + "mean_onset_pulses_above_10" = lubridate::dminutes(1), + "mean_midpoint_pulses_above_10" = lubridate::dminutes(4), + "mean_offset_pulses_above_10" = lubridate::dminutes(7), + ) + ) + expect_equal( + pulses_above_threshold(MEDI,Time,"below",20,"1 mins", "1 mins", 0, as.df = TRUE), + tibble::tibble( + "n_pulses_below_20" = 2, + "mean_level_pulses_below_20" = 7.5, + "mean_duration_pulses_below_20" = lubridate::dminutes(2), + "total_duration_pulses_below_20" = lubridate::dminutes(4), + "mean_onset_pulses_below_20" = lubridate::dminutes(2.5), + "mean_midpoint_pulses_below_20" = lubridate::dminutes(3), + "mean_offset_pulses_below_20" = lubridate::dminutes(3.5), + ) + ) + expect_equal( + pulses_above_threshold(MEDI,Time,"above",c(10,30),"1 mins", "1 mins", 0, as.df = TRUE), + tibble::tibble( + "n_pulses_within_10-30" = 1, + "mean_level_pulses_within_10-30" = 20, + "mean_duration_pulses_within_10-30" = lubridate::dminutes(3), + "total_duration_pulses_within_10-30" = lubridate::dminutes(3), + "mean_onset_pulses_within_10-30" = lubridate::dminutes(1), + "mean_midpoint_pulses_within_10-30" = lubridate::dminutes(2), + "mean_offset_pulses_within_10-30" = lubridate::dminutes(3), + ) + ) +}) + +test_that("Ouput works for different time variables", { + MEDI = c(10,20,30) + Datetime = lubridate::as_datetime(lubridate::dminutes(1:3)) + HMS = hms::as_hms(Datetime) + Duration = lubridate::as.duration(HMS) + Difftime = Datetime-Datetime[1] + x = pulses_above_threshold(MEDI,Datetime,"above",10,"1 mins","1 mins",0,as.df = TRUE) + + + expect_equal( + pulses_above_threshold(MEDI,Datetime,"above",10,"1 mins","1 mins",0,as.df = TRUE), + tibble::tibble( + "n_pulses_above_10" = 1, + "mean_level_pulses_above_10" = 20, + "mean_duration_pulses_above_10" = lubridate::dminutes(3), + "total_duration_pulses_above_10" = lubridate::dminutes(3), + "mean_onset_pulses_above_10" = lubridate::as_datetime(1*60), + "mean_midpoint_pulses_above_10" = lubridate::as_datetime(2*60), + "mean_offset_pulses_above_10" = lubridate::as_datetime(3*60) + ) + ) + + expect_equal( + pulses_above_threshold(MEDI,HMS,"above",10,"1 mins", "1 mins", 0, as.df = TRUE), + tibble::tibble( + "n_pulses_above_10" = 1, + "mean_level_pulses_above_10" = 20, + "mean_duration_pulses_above_10" = lubridate::dminutes(3), + "total_duration_pulses_above_10" = lubridate::dminutes(3), + "mean_onset_pulses_above_10" = hms::as_hms(1*60), + "mean_midpoint_pulses_above_10" = hms::as_hms(2*60), + "mean_offset_pulses_above_10" = hms::as_hms(3*60) + ) + ) + + expect_equal( + pulses_above_threshold(MEDI,Duration,"above",10,"1 mins", "1 mins", 0, as.df = TRUE), + tibble::tibble( + "n_pulses_above_10" = 1, + "mean_level_pulses_above_10" = 20, + "mean_duration_pulses_above_10" = lubridate::dminutes(3), + "total_duration_pulses_above_10" = lubridate::dminutes(3), + "mean_onset_pulses_above_10" = lubridate::dminutes(1), + "mean_midpoint_pulses_above_10" = lubridate::dminutes(2), + "mean_offset_pulses_above_10" = lubridate::dminutes(3), + ) + ) + + expect_equal( + pulses_above_threshold(MEDI,Difftime,"above",10,"1 mins", "1 mins", 0, as.df = TRUE), + tibble::tibble( + "n_pulses_above_10" = 1, + "mean_level_pulses_above_10" = 20, + "mean_duration_pulses_above_10" = lubridate::dminutes(3), + "total_duration_pulses_above_10" = lubridate::dminutes(3), + "mean_onset_pulses_above_10" = lubridate::as.difftime(lubridate::dminutes(0)), + "mean_midpoint_pulses_above_10" = lubridate::as.difftime(lubridate::dminutes(1)), + "mean_offset_pulses_above_10" = lubridate::as.difftime(lubridate::dminutes(2)) + ) + ) +}) + +test_that("calculation works with missing values", { + MEDI = c(10,20,30,NA,NA,50,60) + Time = lubridate::dminutes(1:7) + expect_equal( + pulses_above_threshold(MEDI,Time,"above",10,"1 mins", "1 mins", 0, as.df = TRUE), + tibble::tibble( + "n_pulses_above_10" = 2, + "mean_level_pulses_above_10" = 37.5, + "mean_duration_pulses_above_10" = lubridate::dminutes(2.5), + "total_duration_pulses_above_10" = lubridate::dminutes(5), + "mean_onset_pulses_above_10" = lubridate::dminutes(3.5), + "mean_midpoint_pulses_above_10" = lubridate::dminutes(4.25), + "mean_offset_pulses_above_10" = lubridate::dminutes(5), + ) + ) + expect_equal( + pulses_above_threshold(MEDI,Time,"below",30,"1 mins", "1 mins", 0, as.df = TRUE), + tibble::tibble( + "n_pulses_below_30" = 1, + "mean_level_pulses_below_30" = 20, + "mean_duration_pulses_below_30" = lubridate::dminutes(3), + "total_duration_pulses_below_30" = lubridate::dminutes(3), + "mean_onset_pulses_below_30" = lubridate::dminutes(1), + "mean_midpoint_pulses_below_30" = lubridate::dminutes(2), + "mean_offset_pulses_below_30" = lubridate::dminutes(3), + ) + ) + expect_equal( + pulses_above_threshold(MEDI,Time,"below",60,"1 mins", "2 mins", 1, as.df = TRUE), + tibble::tibble( + "n_pulses_below_60" = 1, + "mean_level_pulses_below_60" = as.double(NA), + "mean_duration_pulses_below_60" = lubridate::dminutes(7), + "total_duration_pulses_below_60" = lubridate::dminutes(7), + "mean_onset_pulses_below_60" = lubridate::dminutes(1), + "mean_midpoint_pulses_below_60" = lubridate::dminutes(4), + "mean_offset_pulses_below_60" = lubridate::dminutes(7), + ) + ) + expect_equal( + pulses_above_threshold(MEDI,Time,"below",60,"1 mins", "2 mins", 1, na.rm=TRUE, as.df = TRUE), + tibble::tibble( + "n_pulses_below_60" = 1, + "mean_level_pulses_below_60" = 34, + "mean_duration_pulses_below_60" = lubridate::dminutes(7), + "total_duration_pulses_below_60" = lubridate::dminutes(7), + "mean_onset_pulses_below_60" = lubridate::dminutes(1), + "mean_midpoint_pulses_below_60" = lubridate::dminutes(4), + "mean_offset_pulses_below_60" = lubridate::dminutes(7), + ) + ) +}) diff --git a/tests/testthat/test-metric_threshold_for_duration.R b/tests/testthat/test-metric_threshold_for_duration.R new file mode 100644 index 0000000..1c7fba3 --- /dev/null +++ b/tests/testthat/test-metric_threshold_for_duration.R @@ -0,0 +1,58 @@ +test_that("Calculation works", { + MEDI <- c(0,100,99,0,0,101,100,0) + datetime <- lubridate::as_datetime((1:8)*60) + expect_equal( + threshold_for_duration(MEDI, datetime, duration = lubridate::dminutes(3)), 100 + ) + expect_equal( + threshold_for_duration(MEDI, datetime, lubridate::dminutes(7), "below"), 100 + ) +}) + +test_that("Handling of missing values", { + MEDI <- c(0,100,99,0,NA,101,100,0) + datetime <- lubridate::as_datetime((1:8)*60) + expect_equal( + threshold_for_duration(MEDI, datetime, lubridate::dminutes(3)), as.double(NA) + ) + expect_equal( + threshold_for_duration(MEDI, datetime, lubridate::dminutes(3), na.rm = TRUE), 100 + ) +}) + +test_that("Return data frame", { + MEDI <- c(0,100,99,0,0,101,100,0) + datetime <- lubridate::as_datetime((1:8)*60) + expect_equal( + threshold_for_duration(MEDI, datetime, lubridate::dminutes(3), as.df = TRUE), + tibble::tibble("threshold_above_for_3_minutes" = 100) + ) + expect_equal( + threshold_for_duration(MEDI, datetime, lubridate::dminutes(7), "below", as.df = TRUE), + tibble::tibble("threshold_below_for_7_minutes" = 100) + ) +}) + +test_that("Input checks", { + expect_error( + threshold_for_duration("100", lubridate::as_datetime(1), lubridate::dminutes(3)), + "`Light.vector` must be numeric!" + ) + expect_error( + threshold_for_duration(100, 1, lubridate::dminutes(3)), + "`Time.vector` must be POSIXct, hms, duration, or difftime!" + ) + expect_error( + threshold_for_duration(100, lubridate::as_datetime(1), 3), + "`duration` must either be duration or a string!" + ) + expect_error( + threshold_for_duration(100, lubridate::as_datetime(1), lubridate::dminutes(3), epoch = 60), + "`epoch` must either be a duration or a string" + ) + expect_error( + threshold_for_duration(100, lubridate::as_datetime(1), lubridate::dminutes(3), as.df = "TRUE"), + "`as.df` must be logical!" + ) +}) + diff --git a/tests/testthat/test-metric_timing_above_threshold.R b/tests/testthat/test-metric_timing_above_threshold.R new file mode 100644 index 0000000..87059ee --- /dev/null +++ b/tests/testthat/test-metric_timing_above_threshold.R @@ -0,0 +1,94 @@ +test_that("Calculation works", { + MEDI = c(rep(1, 6), rep(250, 4), rep(500, 5), rep(750, 4), rep(1, 5)) + Datetime = lubridate::as_datetime(lubridate::dhours(0:23), tz = "UTC") + expect_equal( + timing_above_threshold(MEDI, Datetime, threshold = 250), + list( + "mean" = lubridate::as_datetime(lubridate::dhours(12), tz = "UTC"), + "first" = lubridate::as_datetime(lubridate::dhours(6), tz = "UTC"), + "last" = lubridate::as_datetime(lubridate::dhours(18), tz = "UTC") + ) + ) + expect_equal( + timing_above_threshold(MEDI, Datetime, comparison = "below", threshold = 250), + list( + "mean" = lubridate::as_datetime(lubridate::dhours(10), tz = "UTC"), + "first" = lubridate::as_datetime(lubridate::dhours(0), tz = "UTC"), + "last" = lubridate::as_datetime(lubridate::dhours(23), tz = "UTC") + ) + ) + expect_equal( + timing_above_threshold(MEDI, Datetime, threshold = c(250, 500)), + list( + "mean" = lubridate::as_datetime(lubridate::dhours(10), tz = "UTC"), + "first" = lubridate::as_datetime(lubridate::dhours(6), tz = "UTC"), + "last" = lubridate::as_datetime(lubridate::dhours(14), tz = "UTC") + ) + ) +}) + +test_that("Works with different time representations", { + MEDI = c(rep(1, 6), rep(250, 13), rep(1, 5)) + Datetime = lubridate::as_datetime(lubridate::dhours(0:23), tz = "UTC") + HMS = hms::as_hms(Datetime) + Duration = lubridate::as.duration(HMS) + Difftime = Datetime-Datetime[1] + expect_equal( + timing_above_threshold(MEDI, HMS, threshold = 250), + list( + "mean" = hms::hms(0,0,12), + "first" = hms::hms(0,0,6), + "last" = hms::hms(0,0,18) + ) + ) + expect_equal( + timing_above_threshold(MEDI, Duration, threshold = 250), + list( + "mean" = lubridate::dhours(12), + "first" = lubridate::dhours(6), + "last" = lubridate::dhours(18) + ) + ) + expect_equal( + timing_above_threshold(MEDI, Difftime, threshold = 250), + list( + "mean" = lubridate::as.difftime(lubridate::dhours(12)), + "first" = lubridate::as.difftime(lubridate::dhours(6)), + "last" = lubridate::as.difftime(lubridate::dhours(18)) + ) + ) +}) + +test_that("Handling of missing values works", { + MEDI = c(rep(1, 6), rep(250, 12), NA, rep(1, 5)) + Datetime = lubridate::as_datetime(lubridate::dhours(0:23), tz = "UTC") + expect_equal( + timing_above_threshold(MEDI, Datetime, threshold = 250), + list( + "mean" = lubridate::as_datetime(lubridate::dhours(NA), tz = "UTC"), + "first" = lubridate::as_datetime(lubridate::dhours(6), tz = "UTC"), + "last" = lubridate::as_datetime(lubridate::dhours(NA), tz = "UTC") + ) + ) + expect_equal( + timing_above_threshold(MEDI, Datetime, threshold = 250, na.rm = TRUE), + list( + "mean" = lubridate::as_datetime(lubridate::dhours(11.5), tz = "UTC"), + "first" = lubridate::as_datetime(lubridate::dhours(6), tz = "UTC"), + "last" = lubridate::as_datetime(lubridate::dhours(17), tz = "UTC") + ) + ) +}) + +test_that("Return data frame works", { + MEDI = c(rep(1, 6), rep(250, 13), rep(1, 5)) + Datetime = lubridate::as_datetime(lubridate::dhours(0:23), tz = "UTC") + expect_equal( + timing_above_threshold(MEDI, Datetime, threshold = 250, as.df = TRUE), + tibble::tibble( + "mean_timing_above_250" = lubridate::as_datetime(lubridate::dhours(12), tz = "UTC"), + "first_timing_above_250" = lubridate::as_datetime(lubridate::dhours(6), tz = "UTC"), + "last_timing_above_250" = lubridate::as_datetime(lubridate::dhours(18), tz = "UTC") + ) + ) +})