SARS-CoV-2感染的雪貂支气管肺泡灌洗液单细胞转录组数据挖掘（2） 细分NK和CD8的T单细胞亚群

2022-03-03 13:18:41 浏览数 (1)

给学徒们收集整理了几套带GitHub源代码的文献图表合辑，让优秀者一点一滴拆解开来分享给大家。（全部的代码复制粘贴即可运行，欢迎尝试以及批评指正）

现在是雪貂支气管肺泡灌洗液单细胞转录组显示SARS-CoV-2感染期间巨噬细胞的顺序变化专辑第2讲：主要是对 NK cells and CD8 T 进行细分亚群

下面是前年实习生（日行一膳）的分享

1643459432584

本次复现的是于2021年7月27日发表在Nature Communications上的”Single-cell transcriptome of bronchoalveolar lavage fluid reveals sequential change of macrophages during SARS-CoV-2 infection in ferrets“中的Figure2

Fig. 2 Subpopulation analysis of NK cells and CD8 T cells.

1643459512407

代码语言：javascript复制
## 安装所需要的R包
install.packages('Seurat')
install.packages("tidyverse")
install.packages("ggpubr")
install.packages("cowplot")
install.packages("dplyr")
install.packages("ggplot2")
install.packages("tidyr")

## 载入R包
library(Seurat)
library(tidyverse)
library(ggpubr)
library(cowplot)
library(dplyr)
library(ggplot2)
library(tidyr)
## 分析NK细胞与CD8 T细胞
MPcov <- readRDS("Seurat_object_total_cells.Rds")
#常规降维聚类分群(NK细胞与CD8 T细胞)
Idents(MPcov) <- "Annotation"
cov.NK <- subset(MPcov, idents = c("NK cell"))
cov.NK <- RunPCA(cov.NK)
ElbowPlot(cov.NK)
cov.NK <- FindNeighbors(cov.NK, dims = 1:5, force.recalc = T)
cov.NK <- FindClusters(cov.NK, resolution = 0.3)
cov.NK <- RunUMAP(object = cov.NK, dims = 1:5)

cov.CD8 <- subset(MPcov, idents = c("CD8 T cell"))
cov.CD8 <- RunPCA(cov.CD8)
ElbowPlot(cov.CD8)
cov.CD8 <- FindNeighbors(cov.CD8, dims = 1:5, force.recalc = T)
cov.CD8 <- FindClusters(cov.CD8, resolution = 0.2)
cov.CD8 <- RunUMAP(object = cov.CD8, dims = 1:5)
#NK细胞的UMAP
pdf("Fig2a.NK_umap_seurat_clusters.pdf",900/72*0.8, 688/72*0.8) 
DimPlot(cov.NK,group.by = "seurat_clusters", label=T)
dev.off()

1643460147179

代码语言：javascript复制
#未感染对照组(n=3)、2 dpi (n=3)和5 dpi (n=4)各细胞类型在NK细胞群中的比例
pdf("Fig2b.NK_stacked_bar_composition_celltype.pdf",3.8,5)
x = table(cov.NK$seurat_clusters,cov.NK$Experimental_group)
for(i in 1:ncol(x)){
  x[,i] = x[,i] /sum(x[,i])
}
colSums(x)
x = x[nrow(x):1,c("Ct","C2","C5")]
colnames(x) = c("dpi 0","dpi 2","dpi 5")

x <- x[nrow(x):1,]

par(mar=c(11,4,1,1))
barplot(x,las=2,main = "Relative proportion",horiz=T,col = row(x))
dev.off()

1643463486312

代码语言：javascript复制
#STAT1、OAS1、ISG15、GZMB、GZMK、PRF1在NK细胞中的表达水平
pdf("Fig2c.vln_plot_NK.pdf",3.8,5)
VlnPlot(cov.NK, features = "STAT1", pt.size = 0)
VlnPlot(cov.NK, features = "OAS1", pt.size = 0)
VlnPlot(cov.NK, features = "ISG15", pt.size = 0)
VlnPlot(cov.NK, features = "GZMB", pt.size = 0)
VlnPlot(cov.NK, features = "GZMK", pt.size = 0)
VlnPlot(cov.NK, features = "PRF1", pt.size = 0)
dev.off()

1643463573291

代码语言：javascript复制
#CD8  T细胞的UMAP
pdf("Fig2d.CD8_umap_seurat_clusters.pdf",900/72*0.8, 688/72*0.8)
DimPlot(cov.CD8,label=T,group.by="seurat_clusters", pt.size = 2)   coord_fixed()
dev.off()

1643464091561

代码语言：javascript复制
#CD69、S1PR1、ITGAE在CD8 T细胞的表达水平的表达水平
pdf("Fig2e.vln_plot_resident_CD8.pdf",3.8,5)
VlnPlot(cov.CD8, features = "CD69", pt.size = 0)
VlnPlot(cov.CD8, features = "S1PR1", pt.size = 0)
VlnPlot(cov.CD8, features = "ITGAE", pt.size = 0)
dev.off()

1643464166815

代码语言：javascript复制
#OAS1、ISG15、IFNG、GZMB和PRF1在CD8 T细胞的表达水平
pdf("Fig2f.vln_plot_activation_CD8.pdf",3.8,5)
VlnPlot(cov.CD8, features = "OAS1", pt.size = 0)
VlnPlot(cov.CD8, features = "ISG15", pt.size = 0)
VlnPlot(cov.CD8, features = "IFNG", pt.size = 0)
VlnPlot(cov.CD8, features = "GZMB", pt.size = 0)
VlnPlot(cov.CD8, features = "PRF1", pt.size = 0)
dev.off()

1643464217936

代码语言：javascript复制
#SARS-CoV-2雪貂CD8 T细胞在阴性对照，2 dpi和5 dpi处的分布
pdf("Fig2g.condition_umap_CD8.pdf",9,15)
par(mfrow=c(1,3), mar = c(1,1,3,1),pty="s", oma=c(0,0,0,0))
for(gr in c("Ct","C2","C5")){
  plot(Embeddings(cov.CD8,reduction="umap")[cov.CD8$Experimental_group!=gr,],col="#00000020",
       pch=16,cex=0.8,bty="n",xaxt="n",yaxt="n",xlab="",ylab="")
  points(Embeddings(cov.CD8,reduction="umap")[cov.CD8$Experimental_group==gr,],col="#168C4060",pch=16,cex=1)
  title(c("Ct"="dpi 0","C2"="dpi 2", "C5"="dpi 5")[gr])
}
dev.off()

1643464265434

代码语言：javascript复制
#CD8 T细胞中OAS1和ISG15的表达
pdf("Fig2h.featureplot.pdf",12,8*0.8)
goi <- c("OAS1", "ISG15")#,
p = FeaturePlot(cov.CD8, features = goi,reduction="umap",ncol=2,combine = FALSE,cols = c("grey80", "#0A64A4"))
for(i in 1:length(p)) {
  p[[i]] <- p[[i]]   NoLegend()   NoAxes()
}
cowplot::plot_grid(plotlist = p,ncol = 2)
dev.off()

这样的基础认知，也可以看单细胞的基础10讲：

  • 01. 上游分析流程
  • 02.课题多少个样品，测序数据量如何
  • 03. 过滤不合格细胞和基因（数据质控很重要）
  • 04. 过滤线粒体核糖体基因
  • 05. 去除细胞效应和基因效应
  • 06.单细胞转录组数据的降维聚类分群
  • 07.单细胞转录组数据处理之细胞亚群注释
  • 08.把拿到的亚群进行更细致的分群
  • 09.单细胞转录组数据处理之细胞亚群比例比较
数据处理

0 人点赞