2024-09-05 | ggplot绘制曼哈顿图----受够了有些包的自定义限制

整体脚本如下，需要你的文件包含CHR BP P三列，无所谓顺序，也无所谓有啥其他列，只要有这三列即可

library(data.table)
library(ggplot2)
library(dplyr)

dt = fread(infile)[, c("CHR", "BP", "P")]
# 1. 计算每个染色体的长度（BP最大值）
chr_lengths <- dt %>%
  group_by(CHR) %>%
  summarise(chr_len = max(BP))

# 2. 累计基因组位置
data <- dt %>%
  arrange(CHR, BP) %>%
  mutate(chr_cumsum = cumsum(c(0, head(chr_lengths$chr_len, -1)))[CHR]) %>%
  mutate(BPcum = BP + chr_cumsum)

# 3. 计算每个染色体在x轴上的中间点，用来显示染色体名称
axis_set <- data %>%
  group_by(CHR) %>%
  summarise(center = (max(BPcum) + min(BPcum)) / 2)

# 4. 绘制曼哈顿图
p = ggplot(data, aes(x=BPcum, y=-log10(P), color=factor(CHR))) +
  geom_point(alpha=0.9) +
  geom_hline(yintercept=-log10(5e-8), color="red", linetype="dashed") +
  scale_color_manual(values=rep(c("#1B2C62", "#4695BC"), 22)) + 
  scale_x_continuous(label=axis_set$CHR, 
                     breaks=axis_set$center, 
                     expand=c(0, 0)) +
  scale_y_continuous(expand=c(0, 0)) + 
  labs(x="Chromosome", y=expression(-log[10](P))) +
  theme_classic() +
  theme(legend.position="none",
        axis.text.x = element_text(face="bold", size=20, color="black"),
        axis.text.y = element_text(face="bold", size=20, color="black"),
        axis.title.x = element_text(face="bold", size=20, margin=margin(t=10), color="black"),
        axis.title.y = element_text(face="bold", size=20, margin=margin(t=10), color="black"),
        axis.line = element_line(color="black"), 
        axis.ticks = element_line(color="black"),
        axis.ticks.length = unit(0.3, "cm"))
ggsave("my.png", p, height = 8, width = 23, dpi=300)

结果如下

my.png