Contents

1 Description

This vignette demonstrates how MEFISTO can be used to integrate single-cell multi-omics data sets. Here we used the scNMT-seq gastrulation data set, which consists of ~1,500 cells collected across four stages of mouse development, and can be downloaded from here. For roughly half of the cells we profiled RNA expression, DNA methylation and chromatin accessibility from the same cell, and for the other half only RNA expression was profiled. Thus, a significant challenge of this data set is the sparsity of the epigenetic readouts.

You may wonder, how can we use MEFISTO here if we don’t have continuous spatial or temporal covariates? Here we’ll use the RNA expression to derive a multivariate developmental trajectory, which we can then take as a proxy for real time, and use this manifold as a covariate in MEFISTO. This approach wil allow us to identify coordinated patterns of variation between the transcriptome and epigenome during development.

2 Load libraries

library(data.table)
library(purrr)
library(ggplot2)
library(ggpubr)
library(ggrepel)

# MOFA
library(MOFA2)

3 Load data

Data is available for download here

The input data consists of three modalities: RNA expression quantified over genes; and DNA methylation and chromatin accessibility quantified over TF motifs. The script to do this step can be found here. Then, it is just a matter of doing some feature selection and preparing the input data for MEFISTO. The script for this step can be found here. For the purpose of this vignette we will load the multi-omics data already processed in a long data.table format

data <- fread("/Users/ricard/data/gastrulation/metaccrna/mefisto/vignette/scnmt_data.txt.gz")
head(data,n=3)
##            sample feature   value view
## 1: E7.5_Plate1_C9   Sox17  0.5903  RNA
## 2: E7.5_Plate1_C9 Atp6v1h -0.3542  RNA
## 3: E7.5_Plate1_C9   Sulf1 -1.4312  RNA

Print number of features per view

data[,length(unique(feature)),by="view"]
##         view  V1
## 1:       RNA 887
## 2: motif_met 500
## 3: motif_acc 500

Load cell metadata

sample_metadata <- fread("/Users/ricard/data/gastrulation/metaccrna/mefisto/vignette/scnmt_sample_metadata.txt")
colnames(sample_metadata)
##  [1] "sample"       "embryo"       "stage"        "pass_rnaQC"   "pass_metQC"  
##  [6] "pass_accQC"   "lineage10x"   "lineage10x_2" "UMAP1"        "UMAP2"

There are two cell type columns: lineage10x and lineage10x_2. The first one is obtained by mapping the scNMT-seq cells to the RNA gastrulation atlas. The problem with the most this annotation is that we have very few cells per cell type. In the second column these lineages are aggregated into bigger classes

table(sample_metadata$lineage10x)
## 
##      Anterior_Primitive_Streak                Caudal_epiblast 
##                              4                             96 
##                Caudal_mesoderm                  Def._endoderm 
##                              5                             44 
##                       Epiblast                   ExE_mesoderm 
##                            548                              1 
##                            Gut Haematoendothelial_progenitors 
##                             38                              5 
##          Intermediate_mesoderm                     Mesenchyme 
##                             23                             27 
##                 Mixed_mesoderm               Nascent_mesoderm 
##                            121                            222 
##                      Notochord              Paraxial_mesoderm 
##                             49                             26 
##            Pharyngeal_mesoderm               Primitive_Streak 
##                             17                             42 
##           Rostral_neurectoderm               Somitic_mesoderm 
##                            181                             34 
##               Surface_ectoderm              Visceral_endoderm 
##                              9                             26
table(sample_metadata$lineage10x_2)
## 
##         Ectoderm         Endoderm         Epiblast         Mesoderm 
##              188              157              550              481 
## Primitive_Streak 
##              142

For defining the multivariate pseudotime trajectory we applied UMAP on the RNA expression data. Due to the randomness of the algorithm, there is flexibility on how to calculate the latent manifold. For reproducibility purposes we included our precomputed UMAP coordinates in the samples metadata columns UMAP1 and UMAP2

ggscatter(sample_metadata, x="UMAP1", y="UMAP2", color="lineage10x_2") +
  scale_color_manual(values=celltype.colors) +
  theme_classic() +
  ggplot_theme_NoAxes()