继续上面MetaNeighbor 2的学习,我们做第三组场景的测试,来测试不同基因集的选择对于结果的影响。
========测试数据3============
go_sets <- readRDS("go_mouse.rds")
biccn_gaba <- readRDS("biccn_gaba.rds")
dim(biccn_gaba)
下面,我们就从GO里面筛选数据集出现的基因。作者paper中建议:10-100个之间。
known_genes <- rownames(biccn_gaba)
go_sets <- lapply(go_sets,function(gene_set){
gene_set[gene_set %in% known_genes]}
)
min_size=10
max_size=100
go_set_size <- sapply(go_sets,length)
go_sets <- go_sets[go_set_size >= min_size & go_set_size <= max_size]
筛选完基因,下面就是用用户提供的参考基因集,进行MetaNeighbor的运行了,其实和前面类似。
aurocs <- MetaNeighbor(dat = biccn_gaba,
experiment_labels = biccn_gaba$study_id,
celltype_labels = biccn_gaba$joint_subclass_label,
genesets=go_sets,
fast_version=TRUE,
bplot=FALSE,batch_size=50)
write.table(aurocs,"functional_aurocs.txt")
结果就是,利用不同GO里面的基因集,各个study的AUROC。
plotBPlot(head(aurocs,100))
我们把前100个数据画成小提琴图看下。从结果可看出,大多数基因集对可重复性的贡献很小(AUROC~0.7),许多基因集的性能接近随机(AUROC ~0.5–0.6),一些基因集具有极高的性能(AUROC>0.8)。
然后,我们看下排名较高的基因集是那些。
gs_size <- sapply(go_sets,length)
aurocs_df <- data.frame(go_term=rownames(aurocs),aurocs)
aurocs_df$average <- rowMeans(aurocs)
aurocs_df$n_genes <- gs_size[rownames(aurocs)]
从排名较高的前10个GO可以看出,主要是那些和inhibitory neruons相关的基因。
也可以像单细胞数据一样查看GO下面的每个基因在不同study里面的平均表达值,从而观察这些基因在每个组间的特异性。
网友评论