R/circular_time_mean_sd.R

In mMARCH.AC: Processing of Accelerometry Data with 'GGIR' in mMARCH

#' @title Circular mean, SD, and center adjustment for clock time
#' @description
#' A function for calculating the circular average timing and circular standard deviation
#' of time variables measured on a 24-hour clock. Times larger than 24 hours are wrapped
#' onto the 24-hour scale using modulo arithmetic. The function also returns signed
#' deviations from a specified center, or from the circular mean if center is not provided.
#'
#' @param times \code{numeric} A vector of clock times or elapsed times. Values larger than
#' 24 hours are converted to clock time using \code{x \%/\% 24}. Missing values are ignored.
#' @param unit_to_minute \code{numeric} The ratio of the unit of time and minute. For example, the input unit is hour, the unit2minute = 60. 
#' @param out \code{character} Specify get the mean or sd of the time variables. Default=c("mean","sd") when both mean and sd are calculated.
#'
#' @return A list with:
#' \describe{
#'   \item{\code{mean}}{Circular mean of the input times, in hours.}
#'   \item{\code{sd}}{Circular standard deviation of the input times, in hours.} 
#' }
#'
#' @details
#' The circular mean is computed by minimizing the total squared circular distance on a
#' 24-hour clock. The signed difference is defined on \code{(-12, 12]} hours, so the
#' shortest clockwise or counterclockwise distance is used.
#'
#' @examples 
#' circular_time_mean_sd(c(1, 2, 3, 21, 20, 23))
#' circular_time_mean_sd(c(23, 1, 2) )
#'
#' @export
 
 

 
circular_time_mean_sd  <- function(times,
                                  unit_to_minute = 60,
                                  out = c("mean", "sd") ) {

  out <- match.arg(out, c("mean", "sd"), several.ok = TRUE)

  # Remove missing values
  times <- times[!is.na(times)]

  # No data
  if (length(times) == 0) {
    result <- c(mean = NA_real_, sd = NA_real_)
    return(result[out])
  }

  # Convert to minutes and wrap onto a 24-hour circle [0,1440)
  day_minutes = 24 * 60
  times <- (times * unit_to_minute) %% day_minutes
 

  # One value
  if (length(times) == 1) {
    result <- c(mean = times[1] / unit_to_minute, sd = 0)
    return(result[out])
  }

  # Candidate means: every minute in a day
  candidate_means <- 0:day_minutes

  # Pairwise differences: rows = observed times, columns = candidate means
  diff_matrix <- outer(times, candidate_means, "-")
  
  # Circular squared distance
  circular_sq_distance <- pmin(
    diff_matrix^2,
    (day_minutes - abs(diff_matrix))^2
  )

  # Total distance for each candidate mean
  total_distance <- colSums(circular_sq_distance)

  # Best mean
  best_index <- which.min(total_distance)
  mean_minutes <- candidate_means[best_index]

  # Circular SD
  sd_minutes <- sqrt(total_distance[best_index] / (length(times) - 1))
 
  result <- c(
    mean = mean_minutes / unit_to_minute,
    sd   = sd_minutes / unit_to_minute
  )

  return(result[out])
}



 




#' Center clock times on a circular 24-hour scale
#'
#' Computes the shortest signed circular distance between clock times
#' and a specified center hour on a 24-hour clock.
#'
#' @param hours Numeric vector of hours.
#' @param center Center hour. If NULL, circular mean is used by minimizing the total squared circular distance on a
#' 24-hour clock.
#' @importFrom dplyr %>% mutate 
#'
#' @return A data frame containing centered circular distances.
#'
#' @details
#' The signed difference is defined on \code{(-12, 12]} hours, so the
#' shortest clockwise or counterclockwise distance is used.

#' @examples 
#' circular_hour_center(c(0:24),8)  
#' circular_hour_center(c(NA,0:24),21) 
#'
#' @export


circular_hour_center <- function(hours,center=NULL){ 

  # adjusted time use adjust_center 
  if (is.null( center)) center <- circular_time_mean_sd(hours,  unit_to_minute = 60, out = c("mean" ) )  

   # case1: 24> hour  > center >0 ,clockwise = hour - center > 0
   # case2:  center >  hour  
       #2a:  clockwise = hour - center = 2-4 or 13-16  < 12 hours gap
       #2b:  clockwise = hour - center (cross midnigt) = 2-22  = -20,  but clockwise distance is just 26-22=4 hours, output=4

  hours<-hours %% 24

   clockwise <- ifelse(
      abs(hours - center) <= abs(24 + hours - center),
      hours - center,
      24 + hours - center
   )
   counterclockwise= 24 - abs(clockwise) 
   df <- data.frame(orig = hours, clockwise, counterclockwise)
   df <- df %>%
        mutate( dir_min = ifelse(abs(clockwise) <=12 & clockwise>=0, "clockwise", "counterclockwise"), 
        time_adjust = ifelse(
                         abs(clockwise) == 12,
                          12,
                         ifelse(
                           abs(clockwise) < 12,
                             clockwise,
                            -abs(counterclockwise)
                           )
                       )
  )  
  df
}