diff --git a/NAMESPACE b/NAMESPACE index 7164f0b..68dabf4 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -1,9 +1,11 @@ # Generated by roxygen2: do not edit by hand +S3method("[",ps_calib) S3method("[",ps_trim) S3method("[",ps_trim_matrix) S3method("[",ps_trunc) S3method("[",ps_trunc_matrix) +S3method(Summary,ps_calib) S3method(Summary,ps_trim) S3method(Summary,ps_trunc) S3method(Summary,psw) @@ -11,10 +13,13 @@ S3method(anyDuplicated,ps_trim) S3method(anyDuplicated,ps_trunc) S3method(anyDuplicated,psw) S3method(as.data.frame,ipw) +S3method(as.double,ps_calib) S3method(diff,ps_trim) S3method(diff,ps_trunc) S3method(diff,psw) S3method(is_ps_calibrated,default) +S3method(is_ps_calibrated,ps_calib) +S3method(is_ps_calibrated,ps_calib_matrix) S3method(is_ps_calibrated,psw) S3method(is_ps_trimmed,default) S3method(is_ps_trimmed,ps_trim) @@ -35,17 +40,22 @@ S3method(is_unit_truncated,default) S3method(is_unit_truncated,ps_trunc) S3method(is_unit_truncated,ps_trunc_matrix) S3method(is_unit_truncated,psw) +S3method(max,ps_calib) S3method(max,ps_trim) S3method(max,ps_trunc) S3method(max,psw) +S3method(median,ps_calib) S3method(median,ps_trim) S3method(median,ps_trunc) S3method(median,psw) +S3method(min,ps_calib) S3method(min,ps_trim) S3method(min,ps_trunc) S3method(min,psw) S3method(na.omit,ps_trim) S3method(print,ipw) +S3method(print,ps_calib) +S3method(print,ps_calib_matrix) S3method(print,ps_trim_matrix) S3method(print,ps_trunc_matrix) S3method(ps_trim,data.frame) @@ -59,6 +69,7 @@ S3method(ps_trunc,ps_trunc) S3method(quantile,ps_trim) S3method(quantile,ps_trunc) S3method(quantile,psw) +S3method(range,ps_calib) S3method(range,ps_trim) S3method(range,ps_trunc) S3method(range,psw) @@ -67,12 +78,18 @@ S3method(sort,ps_trunc) S3method(summary,ps_trim) S3method(summary,ps_trunc) S3method(summary,psw) +S3method(vec_arith,ps_calib) S3method(vec_arith,ps_trim) S3method(vec_arith,ps_trunc) S3method(vec_arith,psw) +S3method(vec_arith.numeric,ps_calib) S3method(vec_arith.numeric,ps_trim) S3method(vec_arith.numeric,ps_trunc) S3method(vec_arith.numeric,psw) +S3method(vec_arith.ps_calib,MISSING) +S3method(vec_arith.ps_calib,default) +S3method(vec_arith.ps_calib,numeric) +S3method(vec_arith.ps_calib,ps_calib) S3method(vec_arith.ps_trim,MISSING) S3method(vec_arith.ps_trim,default) S3method(vec_arith.ps_trim,integer) @@ -92,12 +109,15 @@ S3method(vec_arith.psw,psw) S3method(vec_cast,character.ps_trim) S3method(vec_cast,character.ps_trunc) S3method(vec_cast,character.psw) +S3method(vec_cast,double.ps_calib) S3method(vec_cast,double.ps_trim) S3method(vec_cast,double.ps_trunc) S3method(vec_cast,double.psw) S3method(vec_cast,integer.ps_trim) S3method(vec_cast,integer.ps_trunc) S3method(vec_cast,integer.psw) +S3method(vec_cast,ps_calib.double) +S3method(vec_cast,ps_calib.ps_calib) S3method(vec_cast,ps_trim.double) S3method(vec_cast,ps_trim.integer) S3method(vec_cast,ps_trim.ps_trim) @@ -113,83 +133,103 @@ S3method(vec_cast,psw.integer) S3method(vec_cast,psw.ps_trim) S3method(vec_cast,psw.ps_trunc) S3method(vec_cast,psw.psw) +S3method(vec_math,ps_calib) S3method(vec_math,ps_trim) S3method(vec_math,ps_trunc) S3method(vec_math,psw) S3method(vec_ptype2,character.psw) +S3method(vec_ptype2,double.ps_calib) S3method(vec_ptype2,double.ps_trim) S3method(vec_ptype2,double.ps_trunc) S3method(vec_ptype2,double.psw) S3method(vec_ptype2,integer.ps_trim) S3method(vec_ptype2,integer.ps_trunc) S3method(vec_ptype2,integer.psw) +S3method(vec_ptype2,ps_calib.double) +S3method(vec_ptype2,ps_calib.ps_calib) +S3method(vec_ptype2,ps_calib.ps_trim) +S3method(vec_ptype2,ps_calib.ps_trunc) +S3method(vec_ptype2,ps_calib.psw) S3method(vec_ptype2,ps_trim.double) S3method(vec_ptype2,ps_trim.integer) +S3method(vec_ptype2,ps_trim.ps_calib) S3method(vec_ptype2,ps_trim.ps_trim) S3method(vec_ptype2,ps_trim.ps_trunc) S3method(vec_ptype2,ps_trim.psw) S3method(vec_ptype2,ps_trunc.double) S3method(vec_ptype2,ps_trunc.integer) +S3method(vec_ptype2,ps_trunc.ps_calib) S3method(vec_ptype2,ps_trunc.ps_trim) S3method(vec_ptype2,ps_trunc.ps_trunc) S3method(vec_ptype2,ps_trunc.psw) S3method(vec_ptype2,psw.character) S3method(vec_ptype2,psw.double) S3method(vec_ptype2,psw.integer) +S3method(vec_ptype2,psw.ps_calib) S3method(vec_ptype2,psw.ps_trim) S3method(vec_ptype2,psw.ps_trunc) S3method(vec_ptype2,psw.psw) +S3method(vec_ptype_abbr,ps_calib) S3method(vec_ptype_abbr,ps_trim) S3method(vec_ptype_abbr,ps_trunc) S3method(vec_ptype_abbr,psw) +S3method(vec_ptype_full,ps_calib) S3method(vec_ptype_full,ps_trim) S3method(vec_ptype_full,ps_trunc) S3method(vec_ptype_full,psw) S3method(vec_restore,ps_trim) S3method(vec_restore,ps_trunc) S3method(vec_restore,psw) +S3method(wt_atc,ps_calib) S3method(wt_atc,ps_trim) S3method(wt_atc,ps_trunc) S3method(wt_ate,data.frame) S3method(wt_ate,default) S3method(wt_ate,glm) S3method(wt_ate,numeric) +S3method(wt_ate,ps_calib) S3method(wt_ate,ps_trim) S3method(wt_ate,ps_trunc) S3method(wt_atm,data.frame) S3method(wt_atm,default) S3method(wt_atm,glm) S3method(wt_atm,numeric) +S3method(wt_atm,ps_calib) S3method(wt_atm,ps_trim) S3method(wt_atm,ps_trunc) S3method(wt_ato,data.frame) S3method(wt_ato,default) S3method(wt_ato,glm) S3method(wt_ato,numeric) +S3method(wt_ato,ps_calib) S3method(wt_ato,ps_trim) S3method(wt_ato,ps_trunc) S3method(wt_att,data.frame) S3method(wt_att,default) S3method(wt_att,glm) S3method(wt_att,numeric) +S3method(wt_att,ps_calib) S3method(wt_att,ps_trim) S3method(wt_att,ps_trunc) S3method(wt_atu,data.frame) S3method(wt_atu,default) S3method(wt_atu,glm) S3method(wt_atu,numeric) +S3method(wt_atu,ps_calib) S3method(wt_atu,ps_trim) S3method(wt_atu,ps_trunc) S3method(wt_cens,data.frame) S3method(wt_cens,default) S3method(wt_cens,glm) S3method(wt_cens,numeric) +S3method(wt_cens,ps_calib) S3method(wt_cens,ps_trim) S3method(wt_cens,ps_trunc) S3method(wt_entropy,data.frame) S3method(wt_entropy,default) S3method(wt_entropy,glm) S3method(wt_entropy,numeric) +S3method(wt_entropy,ps_calib) S3method(wt_entropy,ps_trim) S3method(wt_entropy,ps_trunc) export("estimand<-") diff --git a/R/ps_calibration.R b/R/ps_calibration.R index f96d904..b69b5bf 100644 --- a/R/ps_calibration.R +++ b/R/ps_calibration.R @@ -30,9 +30,8 @@ #' If not provided, `ps_calibrate()` will attempt to automatically determine the coding. #' @param .treated `r lifecycle::badge("deprecated")` Use `.focal_level` instead. #' @param .untreated `r lifecycle::badge("deprecated")` Use `.reference_level` instead. -#' @param estimand Character indicating the estimand type. #' -#' @return A calibrated propensity score object (`psw`) +#' @return A calibrated propensity score object (`ps_calib`) #' #' @examples #' # Generate example data @@ -56,7 +55,6 @@ ps_calibrate <- function( smooth = TRUE, .focal_level = NULL, .reference_level = NULL, - estimand = NULL, .treated = NULL, .untreated = NULL ) { @@ -80,15 +78,21 @@ ps_calibrate <- function( ) } - if (is_causal_wt(ps) && is_ps_calibrated(ps)) { + if (is_ps_calibrated(ps)) { abort( "`ps` is already calibrated. Cannot calibrate already calibrated propensity scores.", error_class = "propensity_already_calibrated_error" ) } - # Extract numeric values for comparison if ps is a psw object - ps_numeric <- if (is_psw(ps)) as.numeric(ps) else ps + # Extract numeric values for comparison if ps is a causal weight object + ps_numeric <- if (is_causal_wt(ps)) { + as.numeric(ps) + } else if (inherits(ps, "ps_calib")) { + vec_data(ps) + } else { + ps + } if (any(ps_numeric < 0 | ps_numeric > 1, na.rm = TRUE)) { abort( @@ -111,22 +115,8 @@ ps_calibrate <- function( ) } - # Extract attributes from causal weight objects if applicable - if (is_causal_wt(ps)) { - if (is.null(estimand)) { - estimand <- estimand(ps) - } - stabilized <- is_stabilized(ps) - trimmed <- is_ps_trimmed(ps) - truncated <- is_ps_truncated(ps) - } else { - if (is.null(estimand)) { - estimand <- "unknown" - } - stabilized <- FALSE - trimmed <- FALSE - truncated <- FALSE - } + # ps_calib objects should not have weight-specific attributes like estimand + # Only store calibration-specific metadata # Handle NA values na_idx <- is.na(ps) | is.na(.exposure) @@ -202,6 +192,8 @@ ps_calibrate <- function( calib_ps[na_idx] <- NA } } else if (method == "isoreg") { + # Smoothing is not applicable to isotonic regression + smooth <- FALSE # Isotonic regression calibration if (any(na_idx)) { # Work with non-NA values only @@ -244,17 +236,43 @@ ps_calibrate <- function( calib_ps[!na_idx] <- pmax(0, pmin(1, calib_ps[!na_idx])) } - # Return calibrated scores as a psw object with preserved attributes - psw( + # Return calibrated scores as a ps_calib object + new_ps_calib( x = calib_ps, - estimand = estimand, - stabilized = stabilized, - trimmed = trimmed, - truncated = truncated, - calibrated = TRUE + ps_calib_meta = list( + method = method, + smooth = smooth + ) ) } +#' Constructor for ps_calib objects +#' @noRd +new_ps_calib <- function(x, ps_calib_meta = list()) { + if (is.matrix(x)) { + # For matrices, we don't use vctrs + structure( + x, + ps_calib_meta = ps_calib_meta, + class = c("ps_calib_matrix", "ps_calib", "matrix") + ) + } else { + vec_assert(x, ptype = double()) + new_vctr( + x, + ps_calib_meta = ps_calib_meta, + class = "ps_calib", + inherit_base_type = TRUE + ) + } +} + +#' Extract metadata from ps_calib object +#' @noRd +ps_calib_meta <- function(x) { + attr(x, "ps_calib_meta") +} + #' Check if object is calibrated #' #' @description @@ -277,3 +295,304 @@ is_ps_calibrated.default <- function(x) { is_ps_calibrated.psw <- function(x) { isTRUE(attr(x, "calibrated")) } + +#' @export +is_ps_calibrated.ps_calib <- function(x) { + TRUE +} + +#' @export +is_ps_calibrated.ps_calib_matrix <- function(x) { + TRUE +} + +# vctrs machinery for ps_calib + +#' @export +vec_ptype_abbr.ps_calib <- function(x, ...) "ps_calib" + +#' @export +vec_ptype_full.ps_calib <- function(x, ...) { + meta <- ps_calib_meta(x) + paste("ps_calib", meta$method) +} + +#' @export +#' @method vec_arith ps_calib +vec_arith.ps_calib <- function(op, x, y, ...) { + UseMethod("vec_arith.ps_calib", y) +} + +#' @export +#' @method vec_arith.ps_calib default +vec_arith.ps_calib.default <- function(op, x, y, ...) { + vec_arith_base(op, x, y) +} + +#' @export +#' @method vec_arith.ps_calib ps_calib +vec_arith.ps_calib.ps_calib <- function(op, x, y, ...) { + vec_arith_base(op, x, y) +} + +#' @export +#' @method vec_arith.ps_calib MISSING +vec_arith.ps_calib.MISSING <- function(op, x, y, ...) { + switch( + op, + `-` = -1 * vec_data(x), # Returns numeric + `+` = vec_data(x), # Returns numeric + stop_incompatible_op(op, x, y) + ) +} + +#' @export +#' @method vec_arith.ps_calib numeric +vec_arith.ps_calib.numeric <- function(op, x, y, ...) { + vec_arith_base(op, x, y) +} + +#' @export +#' @method vec_arith.numeric ps_calib +vec_arith.numeric.ps_calib <- function(op, x, y, ...) { + vec_arith_base(op, x, y) +} + +#' @export +vec_ptype2.ps_calib.ps_calib <- function(x, y, ...) { + x_meta <- ps_calib_meta(x) + y_meta <- ps_calib_meta(y) + + # Check if calibration parameters match + if ( + !identical(x_meta$method, y_meta$method) || + !identical(x_meta$smooth, y_meta$smooth) + ) { + warn_incompatible_metadata( + x, + y, + "different calibration parameters" + ) + return(double()) + } + + # If parameters match, return ps_calib prototype + new_ps_calib(double(), ps_calib_meta = x_meta) +} + +#' @export +vec_ptype2.ps_calib.double <- function(x, y, ...) { + warn_class_downgrade("ps_calib") + double() +} + +#' @export +vec_ptype2.double.ps_calib <- function(x, y, ...) { + warn_class_downgrade("ps_calib") + double() +} + +#' @export +vec_cast.ps_calib.ps_calib <- function(x, to, ...) { + # Check if metadata matches + x_meta <- ps_calib_meta(x) + to_meta <- ps_calib_meta(to) + + if ( + !identical(x_meta$method, to_meta$method) || + !identical(x_meta$smooth, to_meta$smooth) + ) { + vctrs::stop_incompatible_cast(x, to, x_arg = "", to_arg = "") + } + + # Return x as-is if metadata matches + x +} + +#' @export +vec_cast.double.ps_calib <- function(x, to, ...) { + # degrade to numeric + vec_data(x) +} + +#' @export +as.double.ps_calib <- function(x, ...) { + vec_data(x) +} + +#' @export +vec_cast.ps_calib.double <- function(x, to, ...) { + # create a default ps_calib + new_ps_calib( + x, + ps_calib_meta = list( + method = "unknown", + smooth = FALSE + ) + ) +} + +#' @export +vec_ptype2.psw.ps_calib <- function(x, y, ...) { + warn_class_downgrade(c("psw", "ps_calib")) + double() +} + +#' @export +vec_ptype2.ps_calib.psw <- function(x, y, ...) { + warn_class_downgrade(c("ps_calib", "psw")) + double() +} + +#' @export +vec_ptype2.ps_calib.ps_trim <- function(x, y, ...) { + warn_class_downgrade(c("ps_calib", "ps_trim")) + double() +} + +#' @export +vec_ptype2.ps_trim.ps_calib <- function(x, y, ...) { + warn_class_downgrade(c("ps_trim", "ps_calib")) + double() +} + +#' @export +vec_ptype2.ps_calib.ps_trunc <- function(x, y, ...) { + warn_class_downgrade(c("ps_calib", "ps_trunc")) + double() +} + +#' @export +vec_ptype2.ps_trunc.ps_calib <- function(x, y, ...) { + warn_class_downgrade(c("ps_trunc", "ps_calib")) + double() +} + +#' @export +vec_math.ps_calib <- function(.fn, .x, ...) { + vec_math_base(.fn, vec_data(.x), ...) +} + +#' @export +Summary.ps_calib <- function(..., na.rm = FALSE) { + .fn <- .Generic + args <- list(...) + numeric_args <- lapply(args, vec_data) + do.call(.fn, c(numeric_args, list(na.rm = na.rm))) +} + +#' @export +min.ps_calib <- function(..., na.rm = FALSE) { + args <- list(...) + numeric_args <- lapply(args, vec_data) + do.call("min", c(numeric_args, list(na.rm = na.rm))) +} + +#' @export +max.ps_calib <- function(..., na.rm = FALSE) { + args <- list(...) + numeric_args <- lapply(args, vec_data) + do.call("max", c(numeric_args, list(na.rm = na.rm))) +} + +#' @export +range.ps_calib <- function(..., na.rm = FALSE) { + args <- list(...) + numeric_args <- lapply(args, vec_data) + do.call("range", c(numeric_args, list(na.rm = na.rm))) +} + +#' @export +median.ps_calib <- function(x, na.rm = FALSE, ...) { + median(vec_data(x), na.rm = na.rm, ...) +} + +#' @export +`[.ps_calib` <- function(x, i, ...) { + # If i is missing, just call NextMethod + if (missing(i)) { + return(NextMethod()) + } + + # Get the subset of data using NextMethod to handle vctrs subsetting + result <- NextMethod() + + # Preserve metadata + attr(result, "ps_calib_meta") <- ps_calib_meta(x) + + result +} + +#' @export +print.ps_calib <- function(x, ..., n = NULL) { + meta <- ps_calib_meta(x) + n_vals <- length(x) + + # Create header + cat(sprintf( + "\n", + n_vals, + meta$method, + if (isTRUE(meta$smooth)) " (smooth)" else "" + )) + + # Determine how many values to show + if (is.null(n)) { + n_print <- getOption("propensity.print_max", default = 10) + } else { + n_print <- n + } + + # Show values + if (is.infinite(n_print) || n_print >= n_vals) { + print(vec_data(x)) + } else { + n_show <- min(n_print, n_vals) + print(vec_data(x)[seq_len(n_show)]) + + if (n_vals > n_show) { + cat("# ... with", n_vals - n_show, "more values\n") + } + } + + invisible(x) +} + +#' @export +print.ps_calib_matrix <- function(x, ..., n = NULL) { + meta <- ps_calib_meta(x) + n_rows <- nrow(x) + k <- ncol(x) + + # Create header + cat(sprintf( + "\n", + n_rows, + k, + meta$method, + if (isTRUE(meta$smooth)) " (smooth)" else "" + )) + + # Determine how many rows to show + if (is.null(n)) { + n_print <- getOption("propensity.print_max", default = 10) + } else { + n_print <- n + } + + # Show all rows if n is Inf or very large + if (is.infinite(n_print) || n_print >= n_rows) { + print(unclass(x)) + } else { + # Show first n_print rows + n_show <- min(n_print, n_rows) + x_sub <- x[seq_len(n_show), , drop = FALSE] + print(unclass(x_sub)) + + if (n_rows > n_show) { + cat("# ... with", n_rows - n_show, "more rows\n") + } + } + + invisible(x) +} diff --git a/R/weights-utils.R b/R/weights-utils.R index 2750022..3f7353d 100644 --- a/R/weights-utils.R +++ b/R/weights-utils.R @@ -430,7 +430,7 @@ calculate_weight_from_modified_ps <- function( .propensity, .exposure, weight_fn, - modification_type = c("trim", "trunc"), + modification_type = c("trim", "trunc", "calib"), ... ) { modification_type <- rlang::arg_match(modification_type) @@ -441,7 +441,12 @@ calculate_weight_from_modified_ps <- function( } # Handle matrix or vector propensity scores - if (inherits(.propensity, c("ps_trim_matrix", "ps_trunc_matrix"))) { + if ( + inherits( + .propensity, + c("ps_trim_matrix", "ps_trunc_matrix", "ps_calib_matrix") + ) + ) { # For matrix propensity scores, pass them directly # The weight function should handle the matrix appropriately base_wt <- weight_fn( @@ -467,10 +472,15 @@ calculate_weight_from_modified_ps <- function( estimand(base_wt) <- paste0(old_est, "; trimmed") attr(base_wt, "trimmed") <- TRUE attr(base_wt, "ps_trim_meta") <- attr(.propensity, "ps_trim_meta") - } else { + } else if (modification_type == "trunc") { estimand(base_wt) <- paste0(estimand(base_wt), "; truncated") attr(base_wt, "truncated") <- TRUE attr(base_wt, "ps_trunc_meta") <- ps_trunc_meta(.propensity) + } else { + # calib + estimand(base_wt) <- paste0(estimand(base_wt), "; calibrated") + attr(base_wt, "calibrated") <- TRUE + attr(base_wt, "ps_calib_meta") <- ps_calib_meta(.propensity) } base_wt diff --git a/R/weights.R b/R/weights.R index 8a6988d..e6b8f83 100644 --- a/R/weights.R +++ b/R/weights.R @@ -2036,3 +2036,195 @@ calculate_categorical_weights <- function( weights } + +# ps_calib methods ---- + +#' @export +wt_ate.ps_calib <- function( + .propensity, + .exposure, + .sigma = NULL, + exposure_type = c("auto", "binary", "categorical", "continuous"), + .focal_level = NULL, + .reference_level = NULL, + stabilize = FALSE, + stabilization_score = NULL, + ..., + .treated = NULL, + .untreated = NULL +) { + calculate_weight_from_modified_ps( + .propensity = .propensity, + .exposure = .exposure, + weight_fn = wt_ate.numeric, + modification_type = "calib", + .sigma = .sigma, + exposure_type = exposure_type, + .focal_level = .focal_level, + .reference_level = .reference_level, + .treated = .treated, + .untreated = .untreated, + stabilize = stabilize, + stabilization_score = stabilization_score, + ... + ) +} + +#' @export +wt_att.ps_calib <- function( + .propensity, + .exposure, + exposure_type = c("auto", "binary", "categorical"), + .focal_level = NULL, + .reference_level = NULL, + ..., + .treated = NULL, + .untreated = NULL +) { + calculate_weight_from_modified_ps( + .propensity = .propensity, + .exposure = .exposure, + weight_fn = wt_att.numeric, + modification_type = "calib", + exposure_type = exposure_type, + .focal_level = .focal_level, + .reference_level = .reference_level, + .treated = .treated, + .untreated = .untreated, + ... + ) +} + +#' @export +wt_atu.ps_calib <- function( + .propensity, + .exposure, + exposure_type = c("auto", "binary", "categorical"), + .focal_level = NULL, + .reference_level = NULL, + ..., + .treated = NULL, + .untreated = NULL +) { + calculate_weight_from_modified_ps( + .propensity = .propensity, + .exposure = .exposure, + weight_fn = wt_atu.numeric, + modification_type = "calib", + exposure_type = exposure_type, + .focal_level = .focal_level, + .reference_level = .reference_level, + .treated = .treated, + .untreated = .untreated, + ... + ) +} + +#' @export +wt_atm.ps_calib <- function( + .propensity, + .exposure, + exposure_type = c("auto", "binary", "categorical"), + .focal_level = NULL, + .reference_level = NULL, + ..., + .treated = NULL, + .untreated = NULL +) { + calculate_weight_from_modified_ps( + .propensity = .propensity, + .exposure = .exposure, + weight_fn = wt_atm.numeric, + modification_type = "calib", + exposure_type = exposure_type, + .focal_level = .focal_level, + .reference_level = .reference_level, + .treated = .treated, + .untreated = .untreated, + ... + ) +} + +#' @export +wt_ato.ps_calib <- function( + .propensity, + .exposure, + exposure_type = c("auto", "binary", "categorical"), + .focal_level = NULL, + .reference_level = NULL, + ..., + .treated = NULL, + .untreated = NULL +) { + calculate_weight_from_modified_ps( + .propensity = .propensity, + .exposure = .exposure, + weight_fn = wt_ato.numeric, + modification_type = "calib", + exposure_type = exposure_type, + .focal_level = .focal_level, + .reference_level = .reference_level, + .treated = .treated, + .untreated = .untreated, + ... + ) +} + +#' @export +wt_entropy.ps_calib <- function( + .propensity, + .exposure, + exposure_type = c("auto", "binary", "categorical"), + .focal_level = NULL, + .reference_level = NULL, + ..., + .treated = NULL, + .untreated = NULL +) { + calculate_weight_from_modified_ps( + .propensity = .propensity, + .exposure = .exposure, + weight_fn = wt_entropy.numeric, + modification_type = "calib", + exposure_type = exposure_type, + .focal_level = .focal_level, + .reference_level = .reference_level, + .treated = .treated, + .untreated = .untreated, + ... + ) +} + +#' @export +wt_atc.ps_calib <- wt_atu.ps_calib + +#' @export +wt_cens.ps_calib <- function( + .propensity, + .exposure, + .sigma = NULL, + exposure_type = c("auto", "binary", "categorical", "continuous"), + .focal_level = NULL, + .reference_level = NULL, + stabilize = FALSE, + stabilization_score = NULL, + ..., + .treated = NULL, + .untreated = NULL +) { + calculate_weight_from_modified_ps( + .propensity = .propensity, + .exposure = .exposure, + weight_fn = wt_cens.numeric, + modification_type = "calib", + .sigma = .sigma, + exposure_type = exposure_type, + .focal_level = .focal_level, + .reference_level = .reference_level, + .treated = .treated, + .untreated = .untreated, + stabilize = stabilize, + stabilization_score = stabilization_score, + ... + ) +} diff --git a/man/ps_calibrate.Rd b/man/ps_calibrate.Rd index 82c8466..ec8d3b5 100644 --- a/man/ps_calibrate.Rd +++ b/man/ps_calibrate.Rd @@ -11,7 +11,6 @@ ps_calibrate( smooth = TRUE, .focal_level = NULL, .reference_level = NULL, - estimand = NULL, .treated = NULL, .untreated = NULL ) @@ -47,14 +46,12 @@ If not provided, \code{ps_calibrate()} will attempt to automatically determine t \item{.reference_level}{The value representing the reference group (typically control). If not provided, \code{ps_calibrate()} will attempt to automatically determine the coding.} -\item{estimand}{Character indicating the estimand type.} - \item{.treated}{\ifelse{html}{\href{https://lifecycle.r-lib.org/articles/stages.html#deprecated}{\figure{lifecycle-deprecated.svg}{options: alt='[Deprecated]'}}}{\strong{[Deprecated]}} Use \code{.focal_level} instead.} \item{.untreated}{\ifelse{html}{\href{https://lifecycle.r-lib.org/articles/stages.html#deprecated}{\figure{lifecycle-deprecated.svg}{options: alt='[Deprecated]'}}}{\strong{[Deprecated]}} Use \code{.reference_level} instead.} } \value{ -A calibrated propensity score object (\code{psw}) +A calibrated propensity score object (\code{ps_calib}) } \description{ This function calibrates propensity scores to improve their accuracy using diff --git a/tests/testthat/test-ps_calibration.R b/tests/testthat/test-ps_calibration.R index 20d7b75..687a86b 100644 --- a/tests/testthat/test-ps_calibration.R +++ b/tests/testthat/test-ps_calibration.R @@ -19,14 +19,14 @@ test_that("errors when ps and .exposure have different lengths", { ) }) -test_that("returns a psw object of correct length and range", { +test_that("returns a ps_calib object of correct length and range", { set.seed(42) ps <- rep(0.5, 100) treat <- rbinom(100, 1, 0.3) out <- ps_calibrate(ps, treat) - expect_s3_class(out, "psw") + expect_s3_class(out, "ps_calib") expect_length(out, 100) expect_true(all(as.numeric(out) >= 0 & as.numeric(out) <= 1)) }) @@ -40,23 +40,27 @@ test_that("constant ps yields calibrated = observed prevalence", { expect_equal(unique(as.numeric(out)), 0.5) }) -test_that("preserves psw attributes from an existing causal-weights object", { - # Create a dummy psw object - ps_orig <- psw( - x = runif(10), - estimand = "ATT", - stabilized = TRUE, - trimmed = TRUE, - truncated = FALSE - ) - treat <- rbinom(10, 1, ps_orig) +test_that("calibration metadata is properly stored", { + ps <- runif(10) + treat <- rbinom(10, 1, ps) + + # Test logistic calibration with smooth = TRUE + out_smooth <- ps_calibrate(ps, treat, method = "logistic", smooth = TRUE) + meta_smooth <- ps_calib_meta(out_smooth) + expect_equal(meta_smooth$method, "logistic") + expect_true(meta_smooth$smooth) - out <- ps_calibrate(ps_orig, treat) + # Test logistic calibration with smooth = FALSE + out_simple <- ps_calibrate(ps, treat, method = "logistic", smooth = FALSE) + meta_simple <- ps_calib_meta(out_simple) + expect_equal(meta_simple$method, "logistic") + expect_false(meta_simple$smooth) - expect_equal(attr(out, "estimand"), "ATT") - expect_true(attr(out, "stabilized")) - expect_true(attr(out, "trimmed")) - expect_false(attr(out, "truncated")) + # Test isotonic regression + out_iso <- ps_calibrate(ps, treat, method = "isoreg") + meta_iso <- ps_calib_meta(out_iso) + expect_equal(meta_iso$method, "isoreg") + expect_false(meta_iso$smooth) }) test_that("calibration changes the distribution", { @@ -151,7 +155,7 @@ test_that(".focal_level and .reference_level parameters work consistently with p # Test automatic detection with 0/1 coding treat_01 <- rbinom(30, 1, ps) calib_auto <- ps_calibrate(ps, treat_01) - expect_s3_class(calib_auto, "psw") + expect_s3_class(calib_auto, "ps_calib") # Test explicit specification with 0/1 coding calib_explicit <- ps_calibrate( @@ -209,7 +213,7 @@ test_that("automatic treatment detection works with binary vectors", { # All should work with automatic detection calib_01 <- ps_calibrate(ps, treat_01) - expect_s3_class(calib_01, "psw") + expect_s3_class(calib_01, "ps_calib") # These require explicit specification calib_12 <- ps_calibrate(ps, treat_12, .focal_level = 2, .reference_level = 1) @@ -278,7 +282,7 @@ test_that("handles NA values appropriately", { suppressWarnings(calibrated <- ps_calibrate(ps, treat)) expect_length(calibrated, 5) expect_true(is.na(calibrated[3])) - expect_s3_class(calibrated, "psw") + expect_s3_class(calibrated, "ps_calib") # Test with isotonic regression too calibrated_iso <- ps_calibrate(ps, treat, method = "isoreg") @@ -294,7 +298,7 @@ test_that("isotonic regression calibration works", { # Should not error calibrated_iso <- ps_calibrate(ps, treat, method = "isoreg") - expect_s3_class(calibrated_iso, "psw") + expect_s3_class(calibrated_iso, "ps_calib") expect_length(calibrated_iso, 100) expect_true(all( as.numeric(calibrated_iso) >= 0 & as.numeric(calibrated_iso) <= 1 @@ -422,8 +426,8 @@ test_that("smooth parameter works correctly for logistic calibration", { calib_simple <- ps_calibrate(ps, treat, smooth = FALSE) # Both should be psw objects - expect_s3_class(calib_smooth, "psw") - expect_s3_class(calib_simple, "psw") + expect_s3_class(calib_smooth, "ps_calib") + expect_s3_class(calib_simple, "ps_calib") # Both should be in valid range expect_true(all(