最近学员提问了一个有意思的数据集,他使用我们授课的示例代码发现跑降维聚类分群是没有问题,但是在后面的特异性基因的可视化的时候就全军覆没了。
我让学员发来一下对应的gse数据集,然后去下载这个这个文件,自己读取看了看;
counts <- data.table::fread('GSE190482_UMIsMatrix.txt.gz',data.table = F)
counts[1:4,1:4]
> counts[1:4,1:4]
V1 AAACCCAAGGTCGAGT-1_KO_B0310 AAACGAAAGGACACTG-1_WT_B0302
1 ENSMUSG00000051951 0 0
2 ENSMUSG00000089699 0 0
3 ENSMUSG00000102331 0 0
4 ENSMUSG00000102343 0 0
发现是因为里面的基因名字问题, 目前人类基因的命名体系主要包括 Ensembl、Gene Symbol 和 RefSeq。这些体系都是用于标识和命名基因的,是可以互相转换的 :
既然是需要转换,所以我使用AnnoProbe包进行这个步骤,然后我测试了一下,转换的成功率还很高的:
> ids=annoGene(ensID,'ENSEMBL','mouse')
Warning message:
In annoGene(ensID, "ENSEMBL", "mouse") :
0.27% of input IDs are fail to annotate...
> head(ids)
SYMBOL biotypes ENSEMBL chr start end
3 Xkr4 protein_coding ENSMUSG00000051951 chr1 3205901 3671498
8 Gm1992 antisense ENSMUSG00000089699 chr1 3466587 3513553
12 Gm19938 sense_intronic ENSMUSG00000102331 chr1 3647309 3658904
16 Gm37381 lincRNA ENSMUSG00000102343 chr1 3905739 3986215
17 Rp1 protein_coding ENSMUSG00000025900 chr1 3999557 4409241
20 Sox17 protein_coding ENSMUSG00000025902 chr1 4490931 4497354
最后的完整的代码如下所示:
counts <- data.table::fread('GSE190482_UMIsMatrix.txt.gz',data.table = F)
counts[1:4,1:4]
ensID = counts$V1
library(AnnoProbe)
ids=annoGene(ensID,'ENSEMBL','mouse')
head(ids)
tail(sort(table(ids$biotypes)))
ids=ids[ids$biotypes=='protein_coding',]
ids=ids[!duplicated(ids$SYMBOL),]
ids=ids[!duplicated(ids$ENSEMBL),]
pos = match(ids$ENSEMBL,counts$V1)
counts=counts[pos,]
rownames(counts)=ids$SYMBOL
counts=counts[,-1]
sce.all <- CreateSeuratObject(counts = counts,
# project = "Humanized vasculitis mice",
min.cells = 3)
sce.all
有了这个seurat的对象,后面就是我们常规的!代码在:(链接: https://pan.baidu.com/s/1pKEnPmWXi-pTab0WZUWzgg?pwd=a7s1) 相信大家很容易跟着去复现一次!
值得一提的是这个GSE190482数据集关联的文章是Kaymak T, et al. Gut 2023;标题是 :《IL-20 subfamily cytokines impair the oesophageal epithelial barrier by diminishing filaggrin in eosinophilic oesophagitis》,但是文章正文里面根本就没有出现过一个单细胞图表,仅仅是在附件出现了一个Supplemental Figure 4. 的子图:(J) Normalised expression levels of Il20rα and Il20rβ in annotated cell types from the scRNA-seq dataset.
看了两个基因在不同单细胞亚群的表达量