Skip to content

Commit ed6eb14

Browse files
gaowgaow
authored
Merge pull request #13 from StatFunGen/dev
update package
2 parents 59727db + 1abfeca commit ed6eb14

14 files changed

+233
-149
lines changed

‎NAMESPACE‎

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -1,9 +1,9 @@
11
# Generated by roxygen2: do not edit by hand
22

33
export(colocboost)
4+
export(colocboost_assemble)
45
export(colocboost_check_update_jk)
56
export(colocboost_get_methods)
67
export(colocboost_inits)
7-
export(colocboost_posthoc)
88
export(colocboost_update)
99
export(colocboost_workhorse)

‎R/colocboost.R‎

Lines changed: 28 additions & 21 deletions
Original file line numberDiff line numberDiff line change
@@ -224,7 +224,14 @@ colocboost <- function(X = NULL, Y = NULL, # individual data
224224
warning("Error: X matrices do not have the same number of variables. Provide variable names to the colnames of X matrix.")
225225
return(NULL)
226226
}
227+
} else {
228+
# --- check if there is only one X, default variable names as X_1, ..., X_p
229+
X <- lapply(X, function(xx){
230+
if (is.null(colnames(xx))){ colnames(xx) <- paste0("X_", 1:ncol(xx)) }
231+
return(xx)
232+
})
227233
}
234+
228235
keep.variable.individual <- lapply(X, colnames)
229236
if (!is.list(X) & !is.list(Y)){
230237
warning("Error: Input X and Y must be the list containing genotype matrics and all phenotype vectors!")
@@ -553,27 +560,27 @@ colocboost <- function(X = NULL, Y = NULL, # individual data
553560
outcome_names = outcome_names)
554561

555562
# --- post-processing of the colocboost updates
556-
message("Starting post-hoc analyses and results summary.")
557-
cb_output <- colocboost_posthoc(cb_obj,
558-
coverage = coverage,
559-
func_intw = func_intw,
560-
alpha = alpha,
561-
check_null = check_null,
562-
check_null_method = check_null_method,
563-
check_null_max = check_null_max,
564-
dedup = dedup,
565-
overlap = overlap,
566-
n_purity = n_purity,
567-
min_abs_corr = min_abs_corr,
568-
coverage_singlew = coverage_singlew,
569-
median_abs_corr = median_abs_corr,
570-
between_cluster = between_cluster,
571-
between_purity = between_purity,
572-
weaker_ucos = weaker_ucos,
573-
merging = merging,
574-
tol = tol,
575-
output_level = output_level)
576-
563+
message("Starting assemble analyses and results summary.")
564+
cb_output <- colocboost_assemble(cb_obj,
565+
coverage = coverage,
566+
func_intw = func_intw,
567+
alpha = alpha,
568+
check_null = check_null,
569+
check_null_method = check_null_method,
570+
check_null_max = check_null_max,
571+
dedup = dedup,
572+
overlap = overlap,
573+
n_purity = n_purity,
574+
min_abs_corr = min_abs_corr,
575+
coverage_singlew = coverage_singlew,
576+
median_abs_corr = median_abs_corr,
577+
between_cluster = between_cluster,
578+
between_purity = between_purity,
579+
weaker_ucos = weaker_ucos,
580+
merging = merging,
581+
tol = tol,
582+
output_level = output_level)
583+
577584
return(cb_output)
578585
}
579586

R/colocboost_posthoc_utils.R renamed to R/colocboost_addhoc_utils.R

Lines changed: 8 additions & 7 deletions
Original file line numberDiff line numberDiff line change
@@ -74,9 +74,9 @@ merge_cos_ucos <- function(cb_obj, out_cos, out_ucos, coverage = 0.95,
7474
out_cos$cos$avWeight[[i]] <- out_cos$cos$avWeight[[i]][,pos_name]
7575
change_obj_coloc_i <- change_obj_coloc[i,]
7676
change_obj_each_j <- change_obj_each[j,]
77-
out_cos$cos$cs_change[i,] <- pmax(change_obj_coloc_i, change_obj_each_j)
78-
coloc_outcome <- c(coloc_outcome, fine_outcome)
79-
out_cos$cos$coloc_outcomes[[i]] <- coloc_outcome
77+
out_cos$cos$cs_change[i,] <- pmax(change_obj_coloc_i, change_obj_each_j)
78+
coloc_outcome <- sort(c(coloc_outcome, fine_outcome))
79+
out_cos$cos$coloc_outcomes[[i]] <- coloc_outcome
8080
}
8181

8282
}
@@ -168,21 +168,22 @@ merge_ucos <- function(cb_obj, past_out,
168168
is_merged <- c()
169169
for (i.m in 1:length(potential_merged)){
170170
temp_set <- as.numeric(potential_merged[[i.m]])
171-
is_merged <- c(is_merged, temp_set)
172-
# define merged set
173-
coloc_sets_merged <- c(coloc_sets_merged, list( unique(unlist(ucos_each[temp_set])) ))
174171
# refine avWeight
175172
merged <- out_ucos$avW_ucos_each[, temp_set]
176173
unique_coloc_outcomes <- as.numeric(gsub(".*Y(\\d+).*", "\\1", colnames(merged)))
174+
if (length(unique(unique_coloc_outcomes))==1) next
175+
# define merged set
176+
coloc_sets_merged <- c(coloc_sets_merged, list( unique(unlist(ucos_each[temp_set])) ))
177177
colnames(merged) <- paste0("outcome", unique_coloc_outcomes)
178178
coloc_outcomes_merged <- c(coloc_outcomes_merged,
179-
list(unique(unique_coloc_outcomes[order(unique_coloc_outcomes)])))
179+
list(unique(sort(unique_coloc_outcomes))))
180180
# colnames(temp) <- unique_coloc_outcomes
181181
avWeight_merged <- c(avWeight_merged, list(merged))
182182
# refine cs_change
183183
change_cs_tmp <- change_obj_each[temp_set, , drop = FALSE]
184184
cs_change_merged <- c(cs_change_merged,
185185
list(apply(change_cs_tmp, 2, max)))
186+
is_merged <- c(is_merged, temp_set)
186187
}
187188

188189
if (length(is_merged) != 0){

R/colocboost_posthoc.R renamed to R/colocboost_assemble.R

Lines changed: 60 additions & 57 deletions
Original file line numberDiff line numberDiff line change
@@ -1,11 +1,11 @@
1-
#' Main function for colocboost post-hoc analysis
1+
#' Main function for colocboost post aggregating analysis
22
#'
33
#' @details
44
#' The following functions are included in the post-hoc analysis:
55
#'
6-
#' Colocalization signal - `colocboost_posthoc_coloc` - identify the colocalized confidence sets and the corresponding causal configurations.
6+
#' Colocalization signal - `colocboost_assemble_cos` - identify the colocalized confidence sets and the corresponding causal configurations.
77
#'
8-
#' Un-colocalization signal - `colocboost_posthoc_noncoloc` - identify the causal confidence sets for each outcome only.
8+
#' Un-colocalization signal - `colocboost_assemble_ucos` - identify the causal confidence sets for each outcome only.
99
#'
1010
#' Add-hoc merging functions including
1111
#'
@@ -19,47 +19,49 @@
1919
#' Summary of the colocboost results and get the output of colocboost (TO-DO-LIST)
2020
#'
2121
#' @export
22-
colocboost_posthoc <- function(cb_obj,
23-
coverage = 0.95,
24-
func_intw = "fun_R",
25-
alpha = 1.5,
26-
check_null = 0.1,
27-
check_null_method = "profile",
28-
check_null_max=2e-5,
29-
dedup = TRUE,
30-
overlap = TRUE,
31-
n_purity = 100,
32-
min_abs_corr = 0.5,
33-
coverage_singlew = 0.8,
34-
median_abs_corr = NULL,
35-
between_cluster = 0.8,
36-
between_purity = 0.8,
37-
weaker_ucos = TRUE,
38-
merging = TRUE,
39-
tol = 1e-9,
40-
output_level = 1){
22+
colocboost_assemble <- function(cb_obj,
23+
coverage = 0.95,
24+
func_intw = "fun_R",
25+
alpha = 1.5,
26+
check_null = 0.1,
27+
check_null_method = "profile",
28+
check_null_max=2e-5,
29+
dedup = TRUE,
30+
overlap = TRUE,
31+
n_purity = 100,
32+
min_abs_corr = 0.5,
33+
coverage_singlew = 0.8,
34+
median_abs_corr = NULL,
35+
between_cluster = 0.8,
36+
between_purity = 0.8,
37+
weaker_ucos = TRUE,
38+
merging = TRUE,
39+
tol = 1e-9,
40+
output_level = 1){
4141

4242
if (class(cb_obj) != "colocboost"){
4343
stop("Input must from colocboost function!")}
4444

4545
# - data information
4646
data_info <- get_data_info(cb_obj)
47+
model_info <- get_model_info(cb_obj, outcome_names = data_info$outcome_info$outcome_names)
4748
if (data_info$n_outcomes == 1 & output_level == 1){ output_level = 2 }
4849
if (cb_obj$cb_model_para$num_updates == 1){
4950
cb_output <- list("cos_summary" = NULL,
5051
"vcp" = NULL,
5152
"cos_details" = NULL,
52-
"data_info" = data_info)
53+
"data_info" = data_info,
54+
"model_info" = model_info)
5355
# - save model and all coloc and single information for diagnostic
5456
if (output_level != 1){
5557
tmp <- get_full_output(cb_obj = cb_obj, past_out = NULL, variables = NULL)
5658
if (output_level == 2){
5759
cb_output$ucos_details = tmp$ucos_detials
58-
cb_output <- cb_output[c("cos_summary", "vcp", "cos_details", "data_info", "ucos_details")]
60+
cb_output <- cb_output[c("cos_summary", "vcp", "cos_details", "data_info", "model_info", "ucos_details")]
5961
} else {
6062
cb_output$ucos_details = tmp$ucos_detials
6163
cb_output$diagnostic_details = tmp[-1]
62-
cb_output <- cb_output[c("cos_summary", "vcp", "cos_details", "data_info", "ucos_details", "diagnostic_details")]
64+
cb_output <- cb_output[c("cos_summary", "vcp", "cos_details", "data_info", "model_info", "ucos_details", "diagnostic_details")]
6365
}
6466
if (data_info$n_outcome == 1){
6567
cb_output <- list("ucos_summary" = NULL, "pip" = NULL,
@@ -75,21 +77,21 @@ colocboost_posthoc <- function(cb_obj,
7577
}
7678
cb_obj <- get_max_profile(cb_obj, check_null_max=check_null_max, check_null_method = check_null_method)
7779
# --------- about colocalized confidence sets ---------------------------------
78-
out_cos <- colocboost_posthoc_cos(cb_obj,
79-
coverage = coverage,
80-
func_intw = func_intw,
81-
alpha = alpha,
82-
check_null = check_null,
83-
check_null_method = check_null_method,
84-
dedup = dedup,
85-
overlap = overlap,
86-
n_purity = n_purity,
87-
min_abs_corr = min_abs_corr,
88-
coverage_singlew = coverage_singlew,
89-
median_abs_corr = median_abs_corr,
90-
between_cluster = between_cluster,
91-
between_purity = between_purity,
92-
tol = tol)
80+
out_cos <- colocboost_assemble_cos(cb_obj,
81+
coverage = coverage,
82+
func_intw = func_intw,
83+
alpha = alpha,
84+
check_null = check_null,
85+
check_null_method = check_null_method,
86+
dedup = dedup,
87+
overlap = overlap,
88+
n_purity = n_purity,
89+
min_abs_corr = min_abs_corr,
90+
coverage_singlew = coverage_singlew,
91+
median_abs_corr = median_abs_corr,
92+
between_cluster = between_cluster,
93+
between_purity = between_purity,
94+
tol = tol)
9395

9496
# --------- about non-colocalized confidence sets ---------------------------------
9597
L <- cb_obj$cb_model_para$L
@@ -124,19 +126,19 @@ colocboost_posthoc <- function(cb_obj,
124126
}
125127
}
126128
class(cb_obj_single) <- "colocboost"
127-
out_ucos_each <- colocboost_posthoc_ucos(cb_obj_single,
128-
coverage = coverage,
129-
check_null = check_null,
130-
check_null_method = check_null_method,
131-
dedup = dedup,
132-
overlap = overlap,
133-
n_purity = n_purity,
134-
min_abs_corr = min_abs_corr,
135-
median_abs_corr = median_abs_corr,
136-
between_cluster = between_cluster,
137-
between_purity = between_purity,
138-
weaker_ucos = weaker_ucos,
139-
tol = tol)
129+
out_ucos_each <- colocboost_assemble_ucos(cb_obj_single,
130+
coverage = coverage,
131+
check_null = check_null,
132+
check_null_method = check_null_method,
133+
dedup = dedup,
134+
overlap = overlap,
135+
n_purity = n_purity,
136+
min_abs_corr = min_abs_corr,
137+
median_abs_corr = median_abs_corr,
138+
between_cluster = between_cluster,
139+
between_purity = between_purity,
140+
weaker_ucos = weaker_ucos,
141+
tol = tol)
140142
aaa <- out_ucos_each$ucos$ucos
141143
if (length(aaa) != 0){
142144
ucos_outcome <- c(ucos_outcome, rep(i, length(aaa)))
@@ -193,24 +195,25 @@ colocboost_posthoc <- function(cb_obj,
193195
cos_results <- get_cos_details(cb_obj, coloc_out = past_out$cos$cos, data_info = data_info)
194196
cb_output <- list("vcp" = cos_results$vcp,
195197
"cos_details" = cos_results$cos_results,
196-
"data_info" = data_info)
198+
"data_info" = data_info,
199+
"model_info" = model_info)
197200

198201
### - extract summary table
199202
target_idx <- cb_obj$cb_model_para$target_idx
200203
summary_table <- get_cos_summary(cb_output, target_outcome = data_info$outcome_info$outcome_names[target_idx])
201204
cb_output <- c(cb_output, list(cos_summary = summary_table))
202-
cb_output <- cb_output[c("cos_summary", "vcp", "cos_details", "data_info")]
205+
cb_output <- cb_output[c("cos_summary", "vcp", "cos_details", "data_info", "model_info")]
203206

204207
# - save model and all coloc and single information for diagnostic
205208
if (output_level != 1){
206209
tmp <- get_full_output(cb_obj = cb_obj, past_out = past_out, variables = data_info$variables, cb_output = cb_output)
207210
if (output_level == 2){
208211
cb_output <- c(cb_output, list("ucos_details" = tmp$ucos_details))
209-
cb_output <- cb_output[c("cos_summary", "vcp", "cos_details", "data_info", "ucos_details")]
212+
cb_output <- cb_output[c("cos_summary", "vcp", "cos_details", "data_info", "model_info", "ucos_details")]
210213
} else {
211214
cb_output <- c(cb_output, list("ucos_details" = tmp$ucos_details))
212215
cb_output$diagnostic_details = tmp[-1]
213-
cb_output <- cb_output[c("cos_summary", "vcp", "cos_details", "data_info", "ucos_details", "diagnostic_details")]
216+
cb_output <- cb_output[c("cos_summary", "vcp", "cos_details", "data_info", "model_info", "ucos_details", "diagnostic_details")]
214217
}
215218
# - if fine-boost, the summary table will be the summary of finemapping
216219
if (data_info$n_outcomes == 1){

R/colocboost_posthoc_cos.R renamed to R/colocboost_assemble_cos.R

Lines changed: 16 additions & 16 deletions
Original file line numberDiff line numberDiff line change
@@ -1,18 +1,18 @@
1-
colocboost_posthoc_cos <- function(cb_obj,
2-
coverage = 0.95,
3-
func_intw = "fun_R",
4-
alpha = 1.5,
5-
check_null = 0.1,
6-
check_null_method = "profile",
7-
dedup = TRUE,
8-
overlap = TRUE,
9-
n_purity = 100,
10-
min_abs_corr = 0.5,
11-
coverage_singlew = 0.8,
12-
median_abs_corr = NULL,
13-
between_cluster = 0.8,
14-
between_purity = 0.8,
15-
tol = 1e-9){
1+
colocboost_assemble_cos <- function(cb_obj,
2+
coverage = 0.95,
3+
func_intw = "fun_R",
4+
alpha = 1.5,
5+
check_null = 0.1,
6+
check_null_method = "profile",
7+
dedup = TRUE,
8+
overlap = TRUE,
9+
n_purity = 100,
10+
min_abs_corr = 0.5,
11+
coverage_singlew = 0.8,
12+
median_abs_corr = NULL,
13+
between_cluster = 0.8,
14+
between_purity = 0.8,
15+
tol = 1e-9){
1616

1717
if (class(cb_obj) != "colocboost"){
1818
stop("Input must from colocboost function!")}
@@ -186,7 +186,7 @@ colocboost_posthoc_cos <- function(cb_obj,
186186
weight_coloc <- do.call(cbind, av)
187187

188188
# Hierachical Clustering iteration based on sequenced weights
189-
cormat = cor(t(weight_coloc))
189+
cormat = get_cormat(t(weight_coloc))
190190
hc = hclust(as.dist(1-cormat))
191191
n_cluster = get_n_cluster(hc, cormat, between_cluster = between_cluster)$n_cluster
192192
index = cutree(hc,n_cluster)

0 commit comments

Comments
 (0)