如何用gggenes与ggplot2调整基因组基因绘图方向至一致

通过Prokka生成.gff文件，使用gggenomes::read_feats函数导入并筛选出需绘图的基因数据，得到如下数据框：

df
     file_id      seq_id  start    end strand        feat_id      locus_tag orientation     gene       cluster position
1  GENOME300 GENOME300_1 199743 201911      - IMEHDJCA_00189 IMEHDJCA_00189          -1    geneE cluster_00365   200827
2  GENOME300 GENOME300_1 201914 203275      - IMEHDJCA_00190 IMEHDJCA_00190          -1    geneD cluster_01313   202594
3  GENOME300 GENOME300_1 203272 205377      - IMEHDJCA_00191 IMEHDJCA_00191          -1    geneB cluster_00403   204324
4  GENOME300 GENOME300_1 205817 206176      + IMEHDJCA_00192 IMEHDJCA_00192           1 gene1033 cluster_06099   205996
5  GENOME300 GENOME300_1 206268 206663      + IMEHDJCA_00193 IMEHDJCA_00193           1    geneC cluster_05858   206465
6  GENOME300 GENOME300_1 206686 222306      + IMEHDJCA_00194 IMEHDJCA_00194           1    geneA cluster_00001   214496
7  GENOME320 GENOME320_8 123699 139310      - ILCJGNBA_01570 ILCJGNBA_01570          -1    geneA cluster_00001   131504
8  GENOME320 GENOME320_8 139333 139728      - ILCJGNBA_01571 ILCJGNBA_01571          -1    geneC cluster_05858   139530
9  GENOME320 GENOME320_8 139820 140179      - ILCJGNBA_01572 ILCJGNBA_01572          -1 gene1033 cluster_06099   139999
10 GENOME320 GENOME320_8 140619 142724      + ILCJGNBA_01573 ILCJGNBA_01573           1    geneB cluster_00403   141671
11 GENOME320 GENOME320_8 142721 144082      + ILCJGNBA_01574 ILCJGNBA_01574           1    geneD cluster_01322   143401
12 GENOME320 GENOME320_8 144085 146253      + ILCJGNBA_01575 ILCJGNBA_01575           1    geneE cluster_00365   145169

使用ggplot2结合gggenes的代码绘图时，GENOME320的基因方向与GENOME300相反，需重新计算GENOME320的基因位置，使二者的基因展示方向一致。绘图代码如下：

ggplot(df, aes(xmin = start, xmax = end, y = file_id, fill = gene)) +
    geom_gene_arrow(aes(forward = orientation), arrow_body_height = unit(7, "mm"), 
                    arrowhead_height = unit(7, "mm"), 
                    arrowhead_width = unit(5.7, "mm")) + 
    scale_fill_manual(values = c("orange", "blue", "red",  "green", "yellow", "#878787")) +
    facet_wrap(~ file_id, scales = "free", ncol = 1) +
    scale_y_discrete(position="right") +
    gggenes::theme_genes() +
    theme(panel.background = element_rect(fill = 'white' , color = 'white' ),
          panel.grid.major.y = ggplot2::element_line(colour = "grey", linewidth = 0.6),
          axis.title.y= element_blank(), 
          axis.text.y = element_text(size = 14, # sample name size
                                     family="Times New Roman", 
                                     face="bold"),
          plot.margin = unit(c(t=0.3, b=0.3, r=0.3, l=0.001), "mm") 
    )

要让两个基因组的基因展示方向一致，核心是对GENOME320的基因坐标进行反向映射，同时调整基因的方向标记。具体步骤如下：

1. 计算GENOME320所在序列的最大坐标值
2. 对GENOME320的每个基因，用总长度重新计算start和end：新start = 总长度 - 原end，新end = 总长度 - 原start
3. 反转GENOME320的orientation值（-1变1，1变-1），保证箭头方向正确

修改后的数据处理代码 #

library(dplyr)

# 计算GENOME320序列的最大坐标
genome320_max <- df %>% 
  filter(file_id == "GENOME320") %>% 
  pull(end) %>% 
  max()

# 调整GENOME320的基因坐标和方向
df_aligned <- df %>%
  mutate(
    start = case_when(
      file_id == "GENOME320" ~ genome320_max - end,
      TRUE ~ start
    ),
    end = case_when(
      file_id == "GENOME320" ~ genome320_max - start,
      TRUE ~ end
    ),
    orientation = case_when(
      file_id == "GENOME320" ~ -orientation,
      TRUE ~ orientation
    )
  ) %>%
  arrange(file_id, start) # 按基因组和起始坐标排序，保证基因顺序正确

修改后的绘图代码 #

使用调整后的df_aligned数据框绘图即可：

ggplot(df_aligned, aes(xmin = start, xmax = end, y = file_id, fill = gene)) +
    geom_gene_arrow(aes(forward = orientation), arrow_body_height = unit(7, "mm"), 
                    arrowhead_height = unit(7, "mm"), 
                    arrowhead_width = unit(5.7, "mm")) + 
    scale_fill_manual(values = c("orange", "blue", "red",  "green", "yellow", "#878787")) +
    facet_wrap(~ file_id, scales = "free", ncol = 1) +
    scale_y_discrete(position="right") +
    gggenes::theme_genes() +
    theme(panel.background = element_rect(fill = 'white' , color = 'white' ),
          panel.grid.major.y = ggplot2::element_line(colour = "grey", linewidth = 0.6),
          axis.title.y= element_blank(), 
          axis.text.y = element_text(size = 14, 
                                     family="Times New Roman", 
                                     face="bold"),
          plot.margin = unit(c(t=0.3, b=0.3, r=0.3, l=0.001), "mm") 
    )

原理说明 #

• 坐标反向映射相当于把GENOME320的序列“翻转”，让基因的排列顺序和GENOME300一致
• 反转orientation是因为坐标翻转后，基因的转录方向对应的箭头方向也需要反向，才能保证箭头指向正确的转录方向

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