在RNA-seq或芯片數(shù)據(jù)下游分析中經常遇到需要將基因表達矩陣行名的ensembl_id ( gene_id ) 轉換為gene symbol(gene_name)的情況，而在轉換時經常會出現(xiàn)多個ensembl_id對應與一個gene symbol的情形，此時就出現(xiàn)了重復的gene symbol。
重復的gene symbol當然是不能作為基因表達矩陣行名的，此時就需要我們去除重復的gene symbol。

gene symbol去重復有一般有兩種思路：
1.只保留平均表達量最高的gene symbol
2.合并所有gene symbol（基因表達量進行加和或者取平均）

一、獲取ensembl_id與gene symbol的對應文件

1. 首先需要得到所需的gtf文件（最好是上游基因計數(shù)時所用文件），一般在gencode下載GENCODE - The Mouse GENCODE Release History，本次示例選取小鼠mm10（grcm38）基因組版本，因此選取GENCODE 對應版本M25，選擇regions為ALL的GTF文件進行下載即可
   image.png
   image.png
2. 接著需要提取gtf文件中ensembl_id(gene_id)和gene symbol(gene_name)的對應關系，此步在linux或者R中操作都可以，我個人比較喜歡用linux命令，因此示范一下在linux中的操作，最后會得到g2s_vm25_gencode.txt文件

gunzip  gencode.vM25.chr_patch_hapl_scaff.annotation.gtf.gz 
vim gtf_geneid2symbol_gencode.sh

##################以下為.sh文件內容
gtf="gencode.vM25.chr_patch_hapl_scaff.annotation.gtf"
### gene_id to gene_name
grep 'gene_id' $gtf | awk -F 'gene_id \"' '{print $2}' |awk -F '\"' '{print $1}' >gene_id_tmp
grep 'gene_id' $gtf | awk -F 'gene_name \"' '{print $2}' |awk -F '\"' '{print $1}' >gene_name_tmp
paste gene_id_tmp gene_name_tmp >last_tmp
uniq last_tmp >g2s_vm25_gencode.txt
rm *_tmp

bash gtf_geneid2symbol_gencode.sh

image.png

二、取最高表達量的一個重復gene symbol

以下代碼參考修改自ensembl_id轉換之兩種方法比較 - 簡書 (jianshu.com)
整體思路：
構建包含ensembl_id、gene symbol和基因表達中位數(shù)的ids對象——將gene symbol按照基因表達從大到小排列——去重復gene symbol行——根據(jù)ids的行名保留表達矩陣并將行名轉為gene symbol

###環(huán)境設置
library(tidyverse) # ggplot2 stringer dplyr tidyr readr purrr  tibble forcats
library(data.table) #多核讀取文件
head(counts)  #counts是需要轉換ensembl_id的表達矩陣

##從gtf文件提取信息,獲得gencode的基因id對應symbol的ids矩陣
ids <- data.frame(geneid=rownames(counts),
                    median=apply(counts,1,median)) #計算基因表達中位數(shù)，用于之后排序
g2s <- fread('g2s_vm25_gencode.txt',header = F,data.table = F) #載入從gencode的gtf文件中提取的信息文件
colnames(g2s) <- c("geneid","symbol")
  
table(ids$geneid %in% g2s$geneid) #查看需要轉化的geneid在g2s的匹配情況
ids <- ids[ids$geneid %in% g2s$geneid,] #取出在gencode數(shù)據(jù)庫的gtf注釋中能找到的geneid
ids$symbol <- g2s[match(ids$geneid,g2s$geneid),2] #match返回其第二個參數(shù)中第一個參數(shù)匹配的位置,把g2s的geneid按照ids$geneid的順序一個個取出來，從而得到ids$symbol這一列
ids <- ids[order(ids$symbol,ids$median,decreasing = T),] #把ids$symbol按照ids$median由大到小排序
 
##去重復
dim(ids); table(duplicated(ids$symbol)) #統(tǒng)計查看重復的symbol
#[1] 56262     3
#FALSE  TRUE 
#55492   770
ids <- ids[!duplicated(ids$symbol),]#取出不重復的ids$symbol
dim(ids); table(duplicated(ids$symbol))
#[1] 55492     3
#FALSE 
#55492 
  
##轉化geneid為symbol 
counts <- counts[rownames(ids),] #取出表達矩陣中ids有的行  
rownames(counts) <- ids[match(rownames(counts),ids$geneid),"symbol"]  #根據(jù)geneid和symbol進行匹配
head(counts)

三、合并所有重復的gene symbol(推薦)

主要思路為利用aggregate函數(shù)根據(jù)symbol列中的相同基因合并基因表達矩陣矩陣

library(tidyverse) # ggplot2 stringer dplyr tidyr readr purrr  tibble forcats
library(data.table) #多核讀取文件

g2s <- fread('g2s_vm25_gencode.txt',header = F,data.table = F) #載入從gencode的gtf文件中提取的信息文件
colnames(g2s) <- c("geneid","symbol")
  
symbol <- g2s[match(rownames(counts),g2s$geneid),"symbol"] #匹配counts行名對應的symbol
table(duplicated(symbol))
#FALSE  TRUE 
#55492   770 
  
##使用aggregate根據(jù)symbol列中的相同基因進行合并 
counts <- aggregate(counts, by=list(symbol), FUN=sum)
counts <- column_to_rownames(counts,'Group.1')
dim(counts)
#[1] 55492     4

id轉換前

id轉換后

基因名去重復的一點思考：
這兩種思路的差別在于,第一種只取表達量最高的基因，認為只有這個基因有意義，其余表達量靠后的相同基因不重要。第二種則是合并所有具有相同基因名的基因，考慮到了所有基因的表達情況，考慮更全面，因此個人更推薦。
實際分析中，由于我們一般差異分析是對不同樣品中的同一個基因進行的比較，因此這兩種方法實際差別并不大，按自己需求選擇即可。
（其實分析時有些人嫌麻煩，還會直接合并symbol和表達矩陣后隨緣保留一個重復基因，這種方法這就見仁見智了...）