ggplot2如何在已有基因覆盖度图下方绘制对应坐标外显子内含子

你现有代码已经完成了覆盖度平滑面积图的绘制，只需要额外准备基因结构注释数据，新增对应几何图层即可实现参考图效果，具体操作如下：

前置准备 #

先提取TEX15在目标区间（30689060~30748122）的结构注释：

• 从你使用的参考基因组对应GTF/GFF注释文件中，筛选出gene_name为TEX15、且落在绘图x轴区间内的所有外显子记录
• 整理得到外显子坐标表，至少包含start（外显子起始坐标）、end（外显子终止坐标）两列
• 内含子坐标不需要手动整理，可通过相邻外显子的间隔区间自动计算

方案1：原生ggplot2实现（无额外生信包依赖） #

核心思路是将基因结构固定在y轴底部的空白区域，通过线段绘制内含子、矩形绘制外显子，完全适配你现有代码逻辑。

1. 先整理注释数据，示例代码如下（替换成你提取的真实外显子坐标即可）：

# 外显子坐标表，替换为你从GTF提取的TEX15真实值
tex15_exon <- data.frame(
  start = c(30690002, 30694221, 30697105, 30701298, 30747892),
  end = c(30691024, 30695430, 30698221, 30702401, 30748098)
)
# 自动计算内含子区间
tex15_intron <- data.frame(
  start = tex15_exon$end[-nrow(tex15_exon)],
  end = tex15_exon$start[-1]
)
# 配置基因结构的y轴位置，根据你实际的覆盖度数值范围调整，避免和覆盖度区域重叠
gene_base_y <- -8
exon_height <- 4

2. 在你原有绘图代码基础上叠加基因结构图层，注意新增inherit.aes = F参数避免继承全局映射报错：

ggplot(z, aes(x=inicio, y=promedio, fill=Technology, group=Technology, color=Technology))+ 
  stat_smooth(
    geom = 'area', method = 'loess', span = 1/3,
    alpha = 1/2) + 
  # 绘制内含子横线
  geom_segment(
    data = tex15_intron,
    aes(x = start, xend = end, y = gene_base_y, yend = gene_base_y),
    color = "grey30", linewidth = 0.8, inherit.aes = F
  ) +
  # 可选：绘制内含子链方向箭头，正链箭头朝右、负链箭头朝左
  geom_segment(
    data = tex15_intron,
    aes(x = start, xend = end - (end-start)*0.1, y = gene_base_y, yend = gene_base_y),
    linewidth = 1, arrow = arrow(length = unit(0.1, "inches"), type = "closed"),
    color = "grey30", inherit.aes = F
  ) +
  # 绘制外显子矩形
  geom_rect(
    data = tex15_exon,
    aes(xmin = start, xmax = end, ymin = gene_base_y - exon_height/2, ymax = gene_base_y + exon_height/2),
    fill = "steelblue", color = "grey20", inherit.aes = F
  ) +
  scale_x_continuous(limits=c(30689060,30748122))+
  # 调整y轴范围，底部预留15%空间放置基因结构
  scale_y_continuous(expand = expansion(mult = c(0.15, 0.05))) +
  theme_bw()+ 
  theme(legend.text = element_text(size = 12))+ 
  xlab("Coordinates")+
  ylab("Depth")+
  ggtitle("TEX15")

方案2：用ggbio包快速绘制（适合多转录本场景） #

如果需要展示多个转录本的可变剪接结构，可以直接用生信专用绘图包ggbio自动解析GTF文件生成基因模型，不需要手动拆分内外显子：

library(ggbio)
library(GenomicRanges)
# 读入参考基因组注释文件
gtf <- import("your_reference.gtf")
# 筛选TEX15在目标区间的注释
target_region <- GRanges(seqnames = "TEX15所在染色体编号", IRanges(30689060, 30748122))
tex15_anno <- subsetByOverlaps(gtf, target_region)
tex15_anno <- tex15_anno[tex15_anno$gene_name == "TEX15"]

# 分别绘制覆盖度图和基因模型图，按轨道拼接
p_depth <- # 粘贴你原来的覆盖度绘图代码
p_gene_model <- ggplot() + geom_alignment(tex15_anno, fill = "steelblue", color = "grey20")
# 组合上下两个绘图轨道，高度比4:1
tracks(Depth = p_depth, `Gene Model` = p_gene_model, heights = c(4,1), xlim = c(30689060,30748122))

注意事项 #

• 所有新增的基因结构图层必须加inherit.aes = FALSE，否则会继承全局的fill/color映射，导致报错或颜色错乱
• 基因结构的y轴位置、y轴扩展比例需要根据你实际的测序深度范围调整，确保外显子/内含子不会和上方的覆盖度面积图重叠
• 提取注释时要保证GTF的基因组版本和测序比对用的参考基因组版本完全一致，否则会出现坐标错位

内容的提问来源于stack exchange，提问作者Criss D