diff --git a/analysis/MOFAaIg.Rmd b/analysis/MOFAaIg.Rmd index bdba9c6..31a35a5 100644 --- a/analysis/MOFAaIg.Rmd +++ b/analysis/MOFAaIg.Rmd @@ -9,6 +9,11 @@ editor_options: --- ```{r setup, warning=FALSE, message=FALSE} + +knitr::opts_chunk$set( + message = F, warning = F, echo = T, eval = T +) + source("code/load_packages.R") seu_combined_selectsamples <- readRDS("output/seu_aIG_samples.rds") library(ggrepel) @@ -48,6 +53,8 @@ meta.allcells <- seu_combined_selectsamples@meta.data %>% #### MOFA model +Note, if mofa object was generated before, it will read the generated rds file. (this will speed-up the process of generating this html document if edits are required) + ```{r Mofa object} ## Note, if mofa object was generated before, it will read the generated rds file. (this will speed-up the process of generating this html document if edits are required) @@ -84,9 +91,9 @@ mofa ``` ```{r} -## Rename protein names -featurenamesmofa <- features_names(mofa) +## Hack to rename protein names for visualization +featurenamesmofa <- features_names(mofa) ## todo cleanup/more efficient featurenamesmofa <- rapply(featurenamesmofa,function(x) ifelse(x=="p-Src" ,"p-SRC",x), how = "replace") @@ -102,6 +109,12 @@ featurenamesmofa <- rapply(featurenamesmofa,function(x) ifelse(x=="p-Myc" ,"p-MY featurenamesmofa <- rapply(featurenamesmofa,function(x) ifelse(x=="p-Btk" ,"p-BTK",x), how = "replace") featurenamesmofa <- rapply(featurenamesmofa,function(x) ifelse(x=="Bcl-6" ,"BCL6",x), how = "replace") +featurenamesmofa <- rapply(featurenamesmofa,function(x) ifelse(x=="CD53-PROT" ,"CD53",x), how = "replace") +featurenamesmofa <- rapply(featurenamesmofa,function(x) ifelse(x=="CD70-PROT" ,"CD70",x), how = "replace") +featurenamesmofa <- rapply(featurenamesmofa,function(x) ifelse(x=="KLF6-PROT" ,"KLF6",x), how = "replace") +featurenamesmofa <- rapply(featurenamesmofa,function(x) ifelse(x=="XBP1-PROT" ,"XBP1",x), how = "replace") + + features_names(mofa) <- featurenamesmofa @@ -281,7 +294,6 @@ WNN.Dimplot <- DimPlot(seu_combined_selectsamples, reduction = 'wnn.umap', label ```{r} Factorvalues <- data.frame(get_factors(mofa)[1])[,1:4] colnames(Factorvalues) <- c(paste0("Factor", 1:4)) -Factorvalues seu_combined_selectsamples <- AddMetaData(seu_combined_selectsamples, Factorvalues) plot.wnn.Factor1 <- FeaturePlot(seu_combined_selectsamples, reduction = "wnn.umap", features = "Factor1") + diff --git a/analysis/QC.Rmd b/analysis/QC.Rmd index a0bb683..a39a720 100644 --- a/analysis/QC.Rmd +++ b/analysis/QC.Rmd @@ -148,6 +148,18 @@ seu_combined_filtered <- subset(seu_combined, subset = nFeature_RNA > 150 & nCou ``` +```{r} +## Small hack to prevent mofa error with duplicate gene or protein names +double <- c("CD53", "CD70", "KLF6", "XBP1") +PROT.counts <-as.data.frame(seu_combined_filtered[["PROT"]]@counts) +prot.rownames <- rownames(PROT.counts) +row.names(PROT.counts) <- ifelse(prot.rownames %in% double, paste0(prot.rownames, "-PROT"), paste0(prot.rownames)) + +seu_combined_filtered[["PROT"]] <- CreateAssayObject(counts = PROT.counts) + +``` + + ## Normalize and scale