「ggplot2练习」画基因结构图

基因结构图从本质上就可以看成方块，直线，箭头的组合

解析GFF文件

bioconductor上有一个GenomicFeatures包，里面有一个makeTxDbFromGFF()函数可以解析GFF文件并构建TxDb对象，如何操作该对象见【Bioconductor系列】如何用Bioconductor对基因组注释.

以拟南芥TAIR10上的GFF为例(文件可以从TAIR10上下载),

library(GenomicFeatures)
GFF_file <- "C:/Users/DELL/Desktop/TAIR10_GFF3_genes.gff"
txdb <- makeTxDbFromGFF(GFF_file) 
genes_df <- as.data.frame(genes(txdb))
exons_df <- as.data.frame(exons(txdb))

genes()和exons()函数分别用于获取所有基因，所有外显子的GRanges对象，之后用as.data.frame转成ggplot2使用的数据框格式

提取目标区间的基因结构信息

ggplot2作图

第一步，调整主题。如下的代码就是将画布清空

theme_syntenty <- theme_classic() +
  theme(axis.line = element_blank(),
        axis.text = element_blank(),
        axis.ticks = element_blank(),
        axis.title = element_blank())
theme)set()

第二步，画外显子的结构

外显子就是一个一个方块，我们可以用有一定厚度的线段来表示，geom_segment()

p1 <- ggplot(exons_in_gene_df,
             aes(x=start, xend=end,y=0.5,yend=0.5)) + 
  geom_segment(size=4) 
p1

外显子1

也可以用geom_rect()画框框

p1 <- ggplot(exons_in_gene_df) + 
  geom_rect(aes(xmin=start, xmax=end,ymin=-0.1,ymax=0.1),
            colour="black", fill="white") +
  ylim(c(-1,1))
p1

外显子2

外显子之间一般都是有连线的，所以下一步就把这些线加上去

p1 + geom_segment(data=gene_df, aes(x=start,xend=end,y=0,yend=0))

外显子3

然而，这个连线居然是在外显子的上面，一点都不美观。一种解决方式，就是分别在外显子间加线段，另一种解决方法就是先画线，然后用长方形的图层覆盖线。

p1 <- ggplot(exons_in_gene_df) + 
  geom_segment(data=gene_df, aes(x=start,xend=end,y=0,yend=0)) +
  geom_rect(aes(xmin=start, xmax=end,ymin=-0.1,ymax=0.1),
            colour="black", fill="white") +
  ylim(c(-1,1))

p1 

更好看的外显子

最后我们在加上这个基因的名字，以及把填一个比较好看的颜色

p1 <- ggplot(exons_in_gene_df) + 
  geom_segment(data=gene_df, aes(x=start,xend=end,y=0,yend=0)) +
  geom_rect(aes(xmin=start, xmax=end,ymin=-0.1,ymax=0.1),
            colour="#282a73", fill="#282a73") +
  ylim(c(-1,1)) +
  geom_text(data=gene_df, aes(x=(start + end)/2, y = 0.2,label=gene_id))

最后结果