01 小提琴图

# 带显著性标注的分面小提琴图
# 对应 Cell Figure S4-C 风格复现

required_packages <- c("ggplot2", "readxl", "ggh4x", "ggpubr")
missing_packages <- required_packages[!vapply(
  required_packages,
  requireNamespace,
  logical(1),
  quietly = TRUE
)]

if (length(missing_packages) > 0) {
  stop(
    "请先安装依赖包: ",
    paste(missing_packages, collapse = ", "),
    call. = FALSE
  )
}

library(ggplot2)
library(readxl)
library(ggh4x)
library(ggpubr)

input_file <- "data1.xlsx"
output_file <- "violin_significance.png"

required_columns <- c("Compartment", "CellType", "AT_content")

if (!file.exists(input_file)) {
  stop(
    "找不到输入文件 ", input_file,
    "，请将 data1.xlsx 放在脚本同级目录，或修改 input_file。",
    call. = FALSE
  )
}

data <- read_xlsx(input_file)
missing_columns <- setdiff(required_columns, names(data))
if (length(missing_columns) > 0) {
  stop(
    "输入文件缺少必要列: ",
    paste(missing_columns, collapse = ", "),
    call. = FALSE
  )
}

compartment_levels <- c("Active", "Repressive")
celltype_levels <- c("Exc", "Inh", "Ast", "Oli", "OPC", "Mic/Immune", "Vasc")

unknown_compartments <- setdiff(unique(na.omit(data$Compartment)), compartment_levels)
unknown_celltypes <- setdiff(unique(na.omit(data$CellType)), celltype_levels)
if (length(unknown_compartments) > 0) {
  stop("Compartment 含未定义取值: ", paste(unknown_compartments, collapse = ", "), call. = FALSE)
}
if (length(unknown_celltypes) > 0) {
  stop("CellType 含未定义取值: ", paste(unknown_celltypes, collapse = ", "), call. = FALSE)
}

if (!is.numeric(data$AT_content)) {
  stop("AT_content 必须是数值列。", call. = FALSE)
}

# 固定横轴和分面顺序，确保不同数据批次的图形结构一致。
data$Compartment <- factor(data$Compartment, levels = compartment_levels)
data$CellType <- factor(data$CellType, levels = celltype_levels)

data <- data[complete.cases(data[, required_columns]), ]
if (nrow(data) == 0) {
  stop("去除缺失值后没有可绘图的数据。", call. = FALSE)
}

fill_cols <- c("#F1BED4", "#143E96")
strip_cols <- c(
  "#2FA83B", "#A8D98B", "#E53935", "#FDB36A",
  "#9C4F1E", "#C9B0D9", "#F7C4D7"
)

plot <- ggplot(data, aes(x = Compartment, y = AT_content, fill = Compartment)) +
  geom_violin(width = 0.7, trim = FALSE, color = NA, alpha = 1) +
  stat_summary(
    fun.data = function(x) {
      data.frame(
        y = median(x),
        ymin = quantile(x, 0.25),
        ymax = quantile(x, 0.75)
      )
    },
    geom = "linerange",
    color = "white",
    linewidth = 0.6
  ) +
  stat_summary(
    fun = median,
    geom = "point",
    size = 3,
    color = "white",
    show.legend = FALSE
  ) +
  facet_wrap2(
    ~CellType,
    nrow = 1,
    scales = "free_x",
    strip = strip_nested(
      background_x = elem_list_rect(fill = alpha(strip_cols, alpha = 1))
    )
  ) +
  scale_fill_manual(values = fill_cols, drop = FALSE) +
  stat_compare_means(
    comparisons = list(compartment_levels),
    method = "wilcox.test",
    label = "p.format",
    label.y = 0.72,
    size = 3
  ) +
  labs(x = NULL, y = "AT content", fill = "Compartment") +
  theme_classic() +
  theme(
    strip.background = element_blank(),
    strip.text = element_text(color = "white", face = "bold"),
    panel.spacing = unit(0.6, "lines"),
    legend.position = "none"
  )

ggsave(
  filename = output_file,
  plot = plot,
  width = 12,
  height = 3.6,
  units = "in",
  dpi = 300,
  bg = "white"
)

message("已生成: ", normalizePath(output_file, mustWork = FALSE))

library(tidyverse)
install.packages("glossary")
devtools::install_github("debruine/faux")
devtools::install_github("psyteachr/introdataviz")
library(introdataviz) # Version :0.0.0.9.003
install.packages("rstatix")
library(rstatix)
library(ggpubr)
library(ggnewscale)


dd2 <- read_tsv("Split_Violin.tsv") %>% pivot_longer(- lncRNA_signature)

dd2$lncRNA_signature <- factor(dd2$lncRNA_signature,
                               levels = unique(dd2$lncRNA_signature))
dd2$name <- factor(dd2$name,levels = dd2$name %>% unique())

df_p_val1 <- dd2 %>% group_by(name) %>%
  wilcox_test(value  ~ lncRNA_signature) %>%
  adjust_pvalue(p.col = "p", method = "bonferroni") %>%
  add_significance(p.col = "p.adj") %>% 
  mutate(x=row_number(),y=12.5,
         type=symnum(p.adj,cutpoints = c(0,0.01,0.05,1),
                     symbols = c("Deepening","Shallowing","Non-significant trend"))) %>% 
  mutate(type=as.factor(type),
         xmin=x-0.2,xmax=x+0.2)
df_p_val1$type <- factor(df_p_val1$type,unique(df_p_val1$type))

ggplot(dd2,aes(x = name, y = value,fill= lncRNA_signature))+
  geom_split_violin(trim =F,color = NA,adjust = 1.5)+
  scale_fill_manual(values = c("#C6CDF7","#709AE1FF")) +
  new_scale_fill() +
  stat_summary(fun.data = "mean_sd",position=position_dodge(0.15),
               geom = "errorbar",width = 0.1,color="black") +
  stat_summary(fun = "mean", geom = "point",position=position_dodge(0.15),color="black") +
  geom_point(data=df_p_val1,aes(x=x,y=y,color=type,fill=type),
             inherit.aes = F,shape=19,size=4) +
  geom_segment(data=df_p_val1,aes(x=xmin,xend=xmax,y=12,yend=12),inherit.aes = F) +
  geom_segment(data=df_p_val1,aes(x=xmin,xend=xmin,y=11.7,yend=12),inherit.aes = F) +
  geom_segment(data=df_p_val1,aes(x=xmax,xend=xmax,y=11.7,yend=12),inherit.aes = F) +
  scale_color_manual(values = c("Shallowing"="#709AE1FF","Deepening"="#9986A5",
                                "Non-significant trend"="#F5CDB4")) +
  
  labs(x=NULL,y=NULL)+
  theme(axis.title = element_blank(),
        axis.text.x=element_text(angle =45,hjust =1,vjust =1,
                                 color="black",size = 9,margin = margin(b =2)),
        axis.text.y=element_text(color="black",size = 9,margin = margin(r =1)),
        panel.background = element_rect(fill = NA,color = NA),
        panel.grid.minor= element_line(linewidth=0.2,color="#e5e5e5"),
        panel.grid.major = element_line(linewidth =0.2,color="#e5e5e5"),
        panel.border = element_rect(fill=NA,color="black",size=1,linetype="solid"),
        legend.key=element_blank(),
        legend.title = element_blank(),
        legend.text = element_text(color="black",size=8),
        legend.box = "horizontal",
        legend.key.width=unit(0.5,'cm'),
        legend.key.height=unit(0.5,'cm'),
        legend.position = c(1,0), legend.justification = c(1,0),
        legend.background=element_blank())