biomaRt多个物种同源基因一步转换 2022-08-30

适用背景

之前写了一篇跨物种同源基因转换的文章，里面用到biomaRt这个包，但是这个包只能对2个物种进行同源基因转换。如果有多个物种，使用起来就比较麻烦，因此本文将介绍一步转换多个物种同源基因的函数。

useEnsembl和useMart的区别

上一篇文章提到，用useMart时有时候会出现网络连接错误的问题，可以指定host解决，另外根据一篇Bioconductor文章，建议使用useEnsembl替代useMart，或许可以解决网络问题。
查看帮助文档发现，这两者除了名称没有太大差别，如果要深究需要看源代码了。因此，本文将使用useEnsembl而不用useMart，可以把下面脚本的useEnsembl都替换成useMart使用。

useEnsembl Connects to the selected BioMart database and dataset hosted by Ensembl

useMart Connects to the selected BioMart database and dataset

关于biomaRt未收录的非模式物种

因为有些非模式物种没有被收纳进biomaRt包，因此在找其它同源基因库时发现一个HGNC网站。

HGNC

这个网站有个同源基因转换工具，虽然物种不是很多，但列出了不少同源基因库的信息，感兴趣可以逐个了解一下。

HCOP

2个物种的同源基因转换函数get2sp_orth_gene

ref，参考物种的数据集名称，例如小鼠的是mmusculus_gene_ensembl
que，查询物种的数据集名称
gene.use，指定基因集，默认为空，使用查询物种的所有基因
a.ref，参考物种输出的属性类型向量，默认会输出c("ensembl_gene_id","external_gene_name")
f.ref，参考物种输入的属性类型，默认为基因名"external_gene_name"
a.que，查询物种输出的属性类型向量，默认输出基因名c("external_gene_name")
host，使用的数据库网站版本，默认为"https://dec2021.archive.ensembl.org/"

如果需要查询某个物种的数据集名称以及可用的属性类型，可以查看跨物种同源基因转换。

函数代码如下：

get2sp_orth_gene <- function(ref=NULL,que=NULL,gene.use=NULL,
                          a.ref=NULL,f.ref="external_gene_name",
                          a.que=NULL,uniqueRows=T,host = "https://dec2021.archive.ensembl.org/"){
a.ref <- unique(c("ensembl_gene_id","external_gene_name",a.ref))
a.que <- unique(c("external_gene_name",a.que))

mref <- useEnsembl("ensembl",dataset = ref, host = host)
mque <- useEnsembl("ensembl",dataset = que, host = host)

if (is.null(gene.use)) {
g1 <- getBM(attributes = c("external_gene_name"),mart = useEnsembl("ensembl",dataset = que, host = host),uniqueRows = T)
g2 <- g1$external_gene_name
}else{
g2 <- gene.use
}

gene.out <- getLDS(attributes = a.ref, filters = f.ref, values = g2, mart = mref,
                   attributesL = a.que, martL = mque, uniqueRows = uniqueRows)

nref <- length(a.ref)
nque <- length(a.que)

recol <- colnames(gene.out)
recol[1:2] <- c('ref.gene_ID','ref.gene_name')
recol[nref+1] <- 'que.gene_name'

colnames(gene.out) <- recol
gene.out$ref.specice <- ref
gene.out$que.specice <- que
return(gene.out)
}

使用示例
这个函数返回的是一个数据框，第1列ref.gene_ID是参数物种的基因ID，第2列ref.gene_name是参考物种的基因名，第3列que.gene_name是询问物种的基因名，第4列是ref.specice参考物种的数据集名称，第5列是参考物种的数据集名称que.specice。

> g1 <- get2sp_orth_gene(ref="mmusculus_gene_ensembl",que="mfascicularis_gene_ensembl")
> head(g1)
         ref.gene_ID ref.gene_name que.gene_name            ref.specice
1 ENSMUSG00000038473        Nos1ap        NOS1AP mmusculus_gene_ensembl
2 ENSMUSG00000036885      Arhgef26      ARHGEF26 mmusculus_gene_ensembl
3 ENSMUSG00000033285          Wdr3          WDR3 mmusculus_gene_ensembl
4 ENSMUSG00000074480         Mex3a         MEX3A mmusculus_gene_ensembl
5 ENSMUSG00000062127     Cttnbp2nl     CTTNBP2NL mmusculus_gene_ensembl
6 ENSMUSG00000041355          Ssr2          SSR2 mmusculus_gene_ensembl
                 que.specice
1 mfascicularis_gene_ensembl
2 mfascicularis_gene_ensembl
3 mfascicularis_gene_ensembl
4 mfascicularis_gene_ensembl
5 mfascicularis_gene_ensembl
6 mfascicularis_gene_ensembl

多个个物种的同源基因转换函数getmtsp_orth_gene

getmtsp_orth_gene调用了前面的get2sp_orth_gene函数，可以一步实现多个物种的同源基因转换，大部分参数与get2sp_orth_gene一样，只有que和gene.use不一样，que需要传入查询物种数据集的向量，例如que=c('mmusculus_gene_ensembl','mfascicularis_gene_ensembl','tguttata_gene_ensembl')，而gene.use如果不为空则需输入基因集列表，例如gene.use= list(c("Gad1","Sst"), c("GAD1","SST"))，基因集列表要一一对应查询物种，如果为空则会使用所有基因。