在什么情况下基因ID转换会100%失败？
在什么情况下基因ID转换会100%失败？
生信技能树数据挖掘班的2024年最后一期已经学习完一个多月了，群里有个学员遇到一个报错，他的基因ID在进行不同类型转换的时候居然100% 转换失败了！平时我们转换的时候也可能就10%以内会失败，下面来看看！报错如下：
他的数据截图如下：眼尖的同学肯定一眼就能看出来问题在哪，这个也在我们前面的帖子中提到过：驴的单细胞数据基因ID如何转换？，去这个帖子中看看是怎么回事吧！
答疑解惑
这个ID 的特征是 ENS + MUS小鼠物种缩写 + T转录本特征符号 + 11位唯一数字，很显然就是转录本ID，而不是基因ID。
使用包转换看看：
rm(list = ls())#清空当前的工作环境
options(scipen = 20)#不以科学计数法显示
library(data.table)
library(tinyarray)
library(clusterProfiler)
library(org.Mm.eg.db) #for Mouse
acc <- "GSE270697"
dir.create(acc)
# 下载raw文件夹
# https://ftp.ncbi.nlm.nih.gov/geo/series/GSE163nnn/GSE163558/suppl/GSE163558_RAW.tar
# 使用正则表达式替换：\w{3}$ 匹配每个单词最后的三个字符替换为空字符串，即去掉它们
s <- gsub("(\\w{3}$)", "", acc, perl = TRUE)
s
url <- paste0("https://ftp.ncbi.nlm.nih.gov/geo/series/",s,"nnn/",acc,"/suppl/",acc,"_RAW.tar")
url
file <- paste0(acc,"_RAW.tar")
file
# downloader::download(url, destfile = file)
files <- list.files(path = "GSE270697/", pattern = ".txt.gz")
files
data <- list()
for(i in 1:length(files)) {
# i <- 2
temp <- data.table::fread(paste0("GSE270697/",files[i]), data.table = F)
dim(temp)
head(temp)
data[[i]] <- temp
}
mat <- do.call(cbind, data)
identical(mat[,1], mat[,6])
mat <- mat[,c(1:3,seq(4,ncol(mat),by=5))]
colnames(mat)[4:10] <- stringr::str_split(files,"_",n=2,simplify =T)[,1]
head(mat)
keytypes(org.Mm.eg.db)
gene.trans <- bitr(mat$target_id, fromType = "ENSEMBLTRANS", toType = c("SYMBOL"), OrgDb = org.Mm.eg.db)
head(gene.trans)
mat_symbol <- merge(gene.trans, mat, by.x="ENSEMBLTRANS", by.y="target_id", all.y=T)
head(mat_symbol)
也有 80% 的转换失败~
> gene.trans <- bitr(mat$target_id, fromType = "ENSEMBLTRANS", toType = c("SYMBOL"), OrgDb = org.Mm.eg.db)
'select()' returned 1:many mapping between keys and columns
Warning message:
In bitr(mat$target_id, fromType = "ENSEMBLTRANS", toType = c("SYMBOL"), :
83.09% of input gene IDs are fail to map...
转换成功的：如 ENSMUST00000000001
这个成功的是转录本的来源基因symbol，而不是转录本本身的symbol。
转换失败的：如 ENSMUST00000000003
这个失败了，但是在数据库中也可以查得到。
而且也是有对应的基因的：
还是用gtf文件来转转看看
下载gtf：https://ftp.ensembl.org/pub/release-113/gtf/mus_musculus/Mus_musculus.GRCm39.113.gtf.gz
### gtf
library(tidyverse)
library(rtracklayer)
# 读取gtf文件
gtf <- rtracklayer::import('GSE270697/Mus_musculus.GRCm39.113.gtf.gz')
class(gtf)
str(gtf)
head(gtf)
# 转为数据框，简单探索
gtf_df <- as.data.frame(gtf)
colnames(gtf_df)
# 过滤transcript行
gtf_df <- filter(gtf_df, type=="transcript")
head(gtf_df)
colnames(gtf_df)
gtf_df_trans <- gtf_df[,c("gene_name","gene_biotype","transcript_id","transcript_name","transcript_biotype")]
head(gtf_df_trans)
# 合并
mat_symbol2 <- merge(gtf_df_trans, mat, by.x="transcript_id", by.y="target_id", all.y=T)
head(mat_symbol2)
可以看到转换出来的基本上都有symbol了，不管是基因symbol还是转录本symbol，后面做转录本层面的功能富集分析也好用：
关于est_counts是不是count值
问一下：Kallisto 定量出来的count带有小数吗？
Kallisto是一种用于转录组数据的快速、准确的转录本定量工具，它使用无比对的方法来估计转录本的丰度。根据搜索结果，Kallisto的定量结果输出文件中，
abundance.tsv文件包含了每个基因的表达量，其中est_counts列表示估计的counts，这个值通常是整数，表示映射到特定转录本的reads数量。不过，Kallisto在内部使用了一个更精细的模型来估计转录本的丰度，并且这些估计值可以是小数，以更准确地反映实际的生物学情况。然而，输出的abundance.tsv文件中的est_counts列通常是整数，因为它表示的是在给定条件下估计的reads数目，而不是实际的测量值。在某些情况下，如果需要更精细的定量结果，Kallisto也可以输出小数形式的丰度估计值，这通常在abundance.h5文件中，该文件是HDF5格式，可以包含更复杂的数据结构和更高精度的定量结果。
那这个数据还是自己下载fq去做一下上游的定量吧，顺便得到基因水平的定量~
腾讯云开发者
