从3000左右的外周血单核细胞中拿到CD4,CD8,NK这些T细胞类似的细胞亚群,去看一下文章中那些基因在这些亚群种有没有很好的体现.

    1. #任意文献的单细胞亚群标记基因
    2. rm(list = ls())
    3. library(Seurat)
    4. library(SeuratData)
    5. library(ggplot2)
    6. library(patchwork)
    7. library(dplyr)
    8. load(file = 'basic.sce.pbmc.Rdata')
    9. levels(Idents(pbmc))
    11. #首先提取T细胞子集
    12. DimPlot(pbmc, 
    13.         reduction = 'umap', 
    14.         label = TRUE, 
    15.         pt.size = 0.5) + NoLegend()
    16. FeaturePlot(pbmc, 
    17.             features = c("CD3D","CD3E")) #T细胞常用CD3来分
    18. sce=pbmc
    19. table(Idents(sce))
    20. #取子集
    21. t_sce = sce[, Idents(sce) %in% c( "Naive CD4 T",
    22.                                   "Memory CD4 T",
    23.                                   'CD8 T',
    24.                                   'NK')]
    25. #然后进行标准的降维聚类分群
    26. #代码不要变动
    27. sce=t_sce
    28. sce <- NormalizeData(sce, 
    29.                      normalization.method = "LogNormalize",
    30.                      scale.factor = 1e4) 
    31. sce <- FindVariableFeatures(sce, 
    32.                             selection.method = 'vst',
    33.                             nfeatures = 2000)
    34. sce <- ScaleData(sce,
    35.                  vars.to.regress = "percent.mt")
    36. sce <- RunPCA(sce, 
    37.               features = VariableFeatures(object = sce)) 
    38. sce <- FindNeighbors(sce, 
    39.                      dims = 1:10)
    40. sce <- FindClusters(sce,
    41.                     resolution = 1 )
    43. # Look at cluster IDs of the first 5 cells
    44. head(Idents(sce), 5)
    45. table(sce$seurat_clusters) 
    46. sce <- RunUMAP(sce, 
    47.                dims = 1:10)
    48. DimPlot(sce, 
    49.         reduction = 'umap')
    51. sce.markers <- FindAllMarkers(object = sce, 
    52.                               only.pos = TRUE, 
    53.                               min.pct = 0.25, 
    54.                               thresh.use = 0.25)
    55. DT::datatable(sce.markers)
    56. # write.csv(sce.markers,file=paste0(pro,'_sce.markers.csv'))
    58. library(dplyr) 
    59. # 在seurat V4里面,是 avg_log2FC,是V3版本是avg_logFC
    60. top10 <- sce.markers %>% group_by(cluster) %>% top_n(5, avg_log2FC)
    61. DoHeatmap(sce,
    62.           top10$gene,
    63.           size=3)  
    64. p <- DotPlot(sce,
    65.              features = unique(top10$gene),
    66.              assay='RNA' )  + coord_flip()
    68. p
    69. ggsave('check2-top5-for-all-Tcells.pdf',height = 18)
    71. # 然后查看文献的标记基因
    72. # 2021年5月的文章：《A single-cell map of intratumoral changes during anti-PD1 treatment of patients with breast cancer
    73. # 参考：https://mp.weixin.qq.com/s/FhEASxwTF8e7lNoKxTuxrg 
    75. #文献中的基因
    76. marker_genes= c("LEF1","TCF7","SELL","IL7R","CD40LG","ANXA1","FOS",
    77.                 "JUN","FOXP3","SAT1","IL2RA","CTLA4","PDCD1","CXCL13","CD200",
    78.                 "TNFRSF18","CCR7","NELL2","CD55","KLF2","TOB1","ZNF683","CCL5",
    79.                 "GZMK","EOMES","ITM2C","CX3CR1","GNLY","GZMH","GZMB","LAG3","CCL4L2",
    80.                 "FCGR3A","FGFBP2","TYROBP","AREG","XCL1","KLRC1","TRDV2","TRGV9","MTRNR2L8",
    81.                 "KLRD1","TRDV1","KLRC3",
    82.                 "CTSW","CD7","MKI67","STMN1","TUBA1B","HIST1H4C" )
    85. p <- DotPlot(sce, 
    86.              features = marker_genes,
    87.              assay='RNA'  )  + coord_flip() + theme(axis.text.x = element_text(angle = 90))
    89. p
    90. ggsave('check_g1_markers_by_Tcell-SubType.pdf')

    image.png

    #文献中的基因
marker_genes =c("CD3D","CD4","CD8A","CCR7","LEF1","SELL","TCF7","GNLY","IFNG","NKG7","PRF1",
                "GZMA","GZMB","GZMH","GZMK","HAVCR2","LAG3","PDCD1","CTLA4","TIGIT","BTLA","KLRC1",
                "ANXA1","ANKRD28","CD69","CD40LG","ZNF683","FOXP3","IL2RA","IKZF2","NCR1","NCAM1",
                "TYROBP","KLRD1","KLRF1","KLRB1","CX3CR1","FCGR3A",
                "XCL1","XCL2","TRDV2","TRGV9","TRGC2","MKI67","TOP2A")

p <- DotPlot(sce, 
             features = marker_genes,
             assay='RNA'  )  + coord_flip() + theme(axis.text.x = element_text(angle = 90))

p
ggsave('check_g2_markers_by_Tcell-SubType.pdf')

    image.png

    # naive (LEF1, SELL, TCF7),
# effector (IFNG), 
# cytotoxicity (GZMB, PRF1), 
# early and general exhaustion (PDCD1, CTLA4, ENTPD1 ) .
# antigen presentation (CD74, HLA-DRB1/5, HLA-DQA2)

genes_to_check = c('PTPRC', 'CD3D', 'CD3E', 'FOXP3',
                   'CD4','IL7R','NKG7','CD8A',
                   'LEF1', 'SELL', 'TCF7', # naive marker
                   'IFNG','GZMB', 'PRF1',
                   'PDCD1', 'CTLA4', 'ENTPD1'  )
p <- DotPlot(sce, 
             features = genes_to_check,
             assay='RNA'  )  + coord_flip() + theme(axis.text.x = element_text(angle = 90))

p

    image.png