Add initial files

DESCRIPTION

@@ -1,7 +1,7 @@
 Package: SPARKLE
 Type: Package
 Title: Single-cell Phenotype Association Research Kit for Large-scale Evaluation
-Version: 0.2.4
+Version: 0.3.6
 Authors@R: c(
  person("Chen", "Xi", ,email = "[email protected]", role = c("aut", "cre")),
  person("Ding", "Yang", ,email = "[email protected]", role = "ctb"),

R/f1.1_filter_cell.R

@@ -0,0 +1,22 @@
+filter_cell <- function(sparkle.data){
+
+ celltypeall <- unique(sparkle.data$Celltype)
+
+ filter_cellname <- list()
+ for(i in 1:length(celltypeall)){
+
+ df <- sparkle.data%>% dplyr::select_all() %>% dplyr::filter(Celltype==celltypeall[i])
+ num <- length(unique(df$Phenotype))
+
+ if(num<2){filter_cellname[[i]] <- celltypeall[i]}
+
+ }
+
+ filter_cell <- unlist(filter_cellname)
+
+ print(paste(filter_cell,"is filter due to the incomplete phenotype information"))
+
+ sparkle.data.new <- sparkle.data %>% dplyr::select_all() %>% dplyr::filter(!Celltype%in%filter_cell)
+
+ return(sparkle.data.new)
+}

R/f1.cwas_build_model_data_new.R

@@ -304,8 +304,17 @@ cwas_build_model_data <- function(inputdata,Sample="orig.ident",Phenotype,Cellty
  }
 
 
+ if (length(unique(celldf$Phenotype))==1) {
+
+ print("Warning:There is only one Phenotype in this data,please check your input !")
+
+ }
+
+
  celldf$Celltype <- gsub("-","_",celldf$Celltype)
 
+ celldf <- filter_cell(celldf)
+
  sf.rds(celldf, "01CWAS_data")
 
  return(celldf)

R/f12.2.Mediation_forest_plot.R → R/f12.2.Mediation_forest_plot_v2.R

@@ -10,7 +10,6 @@
 #'
 
 
-
 Mediation_forest_plot <- function(data){
 
 
@@ -57,22 +56,47 @@ Mediation_forest_plot <- function(data){
 
  data$Pvalue <- data$Indirect
  data$FDR <- p.adjust(data$Pvalue, method = "BH")
- show_variable2 <- c("Model","Effect"," ","Pvalue","FDR")
- data$" " <- paste(rep(" ", length(data[,1])), collapse = " ")
+ data$Mediation <-gsub("^.*Mediation:","",data$Model)
+ data$Celltype <-gsub("^.*X:","",data$Model)
+ data$Celltype <-gsub("Mediation.*","", data$Celltype)
+
+ numeric_cols <- sapply(data, is.numeric)
+ data[numeric_cols] <- lapply(data[numeric_cols], function(x) round(x, 2))
+ data$OR<- c(data$CI_lower+data$CI_upper)
+ data$OR<- data$OR/2
+
+ data$sig <- ifelse(data$Pvalue<0.001,"***",
+ ifelse(data$Pvalue<0.01,"**",
+ ifelse(data$Pvalue<0.05,"*","")))
+
+ data$Pvalue <- ifelse(data$Pvalue<0.01,"<0.01",data$Pvalue )
+ data$"P-value" <- paste0(data$Pvalue,data$sig)
+
+ data1 <- data
+ numeric_cols <- sapply(data1, is.numeric)
+ data1[numeric_cols] <- lapply(data[numeric_cols], function(x) format(x, scientific = TRUE, digits = 2))
+
+ data$"OR(95%CI)" <- paste0( data1$OR," ","[",data1$CI_lower,",",data1$CI_upper,"]")
+ data$OR <- ifelse(data$OR > 20, 10, data$OR)
+ data$"OR(95%CI)"[grepl("Inf", data1$OR)] <- "Not reliable"
+
+ show_variable2 <- c("Celltype","Mediation","Effect"," ","P-value")
+
+ data$" " <- paste(rep(" ", length(data[,1])), collapse = " ")
  pp <- forestploter::forest(data[,show_variable2], # 选择要在森林图中使用的数据列，这里包括变量名列、患者数量列、绘图要用的空白列和HR（95%CI）列
  est = data$Effect, # 效应值，也就是HR列
  lower = data$CI_lower, # 置信区间下限
  upper = data$CI_upper, # 置信区间上限
  sizes = 0.8, # 黑框框的大小
- ci_column = 3, # 在第3列（可信区间列）绘制森林图
+ ci_column = 4, # 在第3列（可信区间列）绘制森林图
  ref_line = 0, # 添加参考线
  arrow_lab = c("Down regulate", "Up regulate"), # 箭头标签，用来表示效应方向，如何设置取决于你的样本情况
  #xlim = c(-1, 1), # 设置x轴范围
  #ticks_at = c( -1,0,1 ), # 在指定位置添加刻度
  theme = tm, # 添加自定义主题
  footnote = "Red line:Adjusted.Pvaule<0.05") # 添加脚注信息
 
- pp <- forestploter::edit_plot(pp, row = which(data$FDR<0.0500001), gp = grid::gpar(col = "red4", fontface = "italic"))
+ pp <- forestploter::edit_plot(pp, row = which(data$Pvalue<0.0500001), gp = grid::gpar(col = "red4", fontface = "italic"))
  pp <- forestploter::insert_text(pp,
  text ="CWAS Mediation Analysis",
  col = 1:5,
@@ -81,3 +105,5 @@ Mediation_forest_plot <- function(data){
  pp <- forestploter::add_border(pp, part = "header", where = "bottom")
  return(pp)
 }
+
+

R/f5.cwas_allmodel_cal.R

@@ -73,7 +73,7 @@ if(method=="mix_effect"){
  print("Please select method for [mix_effect],[fix_effect],or[all]")
 
 }
-
+ mdlist[[i]]$inputdata <- cwas.data
 
 }
 

R/f8.4.OR_plot.R

@@ -21,142 +21,19 @@
 
 OR_plot <- function(data,show_variable=c("Celltype","Model","OR(95%CI)","P-value","P-value adjusted"),plot.position=4,Tablename="Tablename"){
 
- cwas.data <- data[["All_models"]][[1]][["Rawdata"]]
- atri_info <- attributes(cwas.data)
- data <- data[[2]]
- show_variable2 <- show_variable[1:plot.position-1]
- show_variable2[plot.position] <- c(" ")
- num <- length(show_variable)
- num1 <- plot.position+1
- num2 <- length(show_variable)+1
- show_variable2[num1:num2] <- show_variable[plot.position:num]
+ raw.sparkle.data <- data[["All_models"]][[1]][["inputdata"]]
 
- # 导入所需的包
- # library(metafor)
- #library(grid)
- #library(forestploter)
+ if(is.null(raw.sparkle.data)){
 
- # 创建数据框，包括研究名称、OR值和95%置信区间的上下限，以及事件数量和总样本量
+  pp <- OR_plot_origin(data = data,show_variable =show_variable,plot.position =plot.position,Tablename =Tablename )
 
- numeric_cols <- sapply(data, is.numeric)
+}else{
 
- data1 <- data
- data[numeric_cols] <- lapply(data[numeric_cols], function(x) round(x, 2))
- data1[numeric_cols] <- lapply(data[numeric_cols], function(x) format(x, scientific = TRUE, digits = 2))
-
- data$" " <- paste(rep(" ", length(data[,1])), collapse = " ")
- data$" " <- paste(rep(" ", length(data[,1])), collapse = " ")
-
-
- #data$CI_upper_table <- ifelse(data$CI_upper > 10, ">10", data$CI_upper)
-
- data$"OR(95%CI)" <- paste0(data1$OR," ","[",data1$CI_lower,",",data1$CI_upper,"]")
- data$OR <- ifelse(data$OR > 20, 10, data$OR)
- data$"OR(95%CI)"[grepl("Inf", data1$OR)] <- "Not reliable"
- # 对所有class是numeric的列保留小数点后两位数字
- data$CI_upper <- ifelse(data$CI_upper > 10, 10, data$CI_upper)
-
- data$"P-value" <-ifelse(data$Pvalue < 0.01, "<0.01", data$Pvalue)
- data$"P-value" <- paste(data$"P-value",data$significance)
- data$"P-value" <- paste0(" ",data$"P-value")
- data$"P-value" [grepl("Inf", data1$OR)] <- "Not reliable"
-
- tryCatch({
- data$"P-value adjusted" <-ifelse(data$adjustedPvalue < 0.01, "<0.01", data$adjustedPvalue)
- data$"P-value adjusted" <- paste(data$"P-value adjusted",data$significance.adjusted)
- data$"P-value adjusted" <- paste0(" ",data$"P-value adjusted")
- data$"P-value adjusted"[grepl("Inf", data1$OR)] <- "Not reliable"
- }, error = function(e) {
- # 捕获到错误时执行的代码
- # 打印出错信息
- print(paste("Warning Message:","No Adjusted p value in this data"))
-
- show_variable2<<- show_variable2[show_variable2 != "P-value adjusted"]
-
- })
-
- data$Pvalue[grepl("Inf", data1$OR)] <- "Not reliable"
-
- if(is.null(data$formula)){ data$formula <- c("Integrated model") }
-
- data$Model <- data$formula
-
- data$Celltype <- data$celltype
- data$Significance <- data$significance
-
-
- tm <- forestploter::forest_theme(
- base_size = 10, # 设置文本的基础大小
-
- # 设置可信区间的外观
- ci_pch = 15, # 可信区间点的形状
- ci_col = "blue4", # 可信区间的边框颜色
- ci_fill = "blue4", # 可信区间的填充颜色
- ci_alpha = 0.8, # 可信区间的透明度
- ci_lty = 1, # 可信区间的线型
- ci_lwd = 1.5, # 可信区间的线宽
- ci_Theight = 0.2, # 设置T字在可信区间末端的高度，默认是NULL
-
- # 设置参考线的外观
- refline_lwd = 2, # 参考线的线宽
- refline_lty = "dashed", # 参考线的线型
- refline_col = "red4", # 参考线的颜色
-
- # 设置垂直线的外观
- vertline_lwd = 1, # 垂直线的线宽，可以添加一条额外的垂直线，如果没有就不显示
- vertline_lty = "dashed", # 垂直线的线型
- vertline_col = "grey20", # 垂直线的颜色
-
- # 设置脚注的字体大小、字体样式和颜色
- footnote_cex = 0.6, # 脚注字体大小
- footnote_fontface = "italic", # 脚注字体样式
- footnote_col = "red4" # 脚注文本的颜色
- )
-
- # 绘制森林图
-
-
- footnote1 <- c("* p<0.05 \ ** p<0.01 \ *** p<0.001 ")
-
- info_parts <- list(
- paste("Phenotype:", ifelse(is.null(atri_info$Phenotype), "Not Provided", atri_info$Phenotype)," ",
- "Control_label:", ifelse(is.null(atri_info$Control_label), "Not Provided", atri_info$Control_label)," ",
- "Disease_label:", ifelse(is.null(atri_info$Disease_label), "Not Provided", atri_info$Disease_label)),
- paste("Group:", ifelse(is.null(atri_info$Group), "Not Provided", atri_info$Group)," ",
- "Subgroup:", ifelse(is.null(atri_info$Subgroup), "Not Provided", atri_info$Subgroup)," ",
- "Covariate 1:", ifelse(is.null(atri_info$Covariate1), "Not Provided", atri_info$Covariate1)," ",
- "Covariate 2:", ifelse(is.null(atri_info$Covariate2), "Not Provided", atri_info$Covariate2) )
- )
-
-
- info <- paste(info_parts, collapse = "\n")
- footnote2 <- paste0(footnote1, info)
-
-
-
- pp <- forestploter::forest(data[,show_variable2], # 选择要在森林图中使用的数据列，这里包括变量名列、患者数量列、绘图要用的空白列和HR（95%CI）列
- est = data$OR, # 效应值，也就是HR列
- lower = data$CI_lower, # 置信区间下限
- upper = data$CI_upper, # 置信区间上限
- sizes = 0.8, # 黑框框的大小
- ci_column = 4, # 在第3列（可信区间列）绘制森林图
- ref_line = 1, # 添加参考线
- arrow_lab = c("Low risk", "High Risk"), # 箭头标签，用来表示效应方向，如何设置取决于你的样本情况
- xlim = c(-1, 10), # 设置x轴范围
- ticks_at = c(0,1, 3, 5), # 在指定位置添加刻度
- theme = tm, # 添加自定义主题
- footnote = footnote2) # 添加脚注信息
-
- pp <- forestploter::edit_plot(pp, row = which(data$Pvalue<0.0500001), gp = grid::gpar(col = "red4", fontface = "italic"))
-
- pp <- forestploter::insert_text(pp,
- text =Tablename,
- col = 1:5,
- part = "header",
- gp = grid::gpar(fontface = "bold"))
- pp <- forestploter::add_border(pp, part = "header", where = "bottom")
- print(pp)
+ pp <- OR_plot_compare(data = data,show_variable =c("Celltype", "Model", "OR(95%CI)",
+ "SPARKLE (p-value)", "SPARKLE (p-value) adjusted","Wilicoxon (p-value)"),
+ plot.position =plot.position,Tablename =Tablename )
 
+ }
  return(pp)
 
 }

R/f8.4.OR_plot_covariant.R

@@ -0,0 +1,85 @@
+or_model <- function (data, show_variable = c("Celltype", "Covariant", "OR(95%CI)",
+ "P-value", "P-value adjusted"), plot.position = 4, Tablename = "Tablename")
+{
+
+
+ cwas.data <- data
+ atri_info <- attributes(cwas.data)
+
+ show_variable2 <- show_variable[1:plot.position - 1]
+ show_variable2[plot.position] <- c(" ")
+ num <- length(show_variable)
+ num1 <- plot.position + 1
+ num2 <- length(show_variable) + 1
+ show_variable2[num1:num2] <- show_variable[plot.position:num]
+ numeric_cols <- sapply(data, is.numeric)
+ data1 <- data
+ data[numeric_cols] <- lapply(data[numeric_cols], function(x) round(x,
+ 2))
+ data1[numeric_cols] <- lapply(data[numeric_cols], function(x) format(x,
+ scientific = TRUE, digits = 2))
+ data$" " <- paste(rep(" ", length(data[, 1])), collapse = " ")
+ data$" " <- paste(rep(" ", length(data[, 1])), collapse = " ")
+ data$"OR(95%CI)" <- paste0(data1$OR, " ", "[", data1$CI_lower,
+ ",", data1$CI_upper, "]")
+ data$OR <- ifelse(data$OR > 20, 10, data$OR)
+ data$"OR(95%CI)"[grepl("Inf", data1$OR)] <- "Not reliable"
+ data$CI_upper <- ifelse(data$CI_upper > 10, 10, data$CI_upper)
+ data$"P-value" <- ifelse(data$Pvalue < 0.01, "<0.01", data$Pvalue)
+ data$"P-value" <- paste(data$"P-value", data$significance)
+ data$"P-value" <- paste0(" ", data$"P-value")
+ data$"P-value"[grepl("Inf", data1$OR)] <- "Not reliable"
+ tryCatch({
+ data$"P-value adjusted" <- ifelse(data$adjustedPvalue <
+ 0.01, "<0.01", data$adjustedPvalue)
+ data$"P-value adjusted" <- paste(data$"P-value adjusted",
+ data$significance.adjusted)
+ data$"P-value adjusted" <- paste0(" ", data$"P-value adjusted")
+ data$"P-value adjusted"[grepl("Inf", data1$OR)] <- "Not reliable"
+ }, error = function(e) {
+ print(paste("Warning Message:", "No Adjusted p value in this data"))
+ show_variable2 <<- show_variable2[show_variable2 !=
+ "P-value adjusted"]
+ })
+ data$Pvalue[grepl("Inf", data1$OR)] <- "Not reliable"
+ if (is.null(data$formula)) {
+ data$formula <- c("Integrated model")
+ }
+ data$Model <- data$formula
+ data$Celltype <- data$celltype
+ data$Significance <- data$significance
+ tm <- forestploter::forest_theme(base_size = 10, ci_pch = 15,
+ ci_col = "blue4", ci_fill = "blue4", ci_alpha = 0.8,
+ ci_lty = 1, ci_lwd = 1.5, ci_Theight = 0.2, refline_lwd = 2,
+ refline_lty = "dashed", refline_col = "red4", vertline_lwd = 1,
+ vertline_lty = "dashed", vertline_col = "grey20", footnote_cex = 0.6,
+ footnote_fontface = "italic", footnote_col = "red4")
+ footnote1 <- c("* p<0.05 ** p<0.01 *** p<0.001 ")
+ info_parts <- list(paste("Phenotype:", ifelse(is.null(atri_info$Phenotype),
+ "Not Provided", atri_info$Phenotype), " ", "Control_label:",
+ ifelse(is.null(atri_info$Control_label), "Not Provided",
+ atri_info$Control_label), " ", "Disease_label:",
+ ifelse(is.null(atri_info$Disease_label), "Not Provided",
+ atri_info$Disease_label)), paste("Group:", ifelse(is.null(atri_info$Group),
+ "Not Provided", atri_info$Group), " ", "Subgroup:",
+ ifelse(is.null(atri_info$Subgroup), "Not Provided",
+ atri_info$Subgroup), " ", "Covariate 1:", ifelse(is.null(atri_info$Covariate1),
+ "Not Provided", atri_info$Covariate1), " ", "Covariate 2:",
+ ifelse(is.null(atri_info$Covariate2), "Not Provided",
+ atri_info$Covariate2)))
+ info <- paste(info_parts, collapse = "\n")
+ #footnote2 <- paste0(footnote1, info)
+ footnote2 <- footnote1
+ pp <- forestploter::forest(data[, show_variable2], est = data$OR,
+ lower = data$CI_lower, upper = data$CI_upper, sizes = 0.8,
+ ci_column = 4, ref_line = 1, arrow_lab = c("Low risk",
+ "High Risk"), xlim = c(-1, 10), ticks_at = c(0,
+ 1, 3, 5), theme = tm, footnote = footnote2)
+ pp <- forestploter::edit_plot(pp, row = which(data$Pvalue <
+ 0.0500001), gp = grid::gpar(col = "red4", fontface = "italic"))
+ pp <- forestploter::insert_text(pp, text = Tablename, col = 1:5,
+ part = "header", gp = grid::gpar(fontface = "bold"))
+ pp <- forestploter::add_border(pp, part = "header", where = "bottom")
+ print(pp)
+ return(pp)
+}