当前位置：首页 > news >正文

跨物种筛选同源基因

news 来源：原创 2025/5/14 11:43:27

工具：R：biomaRt

原始文件：human、mouse、macaque、marmoset四个物种的gene list，有些是用ensembl ID，有的是用gene name来表示。

目的：找到四个物种的gene list之间的1v1同源基因

1. 找到物种间的1v1同源基因

# 加载 biomaRt 包
library(biomaRt)
library(dplyr)
# 连接到 Ensembl 数据库
ensembl <- useMart("ensembl")# 使用read.csv来读取数据
hu <- read.csv(paste0(path,'hu_genelst.csv'), row.names = NULL)
mo <- read.csv(paste0(path,'mo_genelst.csv'), row.names = NULL)head(hu)
str(hu)

这里的hu是我需要筛选的human的基因，我先找到human的gene list对应的其它物种的同源基因。

datasets <- listDatasets(ensembl)
head(datasets)

这里有ensembl中包含的很多物种的dataset

可以用这个命令来筛选你想要的物种的dataset

datasets[grepl("musculus", datasets$dataset), ]

然后可以生成我需要的四个物种的 dataset

# 选择四个物种的数据集
ensembl_human <- useDataset("hsapiens_gene_ensembl", mart = ensembl)
ensembl_mouse <- useDataset("mmusculus_gene_ensembl", mart = ensembl)
ensembl_macaque <- useDataset("mmulatta_gene_ensembl", mart = ensembl)
ensembl_marmoset <- useDataset("cjacchus_gene_ensembl", mart = ensembl)
ensembl_human

ensembl_human这个dataset的attribute属性中会有很多相关的信息，比如human的ensembl_gene_id，以及基因名称等，甚至包括ensembl定义好的与其它物种的同源基因，如with_cjacchus_homolog，，这些信息很多可以根据关键词筛选：

filters = listFilters(ensembl_human) 
filters[grepl('ensembl_gene_id',filters$name),]
h_hom = filters[grepl('homolog',filters$name),]
h_hom[grepl('marmoset',h_hom$description),]
h_hom[grepl('macaque',h_hom$description),]

然后根据getBM函数获取想要的属性

hu_symbol = getBM(
attributes = c('ensembl_gene_id','hgnc_symbol','band'), ## 想要获取的属性
filters = 'ensembl_gene_id',  ## 根据哪个特征来获取这些属性
values = hu$ids[1:5],   ## 想要获取哪些基因的属性
mart = hu_ensembl  ## 从哪个dataset中获取
)
hu_symbol

2. 用getLDS函数进行同源基因 list的转换

也就是用’ensembl_gene_id‘表示的h_genelist，将他从ensembl的human dataset中提取出来，然后对应到mouse dataset的相应信息，提取出mouse dataset中对应的‘ensembl_gene_id’，实现从ENSG00000078808 到 ENSMUSG00000002459 的转换，也就是同源基因的匹配。。（但是！！但是！！我可能是网络连接错误，所以我无法进行这个转换！！！！！【只贴上代码】）

# 将鼠的基因转换为人的，人鼠同源基因-用人的进行表示
h_hvmo <- getLDS(attributes = c("ensembl_gene_id"),   ## 要检索的数据集的属性参数mart = mo_ensembl,    ## 需要查询的基因向量来自哪个mart对象filters = "ensembl_gene_id",   ## 使用哪个属性来查询values = mo$ids[1:3],   ## 需要查询的向量attributesL = c("ensembl_gene_id"),   ## 需要同源转换为人的这些属性的信息martL = hu_ensembl, uniqueRows = TRUE,verbose = TRUE)  # 启用详细的调试日志)# Error: biomaRt has encountered an unexpected server error.
# Consider trying one of the Ensembl mirrors (for more details look at ?useEnsembl)
# Traceback:
# 1. getLDS(attributes = c("ensembl_gene_id"), mart = mo_ensembl, filters = "ensembl_gene_id", values = mo$ids[1:3], attributesL = c("ensembl_gene_id"),martL = hu_ensembl, uniqueRows = TRUE, verbose = TRUE)
# 2. submitQueryXML(host = paste0(martHost(mart), sep), query = xmlQuery,httr_config = martHTTRConfig(mart))
# 3. stop(err_msg, call. = FALSE)

所以我直接尝试另一种方法，

3. 直接使用human dataset中的属性信息，

也就是网站提供的人的基因对应到其它物种上的信息，也就是这个human dataset中会有很多属性，包括ensembl_gene_id，以及从human基因出发与很多物种的同源基因信息，如mmusculus_homolog_ensembl_gene，可以直接提取出来

😊成功！！！

listAttributes(hu_ensembl)[grepl('cjacchus', listAttributes(hu_ensembl)$name), ]

homologs <- getBM(
attributes = c("ensembl_gene_id", "external_gene_name",
"cjacchus_homolog_ensembl_gene", "cjacchus_homolog_associated_gene_name",
'cjacchus_homolog_orthology_type',
"mmusculus_homolog_ensembl_gene",'mmusculus_homolog_associated_gene_name',
'mmusculus_homolog_orthology_type',
"mmulatta_homolog_ensembl_gene", "mmulatta_homolog_associated_gene_name",
'mmulatta_homolog_orthology_type'),
mart = hu_ensembl,
values = hu$ids[1:10],
filters = 'ensembl_gene_id'
)

然后去除为na或空格，以及orthology_type不为‘ortholog_one2one’的基因

# 检查是否存在空字符或 NA
sum(homologs$cjacchus_homolog_ensembl_gene == "")  # 空字符的数量
sum(is.na(homologs$cjacchus_homolog_ensembl_gene))  # NA的数量
sum(homologs$cjacchus_homolog_orthology_type == "ortholog_one2one") ## 1v1同源基因数目sum(homologs$mmusculus_homolog_ensembl_gene == "")  # 空字符的数量
sum(is.na(homologs$mmusculus_homolog_ensembl_gene))  # NA的数量
sum(homologs$mmusculus_homolog_orthology_type == "ortholog_one2one") ## 1v1同源基因数目sum(homologs$mmulatta_homolog_ensembl_gene == "")  # 空字符的数量
sum(is.na(homologs$mmulatta_homolog_ensembl_gene))  # NA的数量
sum(homologs$mmulatta_homolog_orthology_type == "ortholog_one2one") ## 1v1同源基因数目

# 过滤掉空字符和NA
filtered_homologs <- homologs[!(is.na(homologs$cjacchus_homolog_associated_gene_name) & homologs$cjacchus_homolog_associated_gene_name == "") &homologs$cjacchus_homolog_orthology_type == "ortholog_one2one" &!(is.na(homologs$mmusculus_homolog_associated_gene_name) & homologs$mmusculus_homolog_associated_gene_name == "") &homologs$mmusculus_homolog_orthology_type == "ortholog_one2one" &!(is.na(homologs$mmulatta_homolog_associated_gene_name) & homologs$mmulatta_homolog_associated_gene_name == "") &homologs$mmulatta_homolog_orthology_type == "ortholog_one2one", 
]# 查看过滤后的结果
filtered_homologs

🍎如果不筛选 ortholog_one2one 会存在一个人类基因对应某一个物种有多个同源基因的情况，如
13971 ENSG00000280789 PAGR1 ENSCJAG00000007182 PAGR1 ENSMUSG00000030680 Pagr1a ENSMMUG00000061749 PAGR1
13972 ENSG00000280789 PAGR1 ENSCJAG00000007182 PAGR1 ENSMUSG00000107068 Gm42742 ENSMMUG00000061749 PAGR1

🍎因此需要筛选一下1v1同源基因

# 保存为 CSV 文件
write.csv(filtered_homologs, file = paste0(out_dir,'orthologs_1to1_human_mouse_macaque_marmoset_9086.csv'), row.names = FALSE)

1. 找到物种间的1v1同源基因

2. 用getLDS函数进行同源基因 list的转换

3. 直接使用human dataset中的属性信息，

相关文章：