CRUK CI Summer School 2020

Pseudotime Analysis

library(SingleCellExperiment)
library(scran)
library(scater)
library(batchelor)
library(cowplot)
library(pheatmap)
library(tidyverse)
library(SingleR)
library(destiny)
library(gam)
library(viridis)
library(msigdbr)
library(clusterProfiler)
library(cellAlign)

1 Pseudotime Alignment

CellAlign is a tool for quantitative comparison of expression dynamics within or between single-cell trajectories. The input to the CellAlign workflow is any trajectory vector that orders single cell expression with a pseudo-time spacing and the expression matrix for the cells used to define the trajectory. cellAlign has 3 essential steps:

1. Interpolate the data to have N evenly spaced points along the scaled pseudotime vector using a sliding window of Gaussian weights

2. Determine the genes of interest for alignment

3. Align your trajectory among the selected genes either along the whole trajectory or along a partial segment.

The first step is to interpolate the data along the trajectory to represent the data by N (default 200) equally spaced points along the pseudotime trajectory. We included this step because single-cell measurements are often sparse or heterogeneous along the trajectory, leaving gaps that cannot be aligned. Cell-Align interpolates the gene-expression values of equally spaced artificial points using the real single-cell expression data. The expression values of the interpolated points are calculated using all cells, with each single cell assigned a weight given by a Gaussian distribution centered at the interpolated point and a width assigned by a parameter called winSz. The default winSz is 0.1, as this is the range that preserves the dynamics of the trajectory without including excessive noise for standard single cell data sets.

interGlobal_caronPRET1 <- cellAlign::interWeights(expDataBatch = t(caron.PRET1_counts), 
                                                    trajCond = eigenvectors(dm_caron.PRET1)[, 1], 
                                                    winSz = 0.1, numPts=200)

interGlobal_hcaBM1 <- cellAlign::interWeights(expDataBatch = t(tcell_BM1_counts), 
                                                    trajCond = eigenvectors(dm_tcell_BM1)[, 1], 
                                                    winSz = 0.1, numPts=200)

interGlobal_hcaBM2 <- cellAlign::interWeights(expDataBatch = t(tcell_counts_BM2), 
                                                    trajCond = eigenvectors(dm_tcell_BM2)[, 1], 
                                                    winSz = 0.1, numPts=200)

Scale the expression matrix

interGlobal_caronPRET1_scaled = scaleInterpolate(interGlobal_caronPRET1)
interGlobal_hcaBM1_scaled = cellAlign::scaleInterpolate(interGlobal_hcaBM1)
interGlobal_hcaBM2_scaled = cellAlign::scaleInterpolate(interGlobal_hcaBM2)

Identify the shared genes across datasets

sharedMarkers = Reduce(intersect, list(rownames(interGlobal_caronPRET1$interpolatedVals),rownames(interGlobal_hcaBM1$interpolatedVals),rownames(interGlobal_hcaBM2$interpolatedVals)))
length(sharedMarkers)

[1] 107

Finally, there is the alignment step. CellAlign operates much like sequence alignment algorithms, quantifying overall similarity in expression throughout the trajectory (global alignment), or finding areas of highly conserved expression (local alignment). Cell-Align then finds a path through the matrix that minimizes the overall distance while adhering to the following constraints: * for global alignment the alignment must cover the entire extent of both trajectories, always starting in the upper left of the dissimilarity matrix and ending in the lower right.

• for local alignment the alignment is restricted only to highly similar cells, yielding as output regions with conserved expression dynamics

Intuitively, the optimal alignment runs along a “valley” within the dissimilarity matrix.

A=calcDistMat(interGlobal_caronPRET1_scaled$scaledData[sharedMarkers,],interGlobal_hcaBM1_scaled$scaledData[sharedMarkers,], dist.method = 'Euclidean')
alignment = globalAlign(A)
plotAlign(alignment)

[1] NA

B=calcDistMat(interGlobal_hcaBM1_scaled$scaledData[sharedMarkers,],interGlobal_hcaBM2_scaled$scaledData[sharedMarkers,], dist.method = 'Euclidean')
alignment = globalAlign(B)
plotAlign(alignment)

[1] NA

2 Ackowledgements

This notebook uses material from cellAlign vignette.

sessionInfo()

R version 4.0.2 (2020-06-22)
Platform: x86_64-pc-linux-gnu (64-bit)
Running under: Ubuntu 18.04.4 LTS

Matrix products: default
BLAS:   /usr/lib/x86_64-linux-gnu/blas/libblas.so.3.7.1
LAPACK: /usr/lib/x86_64-linux-gnu/lapack/liblapack.so.3.7.1

locale:
 [1] LC_CTYPE=en_US.UTF-8 LC_NUMERIC=C         LC_TIME=C            LC_COLLATE=C         LC_MONETARY=C        LC_MESSAGES=C       
 [7] LC_PAPER=C           LC_NAME=C            LC_ADDRESS=C         LC_TELEPHONE=C       LC_MEASUREMENT=C     LC_IDENTIFICATION=C 

attached base packages:
 [1] splines   parallel  stats4    stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
 [1] scRNAseq_2.2.0              knitr_1.29                  cellAlign_0.1.0             clusterProfiler_3.16.0     
 [5] msigdbr_7.1.1               viridis_0.5.1               viridisLite_0.3.0           gam_1.20                   
 [9] foreach_1.5.0               destiny_3.2.0               SingleR_1.2.4               forcats_0.5.0              
[13] stringr_1.4.0               dplyr_1.0.0                 purrr_0.3.4                 readr_1.3.1                
[17] tidyr_1.1.0                 tibble_3.0.3                tidyverse_1.3.0             pheatmap_1.0.12            
[21] cowplot_1.0.0               batchelor_1.4.0             scater_1.16.2               ggplot2_3.3.2              
[25] scran_1.16.0                SingleCellExperiment_1.10.1 SummarizedExperiment_1.18.2 DelayedArray_0.14.1        
[29] matrixStats_0.56.0          Biobase_2.48.0              GenomicRanges_1.40.0        GenomeInfoDb_1.24.2        
[33] IRanges_2.22.2              S4Vectors_0.26.1            BiocGenerics_0.34.0        

loaded via a namespace (and not attached):
  [1] rappdirs_0.3.1                ggthemes_4.2.0                bit64_0.9-7.1                 irlba_2.3.3                  
  [5] data.table_1.12.8             RCurl_1.98-1.2                generics_0.0.2                RSQLite_2.2.0                
  [9] europepmc_0.4                 proxy_0.4-24                  bit_1.1-15.2                  enrichplot_1.8.1             
 [13] xml2_1.3.2                    lubridate_1.7.9               httpuv_1.5.4                  assertthat_0.2.1             
 [17] xfun_0.15                     hms_0.5.3                     evaluate_0.14                 promises_1.1.1               
 [21] DEoptimR_1.0-8                fansi_0.4.1                   progress_1.2.2                dbplyr_1.4.4                 
 [25] readxl_1.3.1                  igraph_1.2.5                  DBI_1.1.0                     ellipsis_0.3.1               
 [29] RSpectra_0.16-0               backports_1.1.8               bookdown_0.20                 vctrs_0.3.2                  
 [33] TTR_0.23-6                    abind_1.4-5                   RcppEigen_0.3.3.7.0           withr_2.2.0                  
 [37] ggforce_0.3.2                 packrat_0.5.0                 triebeard_0.3.0               robustbase_0.93-6            
 [41] vcd_1.4-7                     xts_0.12-0                    prettyunits_1.1.1             DOSE_3.14.0                  
 [45] ExperimentHub_1.14.0          laeken_0.5.1                  crayon_1.3.4                  labeling_0.3                 
 [49] edgeR_3.30.3                  pkgconfig_2.0.3               tweenr_1.0.1                  nlme_3.1-148                 
 [53] vipor_0.4.5                   nnet_7.3-14                   rlang_0.4.7                   lifecycle_0.2.0              
 [57] downloader_0.4                BiocFileCache_1.12.0          modelr_0.1.8                  rsvd_1.0.3                   
 [61] AnnotationHub_2.20.0          cellranger_1.1.0              polyclip_1.10-0               RcppHNSW_0.2.0               
 [65] lmtest_0.9-37                 Matrix_1.2-18                 urltools_1.7.3                carData_3.0-4                
 [69] boot_1.3-25                   zoo_1.8-8                     reprex_0.3.0                  base64enc_0.1-3              
 [73] beeswarm_0.2.3                ggridges_0.5.2                knn.covertree_1.0             bitops_1.0-6                 
 [77] blob_1.2.1                    DelayedMatrixStats_1.10.1     qvalue_2.20.0                 gridGraphics_0.5-0           
 [81] scales_1.1.1                  memoise_1.1.0                 magrittr_1.5                  plyr_1.8.6                   
 [85] hexbin_1.28.1                 zlibbioc_1.34.0               compiler_4.0.2                scatterpie_0.1.4             
 [89] dqrng_0.2.1                   RColorBrewer_1.1-2            pcaMethods_1.80.0             cli_2.0.2                    
 [93] dtw_1.21-3                    XVector_0.28.0                mgcv_1.8-31                   ggplot.multistats_1.0.0      
 [97] MASS_7.3-51.6                 tidyselect_1.1.0              stringi_1.4.6                 yaml_2.2.1                   
[101] GOSemSim_2.14.0               BiocSingular_1.4.0            locfit_1.5-9.4                ggrepel_0.8.2                
[105] grid_4.0.2                    fastmatch_1.1-0               tools_4.0.2                   rio_0.5.16                   
[109] rstudioapi_0.11               foreign_0.8-80                gridExtra_2.3                 smoother_1.1                 
[113] scatterplot3d_0.3-41          farver_2.0.3                  ggraph_2.0.3                  digest_0.6.25                
[117] rvcheck_0.1.8                 BiocManager_1.30.10           shiny_1.5.0                   Rcpp_1.0.5                   
[121] car_3.0-8                     broom_0.7.0                   BiocVersion_3.11.1            later_1.1.0.1                
[125] httr_1.4.1                    AnnotationDbi_1.50.1          colorspace_1.4-1              rvest_0.3.5                  
[129] fs_1.4.2                      ranger_0.12.1                 statmod_1.4.34                graphlayouts_0.7.0           
[133] sp_1.4-2                      ggplotify_0.0.5               xtable_1.8-4                  jsonlite_1.7.0               
[137] tidygraph_1.2.0               R6_2.4.1                      pillar_1.4.6                  htmltools_0.5.0              
[141] mime_0.9                      glue_1.4.1                    fastmap_1.0.1                 VIM_6.0.0                    
[145] BiocParallel_1.22.0           BiocNeighbors_1.6.0           class_7.3-17                  interactiveDisplayBase_1.26.3
[149] codetools_0.2-16              fgsea_1.14.0                  lattice_0.20-41               curl_4.3                     
[153] ggbeeswarm_0.6.0              zip_2.0.4                     GO.db_3.11.4                  openxlsx_4.1.5               
[157] limma_3.44.3                  rmarkdown_2.3                 munsell_0.5.0                 e1071_1.7-3                  
[161] DO.db_2.9                     GenomeInfoDbData_1.2.3        iterators_1.0.12              haven_2.3.1                  
[165] reshape2_1.4.4                gtable_0.3.0