<p>I' am doing a sc-Seq analysis on a dataset from <a href="https://www.10xgenomics.com/datasets/10-k-pbm-cs-from-a-healthy-donor-v-3-chemistry-3-standard-3-0-0" rel="nofollow noreferrer">10x genomics</a>, using Seurat. I followed the standard workflow and I also did it with the SCTransform workflow. The results are very different when I plot a Feature Scatter with the same 2 genes.</p>
<p>Standard workflow:
<a href="https://i.sstatic.net/PVcQg.jpg" rel="nofollow noreferrer"><img src="https://i.sstatic.net/PVcQg.jpg" alt=" CD3D and CD3E" /></a></p>
<p>SCTransform workflow:
<a href="https://i.sstatic.net/AH7at.png" rel="nofollow noreferrer"><img src="https://i.sstatic.net/AH7at.png" alt="enter image description here" /></a></p>
<p>So I was wondering why they were so different. I know that the perfect lines in the first one happen because of the normalization and that it usually means that I need more data. But I don't understand why the SCTransform one. looks like it has less data and why it looks so perfectly arranged. My goal is to make a plot that looks normal and natural, and that may or may not have a visual correlation; instead of this weirdly looking and unnatural plots. If you also have a tip for this, I would be very grateful.  Thanks in advance.</p>
<p>Here is my code for the normal workflow:</p>
<pre><code>library(Seurat)
library(tidyverse)
library(ggplot2)
library(gridExtra)
library(harmony)
library(future.apply)
library(cowplot)
library(patchwork)
library(&quot;DESeq2&quot;)
library(SeuratData)
library(SeuratDisk)
library(SeuratWrappers)
library(sctransform)
</code></pre>
<pre><code>data &lt;- Read10X(data.dir = &quot;/Users/rodrigohermoza/Desktop/UTEC/2024-0/PoliSia/rstudio/diferential_expression_2/DE2/filtered_gene_bc_matrices/raw_feature_bc_matrix/&quot;)

pbmc &lt;- CreateSeuratObject(counts = data, project = &quot;DE3&quot;)
pbmc &lt;- PercentageFeatureSet(pbmc, pattern = &quot;^MT-&quot;, col.name = &quot;percent.mt&quot;) 
pbmc &lt;- subset(pbmc, subset = nFeature_RNA &gt; 200 &amp; percent.mt &lt;15 &amp; nCount_RNA &gt; 50)
</code></pre>
<pre><code>pbmc &lt;- NormalizeData(pbmc)
pbmc &lt;- FindVariableFeatures(pbmc, nfeatures = 3000)
pbmc &lt;- ScaleData(pbmc, vars.to.regress = &quot;percent.mt&quot;)
pbmc &lt;- RunPCA(pbmc)
</code></pre>
<pre><code>pbmc &lt;- FindNeighbors(pbmc, dims = 1:21)
pbmc &lt;- FindClusters(pbmc, resolution = 1)
pbmc &lt;- RunUMAP(pbmc, dims = 1:21)
DimPlot(pbmc, label = TRUE)
</code></pre>
<pre><code>sweep_res &lt;- paramSweep(pbmc,PCs = 1:20, sct = FALSE)
sweep_stat &lt;- summarizeSweep(sweep_res, GT = F)
bcmv &lt;- find.pK(sweep_stat)

</code></pre>
<pre><code>ggplot(bcmv, aes(pK, BCmetric, group = 1)) +
  geom_point()+
  geom_line()

#0.005
</code></pre>
<pre><code>pK&lt;- bcmv %&gt;%
  filter(BCmetric == max(BCmetric)) %&gt;%
  select(pK)
pK&lt;- as.numeric(as.character(pK[[1]]))
</code></pre>
<pre><code>annotations &lt;- pbmc@meta.data$seurat_clusters
homotipic &lt;- modelHomotypic(annotations)
exp &lt;- round(0.08*nrow(pbmc@meta.data))
exp.adj &lt;- round(exp*(1-homotipic))
</code></pre>
<pre><code>pbmc &lt;- doubletFinder(pbmc,
                          PCs = 1:20,
                          pK = pK,
                          nExp = exp.adj,
                          reuse.pANN = F,
                          sct = F)
DimPlot(pbmc, reduction = &quot;umap&quot;, group.by = &quot;DF.classifications_0.25_0.11_818&quot;)
</code></pre>
<pre><code>table(pbmc@meta.data$DF.classifications_0.25_0.11_818)
pbmc_filtered &lt;- subset(pbmc, subset = DF.classifications_0.25_0.11_818 == &quot;Singlet&quot;)
</code></pre>
<pre><code>pbmc.markers &lt;- FindAllMarkers(pbmc_filtered)
</code></pre>
<pre><code>FeatureScatter(pbmc_filtered, feature1 = &quot;CD3D&quot;, feature2 = &quot;CD3E&quot;, slot = &quot;data&quot;)
</code></pre>
<p>And here is my SCTransform workflow</p>
<pre><code>library(Seurat)
library(tidyverse)
library(ggplot2)
library(gridExtra)
library(harmony)
library(future.apply)
library(cowplot)
library(patchwork)
library(&quot;DESeq2&quot;)
library(SeuratData)
library(SeuratDisk)
library(SeuratWrappers)
library(sctransform)
library(DoubletFinder)
</code></pre>
<pre><code>data &lt;- Read10X(data.dir = &quot;/Users/rodrigohermoza/Desktop/UTEC/2024-0/PoliSia/rstudio/diferential_expression_2/DE2/filtered_gene_bc_matrices/raw_feature_bc_matrix/&quot;)
pbmc &lt;- CreateSeuratObject(counts = data, project = &quot;DE3&quot;)
pbmc &lt;- PercentageFeatureSet(pbmc, pattern = &quot;^MT-&quot;, col.name = &quot;percent.mt&quot;) 
pbmc &lt;- subset(pbmc, subset = nFeature_RNA &gt; 200 &amp; percent.mt &lt;15 &amp; nCount_RNA &gt; 50)
</code></pre>
<pre><code>pbmc &lt;- SCTransform(pbmc, vars.to.regress = &quot;percent.mt&quot;)
pbmc &lt;- RunPCA(pbmc)
ElbowPlot(pbmc,ndims = 50)
</code></pre>
<pre><code>pbmc &lt;- FindNeighbors(pbmc, dims = 1:20)
pbmc &lt;- FindClusters(pbmc, resolution = 1)
pbmc &lt;- RunUMAP(pbmc, dims = 1:20)
DimPlot(pbmc, label = TRUE)
</code></pre>
<pre><code>pbmc.markers &lt;- FindAllMarkers(pbmc)
</code></pre>
<pre><code>FeatureScatter(pbmc, feature1 = &quot;CD3D&quot;, feature2 = &quot;CD3E&quot;, slot = &quot;data&quot;, smooth = T)
</code></pre>

CC BY-SA 4.0
2024-02-15T21:53:46.690Z
Check in detail what is plotted. The "standard" is probably some continuous type of data, like scaled counts or something like this, and the second one is probably some sort of more transformed integers, which is not as tightly continuous as the "standard" one. Like imagine you do log2 of 1:10, that gives defined breaks, while something like a Z-score gives a much more continuous distribution without "hard breaks".
ATpoint
CC BY-SA 4.0
2024-02-23T17:25:02.693Z
Thanks for the clarification :)
Rodrigo Angel Hermoza Sanchez
<r><scrnaseq><seurat><rstudio>
Why my FeatureScatter looks so different, and what does that mean?
# Why my FeatureScatter looks so different, and what does that mean?

(Asked by: username_0 on 2024-02-15)

I' am doing a sc-Seq analysis on a dataset from 10x genomics, using Seurat. I followed the standard workflow and I also did it with the SCTransform workflow. The results are very different when I plot a Feature Scatter with the same 2 genes.
Standard workflow:

SCTransform workflow:

So I was wondering why they were so different. I know that the perfect lines in the first one happen because of the normalization and that it usually means that I need more data. But I don't understand why the SCTransform one. looks like it has less data and why it looks so perfectly arranged. My goal is to make a plot that looks normal and natural, and that may or may not have a visual correlation; instead of this weirdly looking and unnatural plots. If you also have a tip for this, I would be very grateful.  Thanks in advance.
Here is my code for the normal workflow:
library(Seurat)
library(tidyverse)
library(ggplot2)
library(gridExtra)
library(harmony)
library(future.apply)
library(cowplot)
library(patchwork)
library("DESeq2")
library(SeuratData)
library(SeuratDisk)
library(SeuratWrappers)
library(sctransform)

data <- Read10X(data.dir = "/Users/rodrigohermoza/Desktop/UTEC/2024-0/PoliSia/rstudio/diferential_expression_2/DE2/filtered_gene_bc_matrices/raw_feature_bc_matrix/")

pbmc <- CreateSeuratObject(counts = data, project = "DE3")
pbmc <- PercentageFeatureSet(pbmc, pattern = "^MT-", col.name = "percent.mt") 
pbmc <- subset(pbmc, subset = nFeature_RNA > 200 & percent.mt <15 & nCount_RNA > 50)

pbmc <- NormalizeData(pbmc)
pbmc <- FindVariableFeatures(pbmc, nfeatures = 3000)
pbmc <- ScaleData(pbmc, vars.to.regress = "percent.mt")
pbmc <- RunPCA(pbmc)

pbmc <- FindNeighbors(pbmc, dims = 1:21)
pbmc <- FindClusters(pbmc, resolution = 1)
pbmc <- RunUMAP(pbmc, dims = 1:21)
DimPlot(pbmc, label = TRUE)

sweep_res <- paramSweep(pbmc,PCs = 1:20, sct = FALSE)
sweep_stat <- summarizeSweep(sweep_res, GT = F)
bcmv <- find.pK(sweep_stat)


ggplot(bcmv, aes(pK, BCmetric, group = 1)) +
  geom_point()+
  geom_line()

#0.005

pK<- bcmv %>%
  filter(BCmetric == max(BCmetric)) %>%
  select(pK)
pK<- as.numeric(as.character(pK[[1]]))

annotations <- pbmc@meta.data$seurat_clusters
homotipic <- modelHomotypic(annotations)
exp <- round(0.08*nrow(pbmc@meta.data))
exp.adj <- round(exp*(1-homotipic))

pbmc <- doubletFinder(pbmc,
                          PCs = 1:20,
                          pK = pK,
                          nExp = exp.adj,
                          reuse.pANN = F,
                          sct = F)
DimPlot(pbmc, reduction = "umap", group.by = "DF.classifications_0.25_0.11_818")

table(pbmc@meta.data$DF.classifications_0.25_0.11_818)
pbmc_filtered <- subset(pbmc, subset = DF.classifications_0.25_0.11_818 == "Singlet")

pbmc.markers <- FindAllMarkers(pbmc_filtered)

FeatureScatter(pbmc_filtered, feature1 = "CD3D", feature2 = "CD3E", slot = "data")

And here is my SCTransform workflow
library(Seurat)
library(tidyverse)
library(ggplot2)
library(gridExtra)
library(harmony)
library(future.apply)
library(cowplot)
library(patchwork)
library("DESeq2")
library(SeuratData)
library(SeuratDisk)
library(SeuratWrappers)
library(sctransform)
library(DoubletFinder)

data <- Read10X(data.dir = "/Users/rodrigohermoza/Desktop/UTEC/2024-0/PoliSia/rstudio/diferential_expression_2/DE2/filtered_gene_bc_matrices/raw_feature_bc_matrix/")
pbmc <- CreateSeuratObject(counts = data, project = "DE3")
pbmc <- PercentageFeatureSet(pbmc, pattern = "^MT-", col.name = "percent.mt") 
pbmc <- subset(pbmc, subset = nFeature_RNA > 200 & percent.mt <15 & nCount_RNA > 50)

pbmc <- SCTransform(pbmc, vars.to.regress = "percent.mt")
pbmc <- RunPCA(pbmc)
ElbowPlot(pbmc,ndims = 50)

pbmc <- FindNeighbors(pbmc, dims = 1:20)
pbmc <- FindClusters(pbmc, resolution = 1)
pbmc <- RunUMAP(pbmc, dims = 1:20)
DimPlot(pbmc, label = TRUE)

pbmc.markers <- FindAllMarkers(pbmc)

FeatureScatter(pbmc, feature1 = "CD3D", feature2 = "CD3E", slot = "data", smooth = T)

**Tags:** r, scrnaseq, seurat, rstudio
> **Comments:**
> **username_1**: Check in detail what is plotted. The "standard" is probably some continuous type of data, like scaled counts or something like this, and the second one is probably some sort of more transformed integers, which is not as tightly continuous as the "standard" one. Like imagine you do log2 of 1:10, that gives defined breaks, while something like a Z-score gives a much more continuous distribution without "hard breaks".  
> **username_0**: Thanks for the clarification :)  

---

## 0 answers

