link + minor fixes

satijalab · satijalab · Feb 1, 2021 · Jan 1, 2021 · Jan 4, 2021 · Jan 4, 2021
commit 9a198c2d0b45b22487f6262400ebaf029181e59e
diff --git a/vignettes/de_vignette.Rmd b/vignettes/de_vignette.Rmd
@@ -55,7 +55,7 @@ head(monocyte.de.markers)
 The results data frame has the following columns : 
 
  * p_val : p_val (unadjusted)
- * avg_logFC : log fold-chage of the average expression between the two groups. Positive values indicate that the feature is more highly expressed in the first group.
+ * avg_logFC : log fold-change of the average expression between the two groups. Positive values indicate that the feature is more highly expressed in the first group.
  * pct.1 : The percentage of cells where the feature is detected in the first group
  * pct.2 : The percentage of cells where the feature is detected in the second group
  * p_val_adj : Adjusted p-value, based on Bonferroni correction using all features in the dataset.

diff --git a/vignettes/future_vignette.Rmd b/vignettes/future_vignette.Rmd
@@ -108,7 +108,7 @@ ggplot(timing.comparisons, aes(fxn, time)) + geom_bar(aes(fill = strategy), stat
 # Frequently asked questions
 
 1. **Where did my progress bar go?**
-<br>Unfortantely, the when running these functions in any of the parallel plan modes you will lose the progress bar. This is due to some technical limitations in the `future` framework and R generally. If you want to monitor function progress, you'll need to forgo parallelization and use `plan("sequential")`.
+<br>Unfortunately, the when running these functions in any of the parallel plan modes you will lose the progress bar. This is due to some technical limitations in the `future` framework and R generally. If you want to monitor function progress, you'll need to forgo parallelization and use `plan("sequential")`.
 
 2. **What should I do if I keep seeing the following error?**
 ```

diff --git a/vignettes/hashing_vignette.Rmd b/vignettes/hashing_vignette.Rmd
@@ -104,7 +104,7 @@ pbmc.hashtag <- ScaleData(pbmc.hashtag, features = VariableFeatures(pbmc.hashtag
 
 ## Adding HTO data as an independent assay
 
-You can read more about working with multi-modal data [here](https://satijalab.org/seurat/multimodal_vignette.html)
+You can read more about working with multi-modal data [here](multimodal_vignette.html)
 
 ```{r hto_assay}
 # Add HTO data as a new assay independent from RNA

diff --git a/vignettes/integration_introduction.Rmd b/vignettes/integration_introduction.Rmd
@@ -88,8 +88,8 @@ We then identify anchors using the `FindIntegrationAnchors()` function, which ta
 ```{r find.anchors}
 immune.anchors <- FindIntegrationAnchors(object.list = ifnb.list, anchor.features = features)
 ```
-```{r integrate.data}
 
+```{r integrate.data}
 # this command creates an 'integrated' data assay
 immune.combined <- IntegrateData(anchorset = immune.anchors)
 ```
@@ -99,7 +99,6 @@ immune.combined <- IntegrateData(anchorset = immune.anchors)
 Now we can run a single integrated analysis on all cells!
 
 ```{r clustering, results='hide', message=FALSE}
-
 # specify that we will perform downstream analysis on the corrected data
 # note that the original unmodified data still resides in the 'RNA' assay
 DefaultAssay(immune.combined) <- "integrated"
@@ -140,14 +139,12 @@ We can explore these marker genes for each cluster and use them to annotate our
 
 ```{r annotate, results = 'hide', message=FALSE, fig.height = 8}
 FeaturePlot(immune.combined, features = c("CD3D", "SELL", "CREM", "CD8A", "GNLY", "CD79A", "FCGR3A", "CCL2", "PPBP"), min.cutoff = "q9")
-
-immune.combined <- RenameIdents(immune.combined, "0" = "CD14 Mono", "1" = "CD4 Naive T", "2" = "CD4 Memory T", "3" = "CD16 Mono", "4" = "B", "5" = "CD8 T", "6" = "NK" , "7" = "T activated", "8" = "DC", "9" = "B Activated", "10" = "Mk", "11" = "pDC", "12" = "Mono/Mk Doublets", "13" = "Eryth")
+immune.combined <- RenameIdents(immune.combined, "0" = "CD14 Mono", "1" = "CD4 Naive T", "2" = "CD4 Memory T", "3" = "CD16 Mono", "4" = "B", "5" = "CD8 T", "6" = "NK" , "7" = "T activated", "8" = "DC", "9" = "B Activated", "10" = "Mk", "11" = "pDC", "12" = "Eryth", "13" = "Mono/Mk Doublets")
 ifnb.list <- lapply(X = ifnb.list, FUN = function(x) {
   x[['seurat_annotations']] <- Idents(immune.combined)[Cells(x)]
   Idents(x) <- 'seurat_annotations'
   x
 })
-
 DimPlot(immune.combined, label = TRUE)
 ```
 
@@ -163,6 +160,12 @@ DotPlot(immune.combined, features = markers.to.plot, cols = c('blue', 'red'),
                       dot.scale = 8, split.by = "stim") + RotatedAxis()
 ```
 
+```{r save.img, include = FALSE}
+plot <- DotPlot(immune.combined, features = markers.to.plot, cols = c('blue', 'red'),
+                      dot.scale = 6, split.by = "stim") + RotatedAxis() 
+ggsave(filename = "../output/images/pbmc_alignment.png", height = 7, width = 12, plot = plot)
+```
+
 ### Identify differential expressed genes across conditions
 
 Now that we've aligned the stimulated and control cells, we can start to do comparative analyses and look at the differences induced by stimulation. One way to look broadly at these changes is to plot the average expression of both the stimulated and control cells and look for genes that are visual outliers on a scatter plot. Here, we take the average expression of both the stimulated and control naive T cells and CD14 monocyte populations and generate the scatter plots, highlighting genes that exhibit dramatic responses to interferon stimulation. 
@@ -201,11 +204,12 @@ b.interferon.response <- FindMarkers(immune.combined, ident.1 = "B_STIM", ident.
 head(b.interferon.response, n = 15)
 ```
 
-Another useful way to visualize these changes in gene expression is with the `split.by` option to the `FeaturePlot` or `VlnPlot` function. This will display FeaturePlots of the list of given genes, split by a grouping variable (stimulation condition here). Genes such as CD3D and GNLY are canonical cell type markers (for T cells and NK/CD8 T cells) that are virtually unaffected by interferon stimulation and display similar gene expression patterns in the control and stimulated group. IFI6 and ISG15, on the other hand, are core interferon response genes and are upregulated accordingly in all cell types. Finally, CD14 and CXCL10 are genes that show a cell type specific interferon response. CD14 expression decreases after stimulation in CD14 monocytes, which could lead to misclassification in a supervised analysis framework, underscoring the value of integrated analysis. CXCL10 shows a distinct upregulation in monocytes and B cells after interferon stimulation but not in other cell types. 
+Another useful way to visualize these changes in gene expression is with the `split.by` option to the `FeaturePlot()` or `VlnPlot()` function. This will display FeaturePlots of the list of given genes, split by a grouping variable (stimulation condition here). Genes such as CD3D and GNLY are canonical cell type markers (for T cells and NK/CD8 T cells) that are virtually unaffected by interferon stimulation and display similar gene expression patterns in the control and stimulated group. IFI6 and ISG15, on the other hand, are core interferon response genes and are upregulated accordingly in all cell types. Finally, CD14 and CXCL10 are genes that show a cell type specific interferon response. CD14 expression decreases after stimulation in CD14 monocytes, which could lead to misclassification in a supervised analysis framework, underscoring the value of integrated analysis. CXCL10 shows a distinct upregulation in monocytes and B cells after interferon stimulation but not in other cell types. 
 
 ```{r feature.heatmaps, fig.height = 14}
 FeaturePlot(immune.combined, features = c("CD3D", "GNLY", "IFI6"), split.by = "stim", max.cutoff = 3, cols = c("grey", "red"))
 ```
+
 ```{r splitvln, fig.height = 12}
 plots <- VlnPlot(immune.combined, features = c("LYZ", "ISG15", "CXCL10"), split.by = "stim", group.by = "celltype", pt.size = 0, combine = FALSE)
 wrap_plots(plots = plots, ncol = 1)
@@ -215,14 +219,6 @@ wrap_plots(plots = plots, ncol = 1)
 saveRDS(immune.combined, file = "../output/immune.combined.rds")
 ```
 
-```{r save.img, include = FALSE}
-plot <- DimPlot(immune.combined, group.by = "stim") +
-  xlab("UMAP 1") + ylab("UMAP 2") + 
-  theme(axis.title = element_text(size = 18), legend.text = element_text(size = 18)) + 
-  guides(colour = guide_legend(override.aes = list(size = 10)))
-ggsave(filename = "../output/images/pbmc_alignment.png", height = 7, width = 12, plot = plot)
-```
-
 ```{r save.times, include = FALSE}
 write.csv(x = t(as.data.frame(all_times)), file = "../output/timings/immune_alignment_times.csv")
 ```
@@ -233,20 +229,17 @@ In [Hafemeister and Satija, 2019](https://genomebiology.biomedcentral.com/articl
 
 Below, we demonstrate how to modify the Seurat integration workflow for datasets that have been normalized with the sctransform workflow. The commands are largely similar, with a few key differences:
 
-* Normalize datasets individually by `SCTransform`, instead of `NormalizeData` prior to integration
+* Normalize datasets individually by `SCTransform()`, instead of `NormalizeData()` prior to integration
 * As discussed further in our [SCTransform vignette](sctransform_vignette.html), we typically use 3,000 or more features for  analysis downstream of sctransform.
-* Run the `PrepSCTIntegration` function prior to identifying anchors
-* When running `FindIntegrationAnchors`, and `IntegrateData`, set the `normalization.method` parameter to the value `SCT`.
-* When running sctransform-based workflows, including integration, do not run the `ScaleData` function 
+* Run the `PrepSCTIntegration()` function prior to identifying anchors
+* When running `FindIntegrationAnchors()`, and `IntegrateData()`, set the `normalization.method` parameter to the value `SCT`.
+* When running sctransform-based workflows, including integration, do not run the `ScaleData()` function 
 
 
 ```{r panc8.cca.sct.init, results='hide', message=FALSE, fig.keep='none'}
 LoadData('ifnb')
 ifnb.list <- SplitObject(ifnb, split.by = "stim")
-
-ifnb.list <- lapply(X = ifnb.list, FUN = function(x) {
-  x <- SCTransform(x)
-})
+ifnb.list <- lapply(X = ifnb.list, FUN = SCTransform) 
 features <- SelectIntegrationFeatures(object.list = ifnb.list, nfeatures = 3000)
 ifnb.list <- PrepSCTIntegration(object.list = ifnb.list, anchor.features = features)
 ```
@@ -268,3 +261,14 @@ p1 + p2
 ```
 
 Now that the datasets have been integrated, you can follow the previous steps in this vignette identify cell types and cell type-specific responses.
+
+```{r save.times, include = FALSE}
+write.csv(x = t(as.data.frame(all_times)), file = "../output/timings/integration_introduction.csv")
+```
+
+<details>
+  <summary>**Session Info**</summary>
+```{r}
+sessionInfo()
+```
+</details>
diff --git a/vignettes/integration_introduction_rpca.Rmd → vignettes/integration_rpca.Rmd b/vignettes/integration_introduction_rpca.Rmd → vignettes/integration_rpca.Rmd
@@ -1,5 +1,5 @@
 ---
-title: 'Introduction to scRNA-seq integration'
+title: 'Fast integration using reciprocal PCA (RPCA)'
 output:
   html_document:
     theme: united
@@ -30,9 +30,7 @@ knitr::opts_chunk$set(
 )
 ```
 
-# Fast integration using reciprocal PCA (RPCA)
-
-In this vignette, we present a slightly modified workflow for the integration of scRNA-seq datasets. Instead of utilizing canonical correlation analysis (‘CCA’) to identify anchors, we instead utilize reciprocal PCA (‘RPCA’). When determining anchors between any two datasets using RPCA, we project each dataset into the other's PCA space and constrain the anchors by the same mutual neighborhood requirement. The commands for both workflows are largely identical, but the two methods may be applied in different context.
+In this vignette, we present a slightly modified workflow for the integration of scRNA-seq datasets. Instead of utilizing canonical correlation analysis (‘CCA’) to identify anchors, we instead utilize reciprocal PCA (‘RPCA’). When determining anchors between any two datasets using RPCA, we project each dataset into the others PCA space and constrain the anchors by the same mutual neighborhood requirement. The commands for both workflows are largely identical, but the two methods may be applied in different context.
 
 By identifying shared sources of variation between datasets, CCA is well-suited for identifying anchors when cell types are conserved, but there are very substantial differences in gene expression across experiments. CCA-based integration therefore enables integrative analysis when experimental conditions or disease states introduce very strong expression shifts, or when integrating datasets across modalities and species. However, CCA-based integration may also lead to overcorrection, especially when a large proportion of cells are non-overlapping across datasets.
 
@@ -71,25 +69,24 @@ ifnb.list <- lapply(X = ifnb.list, FUN = function(x) {
 })
 ```
 
-## Perform integration
+# Perform integration
 
 We then identify anchors using the `FindIntegrationAnchors()` function, which takes a list of Seurat objects as input, and use these anchors to integrate the two datasets together with `IntegrateData()`.
 
 ```{r find.anchors}
 immune.anchors <- FindIntegrationAnchors(object.list = ifnb.list, anchor.features = features,reduction = 'rpca')
 ```
-```{r integrate.data}
 
+```{r integrate.data}
 # this command creates an 'integrated' data assay
 immune.combined <- IntegrateData(anchorset = immune.anchors)
 ```
 
-## Perform an integrated analysis
+# Perform an integrated analysis
 
 Now we can run a single integrated analysis on all cells!
 
 ```{r clustering, results='hide', message=FALSE}
-
 # specify that we will perform downstream analysis on the corrected data
 # note that the original unmodified data still resides in the 'RNA' assay
 DefaultAssay(immune.combined) <- "integrated"
@@ -122,34 +119,36 @@ immune.combined <- RunPCA(immune.combined, npcs = 30, verbose = FALSE)
 immune.combined <- RunUMAP(immune.combined, reduction = "pca", dims = 1:30)
 immune.combined <- FindNeighbors(immune.combined, reduction = "pca", dims = 1:30)
 immune.combined <- FindClusters(immune.combined, resolution = 0.5)
-
 ```
+
 ```{r viz, results='hide', message=FALSE}
 # Visualization
 p1 <- DimPlot(immune.combined, reduction = "umap", group.by = "stim")
 p2 <- DimPlot(immune.combined, reduction = "umap", label = TRUE, repel = TRUE)
 p1 + p2
 ```
 
+```{r save.img, include = FALSE}
+plot <- DimPlot(immune.combined, group.by = "stim") +
+  xlab("UMAP 1") + ylab("UMAP 2") + 
+  theme(axis.title = element_text(size = 18), legend.text = element_text(size = 18)) + 
+  guides(colour = guide_legend(override.aes = list(size = 10)))
+ggsave(filename = "../output/images/rpca_integration.png", height = 7, width = 12, plot = plot)
+```
+
 Now that the datasets have been integrated, you can follow the previous steps in the [introduction to scRNA-seq integration vignette](integration_introduction.html) to identify cell types and cell type-specific responses.
 
-## Performing integration on datasets normalized with SCTransform
+# Performing integration on datasets normalized with SCTransform
 
-As an additional example, we repeat the analyses performed above, but normalize the datasets using [SCTransform](sctransform_vignette.html). We may choose to set the `method` parameter to `glmGamPoi` (install [here](https://bioconductor.org/packages/release/bioc/html/glmGamPoi.html)) in order to enable faster estimation of regression parameters in `SCTransform`.  
+As an additional example, we repeat the analyses performed above, but normalize the datasets using [SCTransform](sctransform_vignette.html). We may choose to set the `method` parameter to `glmGamPoi` (install [here](https://bioconductor.org/packages/release/bioc/html/glmGamPoi.html)) in order to enable faster estimation of regression parameters in `SCTransform()`.  
 
 ```{r panc8.cca.sct.init, results='hide', message=FALSE, fig.keep='none'}
 LoadData('ifnb')
 ifnb.list <- SplitObject(ifnb, split.by = "stim")
-
-ifnb.list <- lapply(X = ifnb.list, FUN = function(x) {
-  x <- SCTransform(x, method = 'glmGamPoi')
-})
+ifnb.list <- lapply(X = ifnb.list, FUN = SCTransform, method = "glmGamPoi")
 features <- SelectIntegrationFeatures(object.list = ifnb.list, nfeatures = 3000)
 ifnb.list <- PrepSCTIntegration(object.list = ifnb.list, anchor.features = features)
-
-ifnb.list <- lapply(X = ifnb.list, FUN = function(x) {
-  x <- RunPCA(x)
-})
+ifnb.list <- lapply(X = ifnb.list, FUN = RunPCA)
 ```
 
 ```{r ifnb.cca.sct.anchors}
@@ -168,3 +167,15 @@ p1 <- DimPlot(immune.combined.sct, reduction = "umap", group.by = "stim")
 p2 <- DimPlot(immune.combined.sct, reduction = "umap", group.by = 'seurat_annotations',label = TRUE, repel = TRUE)
 p1 + p2
 ```
+
+```{r save.times, include = FALSE}
+write.csv(x = t(as.data.frame(all_times)), file = "../output/timings/integration_rpca.csv")
+```
+
+<details>
+  <summary>**Session Info**</summary>
+```{r}
+sessionInfo()
+```
+</details>
+
diff --git a/vignettes/interaction_vignette.Rmd b/vignettes/interaction_vignette.Rmd
@@ -31,11 +31,13 @@ knitr::opts_chunk$set(
 
 # Load in the data
 
-This vignette demonstrates some useful features for interacting with the Seurat object. For demonstration purposes, we will be using the 2,700 PBMC object that is created in the first guided tutorial. You can download the pre-computed object [here](https://www.dropbox.com/s/63gnlw45jf7cje8/pbmc3k_final.rds?dl=1). To simulate the scenario where we have two replicates, we will randomly assign half the cells in each cluster to be from "rep1" and other half from "rep2".
+This vignette demonstrates some useful features for interacting with the Seurat object. For demonstration purposes, we will be using the 2,700 PBMC object that is created in the first guided tutorial. You can load the data from our [SeuratData](https://github.com/satijalab/seurat-data) package. To simulate the scenario where we have two replicates, we will randomly assign half the cells in each cluster to be from "rep1" and other half from "rep2".
 
 ```{r load_data}
 library(Seurat)
-pbmc <- readRDS(file = "../data/pbmc3k_final.rds")
+library(SeuratData)
+InstallData("pbmc3k")
+pbmc <- LoadData(ds = pbmc3k.final)
 
 # pretend that cells were originally assigned to one of two replicates (we assign randomly here)
 # if your cells do belong to multiple replicates, and you want to add this info to the Seurat object

diff --git a/vignettes/merge_vignette.Rmd b/vignettes/merge_vignette.Rmd
@@ -63,10 +63,12 @@ table(pbmc.combined$orig.ident)
 
 # Merging More Than Two `Seurat` Objects
 
-To merge more than two `Seurat` objects, simply pass a vector of multiple `Seurat` objects to the `y` parameter for `merge`; we'll demonstrate this using the 4K and 8K PBMC datasets as well as our previously computed Seurat object from the 2,700 PBMC tutorial (download [here](https://www.dropbox.com/s/63gnlw45jf7cje8/pbmc3k_final.rds?dl=1)).
+To merge more than two `Seurat` objects, simply pass a vector of multiple `Seurat` objects to the `y` parameter for `merge`; we'll demonstrate this using the 4K and 8K PBMC datasets as well as our previously computed Seurat object from the 2,700 PBMC tutorial (loaded via the [SeuratData](https://github.com/satijalab/seurat-data) package).
 
 ```{r merge_three}
-pbmc3k <- readRDS(file = '../data/pbmc3k_final.rds')
+library(SeuratData)
+InstallData("pbmc3k")
+pbmc3k <- LoadData(ds = "pbmc3k.final")
 pbmc3k
 
 pbmc.big <- merge(pbmc3k, y = c(pbmc4k, pbmc8k), add.cell.ids = c('3K', '4K', '8K'), project = 'PBMC15K')