From d9ad8cd3427990d8601d3b9f82c377dba8a95df2 Mon Sep 17 00:00:00 2001
From: Davis McCarthy <davismcc@gmail.com>
Date: Tue, 1 Oct 2019 18:15:07 +1000
Subject: [PATCH] Bug fix to remove-conf.Rmd
---
course_files/exprs-norm.Rmd | 12 ++++++------
course_files/remove-conf.Rmd | 8 ++++----
2 files changed, 10 insertions(+), 10 deletions(-)
diff --git a/course_files/exprs-norm.Rmd b/course_files/exprs-norm.Rmd
index ce319c4..c4df548 100644
--- a/course_files/exprs-norm.Rmd
+++ b/course_files/exprs-norm.Rmd
@@ -3,7 +3,7 @@ output: html_document
---
```{r setup, echo=FALSE}
-knitr::opts_chunk$set(out.width='90%', fig.align = 'center', eval=TRUE)
+knitr::opts_chunk$set(out.width='90%', fig.align = 'center', eval=TRUE, warning=FALSE, message=FALSE)
knitr::opts_knit$set(root.dir = normalizePath(".."))
```
@@ -459,7 +459,7 @@ Note that (due to what looks like a bug in this version of `sctransform`) we
need to convert the UMI count matrix to a sparse format to apply sctransform.
-```{r sctransform-apply}
+```{r sctransform-apply, warning=FALSE, message=FALSE}
umi_sparse <- as(counts(umi.qc), "dgCMatrix")
### Genes expressed in at least 5 cells will be kept
sctnorm_data <- sctransform::vst(umi = umi_sparse, min_cells = 1,
@@ -475,7 +475,7 @@ assay(umi.qc, "sctrans_norm") <- sctnorm_data$y
Let us look at the NB GLM model parameters estimated by sctransform.
-```{r sctransform-params-plot}
+```{r sctransform-params-plot, warning=FALSE, message=FALSE}
#sce$log10_total_counts
## Matrix of estimated model parameters per gene (theta and regression coefficients)
sctransform::plot_model_pars(sctnorm_data)
@@ -484,14 +484,14 @@ sctransform::plot_model_pars(sctnorm_data)
We can look at the effect of sctransform's normalization on three particular
genes, ACTB, POU5F1 (aka OCT4) and CD74.
-```{r sctransform-genes-plot}
+```{r sctransform-genes-plot, warning=FALSE, message=FALSE}
##c('ACTB', 'Rpl10', 'Cd74')
genes_plot <- c("ENSG00000075624", "ENSG00000204531", "ENSG00000019582")
sctransform::plot_model(sctnorm_data, umi_sparse, genes_plot,
plot_residual = TRUE, cell_attr = as.data.frame(colData(umi.qc)))
```
-```{r norm-pca-sctransform, fig.cap = "PCA plot of the tung data after sctransform normalisation (Pearson residuals)."}
+```{r norm-pca-sctransform, warning=FALSE, message=FALSE, fig.cap = "PCA plot of the tung data after sctransform normalisation (Pearson residuals)."}
reducedDim(umi.qc, "PCA_sctrans_norm") <- reducedDim(
runPCA(umi.qc[endog_genes, ], exprs_values = "sctrans_norm")
)
@@ -504,7 +504,7 @@ plotReducedDim(
) + ggtitle("PCA plot: sctransform normalization")
```
-```{r norm-ours-rle-sctransform, fig.cap = "Cell-wise RLE of the tung data"}
+```{r norm-ours-rle-sctransform, warning=FALSE, message=FALSE, fig.cap = "Cell-wise RLE of the tung data"}
plotRLE(
umi.qc[endog_genes, ],
exprs_values = "sctrans_norm",
diff --git a/course_files/remove-conf.Rmd b/course_files/remove-conf.Rmd
index 69eccad..9486212 100644
--- a/course_files/remove-conf.Rmd
+++ b/course_files/remove-conf.Rmd
@@ -7,7 +7,7 @@ knitr::opts_chunk$set(out.width='90%', fig.align = 'center', eval=TRUE)
knitr::opts_knit$set(root.dir = normalizePath(".."))
```
-## Dealing with confounders
+## Batch effects
### Introduction
@@ -165,19 +165,19 @@ matrix in the `lm_batch_indi` slot.
umi.qc$cdr <- umi.qc$total_features_by_counts_endogenous / nrow(umi.qc)
## fit a model just accounting for batch by individual
lm_design_batch1 <- model.matrix(~batch + cdr,
- data = coldata(umi.qc)[umi.qc$individual == "na19098",])
+ data = colData(umi.qc)[umi.qc$individual == "na19098",])
fit_indi1 <- lmfit(logcounts(umi.qc)[, umi.qc$individual == "na19098"], lm_design_batch1)
fit_indi1$coefficients[,1] <- 0 ## replace intercept with 0 to preserve reference batch
resids_lm_batch1 <- residuals(fit_indi1, logcounts(umi.qc)[, umi.qc$individual == "na19098"])
lm_design_batch2 <- model.matrix(~batch + cdr,
- data = coldata(umi.qc)[umi.qc$individual == "na19101",])
+ data = colData(umi.qc)[umi.qc$individual == "na19101",])
fit_indi2 <- lmfit(logcounts(umi.qc)[, umi.qc$individual == "na19101"], lm_design_batch2)
fit_indi2$coefficients[,1] <- 0 ## replace intercept with 0 to preserve reference batch
resids_lm_batch2 <- residuals(fit_indi2, logcounts(umi.qc)[, umi.qc$individual == "na19101"])
lm_design_batch3 <- model.matrix(~batch + cdr,
- data = coldata(umi.qc)[umi.qc$individual == "na19239",])
+ data = colData(umi.qc)[umi.qc$individual == "na19239",])
fit_indi3 <- lmfit(logcounts(umi.qc)[, umi.qc$individual == "na19239"], lm_design_batch3)
fit_indi3$coefficients[,1] <- 0 ## replace intercept with 0 to preserve reference batch
resids_lm_batch3 <- residuals(fit_indi3, logcounts(umi.qc)[, umi.qc$individual == "na19239"])
--
GitLab