美文网首页单细胞转录组生物信息学R语言源码
单细胞TPM数据的处理及批次效应的去除

单细胞TPM数据的处理及批次效应的去除

作者: 碌碌无为的杰少 | 来源:发表于2021-01-25 23:18 被阅读0次

    harmony真的好用,天才的发明者。复现一篇30多个样本的主图,这个数据我也想了好久,由于作者提供了tpm数据,这个数据相当于标准化之后的,就相当于SCT之后的,所以不需要对样本进行标准化了,主图是这样的。


    图片.png
    exp1 <- data.table::fread('GSE132465_GEO_processed_CRC_10X_natural_log_TPM_matrix.txt',data.table = F)
    rownames(exp1) <- exp1$Index
    exp1 <-exp1[,-1]
    experiment.aggregate <- CreateSeuratObject(
      exp1,
      project = "multi", 
      min.cells = 10,
      min.features = 200,
      names.field = 1,
      names.delim = "_")
    experiment.aggregate[["percent.mt"]] <- PercentageFeatureSet(experiment.aggregate, 
                                                                 pattern = "^MT-")
    view(experiment.aggregate@meta.data)
    experiment.aggregate <- FindVariableFeatures(experiment.aggregate, 
                                                 selection.method = "vst",
                                                 nfeatures = 1000)
    all.genes <- rownames(experiment.aggregate)
    experiment.aggregate <- ScaleData(experiment.aggregate,
                                      features = all.genes, verbose = FALSE
    )
    experiment.aggregate <- RunPCA(object = experiment.aggregate, 
                                   features = VariableFeatures(experiment.aggregate),
                                   verbose = F,npcs = 50)
    experiment.aggregate <- RunHarmony(experiment.aggregate, group.by.vars = "orig.ident")
    A1 <- ElbowPlot(experiment.aggregate, ndims = 50)
    pdf("A1.pdf",width = 6,height = 5)
    
    print(A1)
    dev.off()
    
    dim.use <- 1:20
    experiment.aggregate <- FindNeighbors(experiment.aggregate, dims = dim.use, reduction = "harmony")
    #基于细胞之间的相似性计算具体的cluster
    experiment.aggregate <- FindClusters(experiment.aggregate, resolution = 0.6, reduction = "harmony")
    #TSNE算法聚类(降维)
    experiment.aggregate <- RunTSNE(experiment.aggregate, dims = dim.use, do.fast = TRUE, reduction = "harmony")
    experiment.aggregate@meta.data$celltype3 = "NA"
    
    
    # 更改 celltype 信息,设置细胞群新名称
    experiment.aggregate@meta.data[which(experiment.aggregate@meta.data$seurat_clusters %in% c(1,2,6,8)), "celltype3"] = "T cells"
    experiment.aggregate@meta.data[which(experiment.aggregate@meta.data$seurat_clusters %in% c(7,11,15,19)), "celltype3"] = "Stromal cells"
    experiment.aggregate@meta.data[which(experiment.aggregate@meta.data$seurat_clusters %in% c(3,4,13,18)), "celltype3"] = "B cells"
    experiment.aggregate@meta.data[which(experiment.aggregate@meta.data$seurat_clusters %in% c(0,10,12,14,16,17)), "celltype3"] = "Epithelial cells"
    experiment.aggregate@meta.data[which(experiment.aggregate@meta.data$seurat_clusters %in% c(5,9)), "celltype3"] = "Myloid cells"
    experiment.aggregate@meta.data[which(experiment.aggregate@meta.data$seurat_clusters %in% c(20)), "celltype3"] = "Mast cells"
    
    A2 <- DimPlot(experiment.aggregate, reduction = "tsne", label = T,group.by = 'celltype3') +NoLegend()
    A2
    pdf("A2.pdf",width = 8,height = 6)
    print(A2)
    dev.off()
    
    singleB<-  SubsetData(experiment.aggregate,
                          # 提取数据根据的组名
                          subset.name = "celltype3", 
                          # 提取的组别
                          accept.value = c('B cells'))
    singleT2<-  SubsetData(experiment.aggregate,
                          # 提取数据根据的组名
                          subset.name = "celltype3", 
                          # 提取的组别
                          accept.value = c('T cells'))
    save(singleB, file = "singleB.Rdata")
    save(singleT2, file = "singleT.Rdata")
    

    最后的图是这样的,感觉一点不差于原文的CCA去除批次效应的结果


    图片.png

    相关文章

      网友评论

        本文标题:单细胞TPM数据的处理及批次效应的去除

        本文链接:https://www.haomeiwen.com/subject/fnlwzktx.html