为了演示目的，我们使用通过SeuratData包提供的 2，700个 PBMC 。）
加载数据
library(Seurat)
library(SeuratData)
pbmc <- LoadData("pbmc3k", type = "pbmc3k.final")
执行默认差异表达检测
Seurat 的大部分差异表达功能可以通过FindMarkers()功能访问。作为默认值，Seurat 执行Wilcoxon rank sum test。要测试两个特定细胞组之间的差异表达，可指定ident参数。
# list options for groups to perform differential expression on
levels(pbmc)
## [1] "Naive CD4 T" "Memory CD4 T" "CD14+ Mono" "B" "CD8 T"
## [6] "FCGR3A+ Mono" "NK" "DC" "Platelet"
# Find differentially expressed features between CD14+ and FCGR3A+ Monocytes
monocyte.de.markers <- FindMarkers(pbmc, ident.1 = "CD14+ Mono", ident.2 = "FCGR3A+ Mono")
# view results
head(monocyte.de.markers)
## p_val avg_log2FC pct.1 pct.2 p_val_adj
## FCGR3A 1.193617e-101 -3.776553 0.131 0.975 1.636926e-97
## LYZ 8.134552e-75 2.614275 1.000 0.988 1.115572e-70
## RHOC 4.479768e-68 -2.325013 0.162 0.864 6.143554e-64
## S100A8 7.471811e-65 3.766437 0.975 0.500 1.024684e-60
## S100A9 1.318422e-64 3.299060 0.996 0.870 1.808084e-60
## IFITM2 4.821669e-64 -2.085807 0.677 1.000 6.612437e-60
结果有以下列：
- p_val ： p_val （未经校正）
- avg_log2FC：记录两组之间平均表达的倍数变化。正值表示该基因在第一组中表达更高。
- pct.1 ： 在第一组检测到该基因的细胞百分比
- pct.2 ： 在第二组检测到该基因的细胞百分比
- p_val_adj： 校正后的 p 值，基于使用数据集中的所有基因的Bonferroni校正。
如果参数被省略或设置为 NULL，FindMarkers()将执行所指定的组和所有其他组之间的比较。
# Find differentially expressed features between CD14+ Monocytes and all other cells, only
# search for positive markers
monocyte.de.markers <- FindMarkers(pbmc, ident.1 = "CD14+ Mono", ident.2 = NULL, only.pos = TRUE)
# view results
head(monocyte.de.markers)
## p_val avg_log2FC pct.1 pct.2 p_val_adj
## S100A9 0.000000e+00 5.570063 0.996 0.215 0.000000e+00
## S100A8 0.000000e+00 5.477394 0.975 0.121 0.000000e+00
## FCN1 0.000000e+00 3.394219 0.952 0.151 0.000000e+00
## LGALS2 0.000000e+00 3.800484 0.908 0.059 0.000000e+00
## CD14 2.856582e-294 2.815626 0.667 0.028 3.917516e-290
## TYROBP 3.190467e-284 3.046798 0.994 0.265 4.375406e-280
预过滤基因或细胞，以提高差异基因检测的速度
为了提高marker检测的速度，特别是对于大型数据集，Seurat 允许对基因或细胞进行预过滤。例如，在两组细胞中很少检测到的基因，或在平均水平表达类似的基因，不太可能有差异表达。下面演示了几个参数的使用。
min.pct
logfc.threshold
min.diff.pct
max.cells.per.ident
# Pre-filter features that are detected at <50% frequency in either CD14+ Monocytes or FCGR3A+
# Monocytes
head(FindMarkers(pbmc, ident.1 = "CD14+ Mono", ident.2 = "FCGR3A+ Mono", min.pct = 0.5))
## p_val avg_log2FC pct.1 pct.2 p_val_adj
## FCGR3A 1.193617e-101 -3.776553 0.131 0.975 1.636926e-97
## LYZ 8.134552e-75 2.614275 1.000 0.988 1.115572e-70
## RHOC 4.479768e-68 -2.325013 0.162 0.864 6.143554e-64
## S100A8 7.471811e-65 3.766437 0.975 0.500 1.024684e-60
## S100A9 1.318422e-64 3.299060 0.996 0.870 1.808084e-60
## IFITM2 4.821669e-64 -2.085807 0.677 1.000 6.612437e-60
# Pre-filter features that have less than a two-fold change between the average expression of
# CD14+ Monocytes vs FCGR3A+ Monocytes
head(FindMarkers(pbmc, ident.1 = "CD14+ Mono", ident.2 = "FCGR3A+ Mono", logfc.threshold = log(2)))
## p_val avg_log2FC pct.1 pct.2 p_val_adj
## FCGR3A 1.193617e-101 -3.776553 0.131 0.975 1.636926e-97
## LYZ 8.134552e-75 2.614275 1.000 0.988 1.115572e-70
## RHOC 4.479768e-68 -2.325013 0.162 0.864 6.143554e-64
## S100A8 7.471811e-65 3.766437 0.975 0.500 1.024684e-60
## S100A9 1.318422e-64 3.299060 0.996 0.870 1.808084e-60
## IFITM2 4.821669e-64 -2.085807 0.677 1.000 6.612437e-60
# Pre-filter features whose detection percentages across the two groups are similar (within
# 0.25)
head(FindMarkers(pbmc, ident.1 = "CD14+ Mono", ident.2 = "FCGR3A+ Mono", min.diff.pct = 0.25))
## p_val avg_log2FC pct.1 pct.2 p_val_adj
## FCGR3A 1.193617e-101 -3.776553 0.131 0.975 1.636926e-97
## RHOC 4.479768e-68 -2.325013 0.162 0.864 6.143554e-64
## S100A8 7.471811e-65 3.766437 0.975 0.500 1.024684e-60
## IFITM2 4.821669e-64 -2.085807 0.677 1.000 6.612437e-60
## LGALS2 1.034540e-57 2.956704 0.908 0.265 1.418768e-53
## CDKN1C 2.886353e-48 -1.453845 0.029 0.506 3.958345e-44
# Increasing min.pct, logfc.threshold, and min.diff.pct, will increase the speed of DE testing,
# but could also miss features that are prefiltered
# Subsample each group to a maximum of 200 cells. Can be very useful for large clusters, or
# computationally-intensive DE tests
head(FindMarkers(pbmc, ident.1 = "CD14+ Mono", ident.2 = "FCGR3A+ Mono", max.cells.per.ident = 200))
## p_val avg_log2FC pct.1 pct.2 p_val_adj
## FCGR3A 4.725441e-61 -3.776553 0.131 0.975 6.480470e-57
## LYZ 6.442891e-56 2.614275 1.000 0.988 8.835781e-52
## S100A8 8.983226e-49 3.766437 0.975 0.500 1.231960e-44
## S100A9 1.812278e-47 3.299060 0.996 0.870 2.485358e-43
## IFITM2 1.185202e-45 -2.085807 0.677 1.000 1.625386e-41
## RPS19 1.685374e-44 -1.091150 0.990 1.000 2.311321e-40
使用其他检测方法执行差异基因分析
当前支持以下差异基因检测方法：
- “wilcox”
- “bimod”
- “roc”
- “Student’s t-test”
- “poisson”
- “negbinom”
- “LR”
- “MAST”
- “DESeq2”
# Test for DE features using the MAST package
head(FindMarkers(pbmc, ident.1 = "CD14+ Mono", ident.2 = "FCGR3A+ Mono", test.use = "MAST"))
## p_val avg_log2FC pct.1 pct.2 p_val_adj
## LYZ 7.653066e-145 2.614275 1.000 0.988 1.049541e-140
## FCGR3A 2.897172e-119 -3.776553 0.131 0.975 3.973182e-115
## S100A9 2.123928e-95 3.299060 0.996 0.870 2.912755e-91
## S100A8 3.279521e-92 3.766437 0.975 0.500 4.497535e-88
## IFITM2 5.591175e-87 -2.085807 0.677 1.000 7.667737e-83
## LGALS2 1.132854e-75 2.956704 0.908 0.265 1.553596e-71
# Test for DE features using the DESeq2 package. Throws an error if DESeq2 has not already been
# installed Note that the DESeq2 workflows can be computationally intensive for large datasets,
# but are incompatible with some feature pre-filtering options We therefore suggest initially
# limiting the number of cells used for testing
head(FindMarkers(pbmc, ident.1 = "CD14+ Mono", ident.2 = "FCGR3A+ Mono", test.use = "DESeq2", max.cells.per.ident = 50))
## p_val avg_log2FC pct.1 pct.2 p_val_adj
## S100A9 1.887262e-58 2.538360 0.996 0.870 2.588191e-54
## LYZ 4.198394e-52 1.987962 1.000 0.988 5.757678e-48
## S100A8 5.747352e-49 2.784248 0.975 0.500 7.881918e-45
## FCGR3A 1.315842e-35 -2.949992 0.131 0.975 1.804546e-31
## RPS19 1.514561e-33 -1.614892 0.990 1.000 2.077068e-29
## IFITM2 1.227042e-26 -2.212583 0.677 1.000 1.682765e-22
