R/variance_bounds.R

In wavecov: Overlap-Aware Workforce Wave Comparisons

#' Cauchy-based bounds for the variance of the change in means
#'
#' Compute generic and overlap-aware Cauchy bounds on the variance of the
#' observed change in wave means when the covariance between wave means is
#' unknown.
#'
#' @param N0 Wave 0 workforce size.
#' @param N1 Wave 1 workforce size.
#' @param V0_bar Variance of the wave-0 mean.
#' @param V1_bar Variance of the wave-1 mean.
#' @param nS Number of stayers (optional, for overlap-aware bounds).
#' @param V0S_ind Stayer individual-level variance at wave 0 (optional).
#' @param V1S_ind Stayer individual-level variance at wave 1 (optional).
#' @param sigma2S Pooled individual-level stayer variance (optional).
#' @param cov_sign Covariance sign restriction: "unrestricted" or "nonnegative".
#' @param overlap_mode Whether to use generic bounds or overlap-aware tightening.
#'
#' @return A list with covariance and variance bounds and associated metadata.
#' @keywords internal
variance_bounds <- function(N0, N1,
                            V0_bar, V1_bar,
                            nS = NULL,
                            V0S_ind = NULL,
                            V1S_ind = NULL,
                            sigma2S = NULL,
                            cov_sign = c("unrestricted", "nonnegative"),
                            overlap_mode = c("generic", "overlap_aware")) {
  cov_sign     <- base::match.arg(cov_sign)
  overlap_mode <- base::match.arg(overlap_mode)
  
  if (!is.numeric(V0_bar) || !is.numeric(V1_bar) || V0_bar <= 0 || V1_bar <= 0) {
    stop("V0_bar and V1_bar must be positive variances of the wave means.")
  }
  
  # Generic Cauchy bound on covariance of wave means
  generic_cov_max <- base::sqrt(V0_bar * V1_bar)
  
  overlap_used <- FALSE
  if (identical(overlap_mode, "overlap_aware")) {
    if (is.null(nS) || !is.numeric(nS) || nS < 0) {
      stop("Overlap-aware bounds require nonnegative nS.")
    }
    
    # Allow either V0S_ind/V1S_ind or a pooled sigma2S
    if (!is.null(V0S_ind) && !is.null(V1S_ind)) {
      if (V0S_ind <= 0 || V1S_ind <= 0) {
        stop("V0S_ind and V1S_ind must be positive stayer variances.")
      }
      V0S_use <- V0S_ind
      V1S_use <- V1S_ind
      overlap_label <- "overlap_aware"
    } else if (!is.null(sigma2S)) {
      if (sigma2S <= 0) {
        stop("sigma2S must be a positive pooled stayer variance.")
      }
      V0S_use <- sigma2S
      V1S_use <- sigma2S
      overlap_label <- "overlap_aware_pooled"
    } else {
      stop("Overlap-aware bounds require V0S_ind/V1S_ind or pooled sigma2S.")
    }
    
    overlap_cov_max <- nS * base::sqrt(V0S_use * V1S_use) / (N0 * N1)
    cov_max <- base::min(generic_cov_max, overlap_cov_max)
    overlap_used <- TRUE
  } else {
    cov_max <- generic_cov_max
    overlap_label <- "generic"
  }
  
  # Covariance sign restriction
  if (identical(cov_sign, "unrestricted")) {
    cov_lower <- -cov_max
    cov_upper <-  cov_max
  } else { # nonnegative
    cov_lower <- 0
    cov_upper <- cov_max
  }
  
  var_lower <- V0_bar + V1_bar - 2 * cov_upper
  var_upper <- V0_bar + V1_bar - 2 * cov_lower
  if (var_lower < 0) var_lower <- 0
  
  list(
    cov_lower    = cov_lower,
    cov_upper    = cov_upper,
    var_lower    = var_lower,
    var_upper    = var_upper,
    se_lower     = base::sqrt(var_lower),
    se_upper     = base::sqrt(var_upper),
    cov_sign     = cov_sign,
    overlap_mode = if (overlap_used) overlap_label else "generic"
  )
}