Skip to content

Commit

Permalink
Warning
Browse files Browse the repository at this point in the history
  • Loading branch information
jinseob2kim committed May 7, 2024
1 parent 6b5fc73 commit 7e5aeaa
Showing 1 changed file with 35 additions and 13 deletions.
48 changes: 35 additions & 13 deletions R/forestcox.R
Original file line number Diff line number Diff line change
Expand Up @@ -8,6 +8,7 @@
#' @param decimal.hr Decimal for hazard ratio, Default: 2
#' @param decimal.percent Decimal for percent, Default: 1
#' @param decimal.pvalue Decimal for pvalue, Default: 3
#' @param cluster Cluster variable for coxph, Default: NULL
#' @return Sub-group analysis table.
#' @details This result is used to make forestplot.
#' @examples
Expand Down Expand Up @@ -50,7 +51,7 @@
#' @importFrom utils tail


TableSubgroupCox <- function(formula, var_subgroup = NULL, var_cov = NULL, data, time_eventrate = 3 * 365, decimal.hr = 2, decimal.percent = 1, decimal.pvalue = 3) {
TableSubgroupCox <- function(formula, var_subgroup = NULL, var_cov = NULL, data, time_eventrate = 3 * 365, decimal.hr = 2, decimal.percent = 1, decimal.pvalue = 3, cluster = NULL) {
. <- NULL

### var_subgroup이 categorical variable이 아닌 경우 중단 ###
Expand Down Expand Up @@ -113,7 +114,11 @@ TableSubgroupCox <- function(formula, var_subgroup = NULL, var_cov = NULL, data,
}
} else {
### survey data가 아닌 경우 ###
model <- survival::coxph(formula, data = data, x = TRUE)
if (is.null(cluster)){
model <- survival::coxph(formula, data = data, x = TRUE)
} else {
model <- survival::coxph(formula, data = data, x = TRUE, cluster = get(cluster))
}
# if (!is.null(model$xlevels) & length(model$xlevels[[1]]) != 2) stop("Categorical independent variable must have 2 levels.")

# KM 구하기(categorical인 경우)
Expand Down Expand Up @@ -236,18 +241,32 @@ TableSubgroupCox <- function(formula, var_subgroup = NULL, var_cov = NULL, data,
}
} else {
### survey data가 아닌 경우 ###
if (is.null(cluster)){
model <- data %>%
filter(!is.na(get(var_subgroup))) %>%
split(.[[var_subgroup]]) %>%
purrr::map(~ possible_coxph(formula, data = ., x = T))
} else{
model <- data %>%
filter(!is.na(get(var_subgroup))) %>%
split(.[[var_subgroup]]) %>%
purrr::map(~ tryCatch(coxph(formula, data = ., x = T, cluster = get(cluster)), error = function(e)NA))
}

data %>%
filter(!is.na(get(var_subgroup))) %>%
split(.[[var_subgroup]]) %>%
purrr::map(~ possible_coxph(formula, data = ., x = T)) -> model

data %>%
filter(!is.na(get(var_subgroup))) %>%
select(dplyr::all_of(var_subgroup)) %>%
table() %>%
names() -> label_val
xlev <- survival::coxph(formula, data = data)$xlevels
model.int <- possible_coxph(as.formula(gsub(xlabel, paste0(xlabel, "*", var_subgroup), deparse(formula))), data = data)

if (is.null(cluster)){
model.int <- possible_coxph(as.formula(gsub(xlabel, paste0(xlabel, "*", var_subgroup), deparse(formula))), data = data)
} else{
model.int <- tryCatch(coxph(as.formula(gsub(xlabel, paste0(xlabel, "*", var_subgroup), deparse(formula))), data = data, cluster = get(cluster)), error = function(e)NA)
}



# KM 구하기(categorical인 경우만)
Expand All @@ -273,7 +292,8 @@ TableSubgroupCox <- function(formula, var_subgroup = NULL, var_cov = NULL, data,

for (i in rownames(prop)) {
if (length(sub_xlev[[i]]) == 1) {
prop[i, paste0(xlabel, "=", sub_xlev[[i]])] <- res.kap.times[["0"]][["surv"]]
#prop[i, paste0(xlabel, "=", sub_xlev[[i]])] <- res.kap.times[["0"]][["surv"]]
prop[i, ] <- res.kap.times[["0"]][["surv"]]
} else if (length(sub_xlev[[i]]) > 1) {
surv.df <- data.frame(res.kap.times[[i]][c("strata", "surv")])
for (j in colnames(prop)) {
Expand Down Expand Up @@ -304,8 +324,9 @@ TableSubgroupCox <- function(formula, var_subgroup = NULL, var_cov = NULL, data,
} else {
model.int$call$formula <- as.formula(gsub(xlabel, paste0(xlabel, "*", var_subgroup), deparse(formula)))
model.int$call$data <- as.name("data")
pv_anova <- anova(model.int)
pv_int <- round(pv_anova[nrow(pv_anova), 4], decimal.pvalue)
pv_anova <- tryCatch(anova(model.int), error = function(e)NA)
if (is.na(pv_anova) & !is.null(cluster)) {warning("Warning: Anova test is not available for cluster data. So Interaction P value is NA when 3 and more categorical subgroup variable.")}
pv_int <- tryCatch(round(pv_anova[nrow(pv_anova), 4], decimal.pvalue), error = function(e)NA)
}
}

Expand Down Expand Up @@ -416,6 +437,7 @@ TableSubgroupCox <- function(formula, var_subgroup = NULL, var_cov = NULL, data,
#' @param decimal.percent Decimal for percent, Default: 1
#' @param decimal.pvalue Decimal for pvalue, Default: 3
#' @param line Include new-line between sub-group variables, Default: F
#' @param cluster Cluster variable for coxph, Default: NULL
#' @return Multiple sub-group analysis table.
#' @details This result is used to make forestplot.
#' @examples
Expand Down Expand Up @@ -449,17 +471,17 @@ TableSubgroupCox <- function(formula, var_subgroup = NULL, var_cov = NULL, data,
#' @importFrom magrittr %>%
#' @importFrom dplyr bind_rows

TableSubgroupMultiCox <- function(formula, var_subgroups = NULL, var_cov = NULL, data, time_eventrate = 3 * 365, decimal.hr = 2, decimal.percent = 1, decimal.pvalue = 3, line = F) {
TableSubgroupMultiCox <- function(formula, var_subgroups = NULL, var_cov = NULL, data, time_eventrate = 3 * 365, decimal.hr = 2, decimal.percent = 1, decimal.pvalue = 3, line = F, cluster = NULL) {
. <- NULL
xlabel <- setdiff(as.character(formula)[[3]], "+")[1]

out.all <- TableSubgroupCox(formula, var_subgroup = NULL, var_cov = var_cov, data = data, time_eventrate = time_eventrate, decimal.hr = decimal.hr, decimal.percent = decimal.percent, decimal.pvalue = decimal.pvalue)
out.all <- TableSubgroupCox(formula, var_subgroup = NULL, var_cov = var_cov, data = data, time_eventrate = time_eventrate, decimal.hr = decimal.hr, decimal.percent = decimal.percent, decimal.pvalue = decimal.pvalue, cluster = cluster)
out.all <- dplyr::mutate_all(out.all, as.character)

if (is.null(var_subgroups)) {
return(out.all)
} else {
out.list <- purrr::map(var_subgroups, ~ TableSubgroupCox(formula, var_subgroup = ., var_cov = var_cov, data = data, time_eventrate = time_eventrate, decimal.hr = decimal.hr, decimal.percent = decimal.percent, decimal.pvalue = decimal.pvalue))
out.list <- purrr::map(var_subgroups, ~ TableSubgroupCox(formula, var_subgroup = ., var_cov = var_cov, data = data, time_eventrate = time_eventrate, decimal.hr = decimal.hr, decimal.percent = decimal.percent, decimal.pvalue = decimal.pvalue, cluster = cluster))
if (line) {
out.newline <- out.list %>% purrr::map(~ rbind(NA, .))
result <- bind_rows(out.all, out.newline %>% dplyr::bind_rows() %>% dplyr::mutate_all(as.character))
Expand Down

0 comments on commit 7e5aeaa

Please sign in to comment.