## ----include = FALSE---------------------------------------------------------- knitr::opts_chunk$set( collapse = TRUE, comment = "#>" ) ## ----eval = FALSE------------------------------------------------------------- # library(remotes) # remotes::install_github("ilia-kats/MuData") ## ----setup, message = FALSE--------------------------------------------------- library(MuData) library(SingleCellExperiment) library(MultiAssayExperiment) library(SingleCellMultiModal) library(scater) library(rhdf5) ## ----------------------------------------------------------------------------- mae <- CITEseq( DataType="cord_blood", modes="*", dry.run=FALSE, version="1.0.0" ) mae ## ----------------------------------------------------------------------------- # Define CLR transformation as in the Seurat workflow clr <- function(data) t( apply(data, 1, function(x) log1p( x / (exp(sum(log1p(x[x > 0]), na.rm = TRUE) / length(x))) )) ) ## ----------------------------------------------------------------------------- adt_counts <- mae[["scADT"]] mae[["scADT"]] <- SingleCellExperiment(adt_counts) assay(mae[["scADT"]], "clr") <- clr(adt_counts) ## ----------------------------------------------------------------------------- mae[["scADT"]] <- runPCA( mae[["scADT"]], exprs_values = "clr", ncomponents = 20 ) ## ----------------------------------------------------------------------------- plotReducedDim(mae[["scADT"]], dimred = "PCA", by_exprs_values = "clr", colour_by = "CD3") plotReducedDim(mae[["scADT"]], dimred = "PCA", by_exprs_values = "clr", colour_by = "CD14") ## ----------------------------------------------------------------------------- writeH5MU(mae, "cord_blood_citeseq.h5mu") ## ----------------------------------------------------------------------------- h5 <- rhdf5::H5Fopen("cord_blood_citeseq.h5mu") h5ls(H5Gopen(h5, "mod"), recursive = FALSE) ## ----------------------------------------------------------------------------- h5ls(H5Gopen(h5, "mod/scADT"), recursive = FALSE) h5ls(H5Gopen(h5, "mod/scADT/layers"), recursive = FALSE) ## ----------------------------------------------------------------------------- h5ls(H5Gopen(h5, "mod/scADT/obsm"), recursive = FALSE) # There is an alternative way to access groups: # h5&'mod'&'scADT'&'obsm' rhdf5::H5close() ## ----------------------------------------------------------------------------- sessionInfo()