16.1 Myeloid Trajectory - Monocyte Differentiation

In this section, we will create a cellular trajectory that approximates the differentiation of HSCs into fully differentiated monocytes. To start, lets review the clusters and cell types that we defined previously, stored in cellColData in columns named “Clusters” and “Clusters2”. Overlaying these cell groupings on our UMAP embedding shows the different cell types that we are interested in.

p1 <- plotEmbedding(ArchRProj = projHeme5, colorBy = "cellColData", name = "Clusters", embedding = "UMAP")

## ArchR logging to : ArchRLogs/ArchR-plotEmbedding-15e1326ab29a8-Date-2020-04-15_Time-14-59-47.log
## If there is an issue, please report to github with logFile!
## Getting UMAP Embedding
## ColorBy = cellColData
## Plotting Embedding
## 1
## ArchR logging successful to : ArchRLogs/ArchR-plotEmbedding-15e1326ab29a8-Date-2020-04-15_Time-14-59-47.log

p2 <- plotEmbedding(ArchRProj = projHeme5, colorBy = "cellColData", name = "Clusters2", embedding = "UMAP")

## ArchR logging to : ArchRLogs/ArchR-plotEmbedding-15e1339bfbfcf-Date-2020-04-15_Time-14-59-49.log
## If there is an issue, please report to github with logFile!
## Getting UMAP Embedding
## ColorBy = cellColData
## Plotting Embedding
## 1
## ArchR logging successful to : ArchRLogs/ArchR-plotEmbedding-15e1339bfbfcf-Date-2020-04-15_Time-14-59-49.log

ggAlignPlots(p1, p2, type = "h")

16.1.1 Pseudo-time UMAPs and individual feature plots

We will use the cell type definitions that are stored in “Clusters2”. As mentioned above, we are creating a trajectory from stem cells (“Progenitor”), through the commited myeloid progenitor (“GMP”), to monocytes (“Mono”). The first step in creating a trajectory is to create a trajectory backbone in the form of an ordered vector of cell group labels.

trajectory <- c("Progenitor", "GMP", "Mono")
trajectory

## [1] “Progenitor” “GMP” “Mono”

We use the addTrajectory() function to create a trajectory and we add this to our ArchRProject. We will call this trajectory “MyeloidU”. What this does is it creates a new column in cellColData called “MyeloidU” that stores the pseudo-time value for each cell in the trajectory. Cells that are not part of the trajectory are labeled with NA.

projHeme5 <- addTrajectory(
    ArchRProj = projHeme5, 
    name = "MyeloidU", 
    groupBy = "Clusters2",
    trajectory = trajectory, 
    embedding = "UMAP", 
    force = TRUE
)

## ArchR logging to : ArchRLogs/ArchR-addTrajectory-15e132f14b12b-Date-2020-04-15_Time-15-00-01.log
## If there is an issue, please report to github with logFile!
## Filtering outliers
## Initial Alignment Before Spline Fit
## Spline Fit
## KNN to Spline
## Aligning cells not in trajectory
## ArchR logging successful to : ArchRLogs/ArchR-addTrajectory-15e132f14b12b-Date-2020-04-15_Time-15-00-01.log

We can look at this information and see that each cell has a unique pseudotime value between 0 and 100. We exclude cells with NA values because these are not part of the trajectory.

head(projHeme5$MyeloidU[!is.na(projHeme5$MyeloidU)])

## [1] 46.07479 53.75027 44.82834 43.18828 47.49617 43.21015

To plot this trajectory, we use the plotTrajectory() function which overlays the pseudo-time values on our UMAP embedding and displays an arrow approximating the trajectory path from the spline-fit. Cells that are not part of the trajectory are colored gray. In this example, we use colorBy = "cellColData" to tell ArchR to look within cellColData for the column specified by name - in this case, the “MyeloidU” pseudo-time trajectory. While it may seem conterintuitive to list “MyeloidU” for both trajectory and name, this is done because trajectory tells ArchR which subset of cells we are interested in and name tells ArchR how to color that subset of cells.

p <- plotTrajectory(projHeme5, trajectory = "MyeloidU", colorBy = "cellColData", name = "MyeloidU")

## ArchR logging to : ArchRLogs/ArchR-plotTrajectory-15e136ad2f619-Date-2020-04-15_Time-15-00-03.log
## If there is an issue, please report to github with logFile!
## Plotting
## Warning: Removed 5678 rows containing non-finite values (stat_summary_hex).
## Plotting Trajectory
## Adding Inferred Arrow Trajectory to Plot
## ArchR logging successful to : ArchRLogs/ArchR-plotTrajectory-15e136ad2f619-Date-2020-04-15_Time-15-00-03.log

p[[1]]

## Warning: Removed 5678 rows containing non-finite values (stat_summary_hex).

To save an editable vectorized version of this plot, we use plotPDF().

plotPDF(p, name = "Plot-MyeloidU-Traj-UMAP.pdf", ArchRProj = projHeme5, addDOC = FALSE, width = 5, height = 5)

## [1] “plotting ggplot!”
## Warning: Removed 5678 rows containing non-finite values (stat_summary_hex).
## Warning: Removed 5678 rows containing non-finite values (stat_summary_hex).
## [1] “plotting ggplot!”
## geom_smooth() using method = ‘gam’ and formula ‘y ~ s(x, bs = “cs”)’
## geom_smooth() using method = ‘gam’ and formula ‘y ~ s(x, bs = “cs”)’
## [1] 0

We can also overlay other features on the trajectory within our UMAP embedding. This allows us to display specific features only within the cells that are relevant to our trajectory.

If you have not already added impute weights to your projHeme5 project, lets do that now.

projHeme5 <- addImputeWeights(projHeme5)

## 2020-04-15 15:00:28 : Computing Impute Weights Using Magic (Cell 2018), 0 mins elapsed.
## 2020-04-15 15:00:40 : Completed Getting Magic Weights!, 0.206 mins elapsed.

Then, we can plot the “MyeloidU” trajectory but color the cells by the gene score value of the CEBPB gene, a known regulator of monocyte function that becomes active during the differentiation process. We indicate the matrix to use via the colorBy parameter and the feature to use via the name parameter.

p1 <- plotTrajectory(projHeme5, trajectory = "MyeloidU", colorBy = "GeneScoreMatrix", name = "CEBPB", continuousSet = "horizonExtra")

## Getting ImputeWeights
## ArchR logging to : ArchRLogs/ArchR-plotTrajectory-15e13179524ae-Date-2020-04-15_Time-15-00-40.log
## If there is an issue, please report to github with logFile!
## Getting Matrix Values…
## Getting Matrix Values…
##
## Imputing Matrix
## Using weights on disk
## Using weights on disk
## Plotting
## Warning: Removed 5678 rows containing non-finite values (stat_summary_hex).
## Plotting Trajectory
## Adding Inferred Arrow Trajectory to Plot
## ArchR logging successful to : ArchRLogs/ArchR-plotTrajectory-15e13179524ae-Date-2020-04-15_Time-15-00-40.log

We can repeat this process but color the cells by their linked gene expression via the GeneIntegrationMatrix.

p2 <- plotTrajectory(projHeme5, trajectory = "MyeloidU", colorBy = "GeneIntegrationMatrix", name = "CEBPB", continuousSet = "blueYellow")

## Getting ImputeWeights
## ArchR logging to : ArchRLogs/ArchR-plotTrajectory-15e133644196a-Date-2020-04-15_Time-15-00-47.log
## If there is an issue, please report to github with logFile!
## Getting Matrix Values…
## Getting Matrix Values…
##
## Imputing Matrix
## Using weights on disk
## Using weights on disk
## Plotting
## Warning: Removed 5678 rows containing non-finite values (stat_summary_hex).
## Plotting Trajectory
## Adding Inferred Arrow Trajectory to Plot
## ArchR logging successful to : ArchRLogs/ArchR-plotTrajectory-15e133644196a-Date-2020-04-15_Time-15-00-47.log

The plotTrajectory() function actually returns a list of relevant plots. The first plot in the list is a UMAP embedding, colorized as specified in the function call.

Comparing these UMAP plots side-by-side for the gene score and gene expression, we see that the activity of the CEBPB gene is highly specific to monocyte cells in the later part of the pseudo-time trajectory.

ggAlignPlots(p1[[1]], p2[[1]], type = "h")

## Warning: Removed 5678 rows containing non-finite values (stat_summary_hex).
## Warning: Removed 5678 rows containing non-finite values (stat_summary_hex).

The second plot that is returned by plotTrajectory() is a dot plot of pseudo-time versus the value of the relevant feature, in this case, the gene score or gene expression of CEBPB. In this case, the cells are colored by their pseudo-time.

ggAlignPlots(p1[[2]], p2[[2]], type = "h")

## geom_smooth() using method = ‘gam’ and formula ‘y ~ s(x, bs = “cs”)’
## geom_smooth() using method = ‘gam’ and formula ‘y ~ s(x, bs = “cs”)’

16.1.2 Pseudo-time heatmaps

We can visualize changes in many features across pseudo-time using heatmaps. To do this, we first retrieve the trajectory of interest from the ArchRProject using the getTrajectory() function which returns the trajectory as a SummarizedExperiment object. We will create these pseudo-time heatmaps for motifs, gene scores, gene expression, and peak accessibility by passing the corresponding matrix to the useMatrix parameter.

trajMM  <- getTrajectory(ArchRProj = projHeme5, name = "MyeloidU", useMatrix = "MotifMatrix", log2Norm = FALSE)

## Creating Trajectory Group Matrix..
## Some values are below 0, this could be a DeviationsMatrix in which scaleTo should be set = NULL.
## Continuing without depth normalization!
## Smoothing…

We then pass this SummarizedExperiment to the plotTrajectoryHeatmap() function.

p1 <- plotTrajectoryHeatmap(trajMM, pal = paletteContinuous(set = "solarExtra"))

## ArchR logging to : ArchRLogs/ArchR-plotTrajectoryHeatmap-15e131c8ea2de-Date-2020-04-15_Time-15-01-31.log
## If there is an issue, please report to github with logFile!
## useSeqnames is NULL or greater than 1 with a Sparse.Assays.Matrix trajectory input.
## 2020-04-15 15:01:31 :
## force=FALSE thus continuing with subsetting useSeqnames = z
## 2020-04-15 15:01:31 :
## Removing rows with NA values…
## 2020-04-15 15:01:31 :
## Preparing Main Heatmap..
## ArchR logging successful to : ArchRLogs/ArchR-plotTrajectoryHeatmap-15e131c8ea2de-Date-2020-04-15_Time-15-01-31.log

We can perform the same steps to obtain a pseudo-time heatmap of gene scores by setting useMatrix = "GeneScoreMatrix".

trajGSM <- getTrajectory(ArchRProj = projHeme5, name = "MyeloidU", useMatrix = "GeneScoreMatrix", log2Norm = TRUE)

## Creating Trajectory Group Matrix..
## Smoothing…

p2 <- trajectoryHeatmap(trajGSM,  pal = paletteContinuous(set = "horizonExtra"))

## ArchR logging to : ArchRLogs/ArchR-trajectoryHeatmap-15e132f1f4173-Date-2020-04-15_Time-15-01-51.log
## If there is an issue, please report to github with logFile!
## Preparing Main Heatmap..
## ArchR logging successful to : ArchRLogs/ArchR-trajectoryHeatmap-15e132f1f4173-Date-2020-04-15_Time-15-01-51.log

Similarly, we can obtain a pseudo-time heatmap of gene expression by setting useMatrix = "GeneIntegrationMatrix".

trajGIM <- getTrajectory(ArchRProj = projHeme5, name = "MyeloidU", useMatrix = "GeneIntegrationMatrix", log2Norm = FALSE)

## Creating Trajectory Group Matrix..
## Smoothing…

p3 <- plotTrajectoryHeatmap(trajGIM,  pal = paletteContinuous(set = "blueYellow"))

## ArchR logging to : ArchRLogs/ArchR-plotTrajectoryHeatmap-15e1333158a4-Date-2020-04-15_Time-15-02-18.log
## If there is an issue, please report to github with logFile!
## Preparing Main Heatmap..
## ArchR logging successful to : ArchRLogs/ArchR-plotTrajectoryHeatmap-15e1333158a4-Date-2020-04-15_Time-15-02-18.log

Lastly, we can obtain a pseudo-time heatmap of peak accessibility by setting useMatrix = "PeakMatrix".

trajPM  <- getTrajectory(ArchRProj = projHeme5, name = "MyeloidU", useMatrix = "PeakMatrix", log2Norm = TRUE)

## Creating Trajectory Group Matrix..
## Smoothing…

p4 <- plotTrajectoryHeatmap(trajPM, pal = paletteContinuous(set = "solarExtra"))

## ArchR logging to : ArchRLogs/ArchR-plotTrajectoryHeatmap-15e131ecfb872-Date-2020-04-15_Time-15-03-01.log
## If there is an issue, please report to github with logFile!
## Preparing Main Heatmap..
## ArchR logging successful to : ArchRLogs/ArchR-plotTrajectoryHeatmap-15e131ecfb872-Date-2020-04-15_Time-15-03-01.log

To save an editable vectorized version of this plot, we use plotPDF().

plotPDF(p1, p2, p3, p4, name = "Plot-MyeloidU-Traj-Heatmaps.pdf", ArchRProj = projHeme5, addDOC = FALSE, width = 6, height = 8)

## [1] 0

16.1.3 Integrative pseudo-time analyses

We can also perform integrative analyses, such as identification of positive TF regulators by integration of gene scores / gene expression with motif accessibility across pseudo-time. This can be very powerful, for example in identifying drivers of differentiation. To do this, we use the correlateTrajectories() function which takes two SummarizedExperiment objects retrived from the getTrajectories() function.

First, lets find motifs whose accessibility across pseudo-time is correlated with the gene score of the TF gene.

corGSM_MM <- correlateTrajectories(trajGSM, trajMM)

## ArchR logging to : ArchRLogs/ArchR-correlateTrajectories-15e13b9ca02a-Date-2020-04-15_Time-15-03-37.log
## If there is an issue, please report to github with logFile!
## Found 36 Correlated Pairings!
## 2020-04-15 15:03:56 :

The primary output of correlateTrajectories() is a list object containing a DataFrame object as the first entry in the list. This DataFrame has columns named idx1, matchname1, name1, and VarAssay1 which correspond to the index, match name, unaltered name, and the variance quantile of the features from the first trajectory (gene scores) passed to the correlateTrajectories() function. A “variance quantile” is a normalized measure of the given feature which allows us to derive a correlation across disparate assays. This DataFrame contains all of the features that met the cutoffs specified in the correlateTrajectories() function.

corGSM_MM[[1]]

## DataFrame with 36 rows and 12 columns
## idx1 idx2 matchname1 matchname2 name1 name2
##
## 1 82 1081 PRDM16 PRDM16 chr1:PRDM16 z:PRDM16_211
## 2 731 932 TAL1 TAL1 chr1:TAL1 z:TAL1_62
## 3 818 1013 JUN JUN chr1:JUN z:JUN_143
## 4 2034 1002 ATF3 ATF3 chr1:ATF3 z:ATF3_132
## 5 2369 1254 GATA3 GATA3 chr10:GATA3 z:GATA3_384
## … … … … … … …
## 32 20775 1022 CEBPD CEBPD chr8:CEBPD z:CEBPD_152
## 33 20780 1031 SNAI2 SNAI2 chr8:SNAI2 z:SNAI2_161
## 34 21658 1003 NFIL3 NFIL3 chr9:NFIL3 z:NFIL3_133
## 35 21793 1078 KLF4 KLF4 chr9:KLF4 z:KLF4_208
## 36 22097 1565 RXRA RXRA chr9:RXRA z:RXRA_695
## Correlation VarAssay1 VarAssay2 TStat
##
## 1 0.859588836579798 0.999783802481948 0.860344827586207 16.6530781620713
## 2 0.687938524545333 0.866735849872443 0.987931034482759 9.38348743943718
## 3 0.530335586469274 0.993989708998141 0.963218390804598 6.19265152142139
## 4 0.776280918431044 0.985168850261599 0.936781609195402 12.1905867620732
## 5 0.714559116657296 0.801573917931422 0.989655172413793 10.111530171907
## … … … … …
## 32 0.73809434207349 0.993816750983699 0.998275862068966 10.8296844080467
## 33 0.643337129188769 0.87252994335625 0.983333333333333 8.31877485229719
## 34 0.833740261914099 0.992303368357331 0.893103448275862 14.9479403447508
## 35 0.82265894955927 0.998702814891685 0.905747126436782 14.3243696597426
## 36 0.59509316489084 0.88640982401522 0.890229885057471 7.33039570155334
## Pval FDR
##
## 1 2.47466579579935e-30 4.04855324192774e-27
## 2 2.63756432350363e-15 1.00350121703533e-13
## 3 1.38781951660218e-08 1.55511830764464e-07
## 4 2.36537648993666e-21 3.22479661461365e-19
## 5 6.9322610846337e-17 3.78039304482025e-15
## … … …
## 32 1.93304822797782e-18 1.31769454207155e-16
## 33 5.27319420246801e-13 1.3914428572964e-11
## 34 5.06319064599113e-27 2.70186610322912e-24
## 35 8.9984663184347e-26 2.45358181615986e-23
## 36 6.60626612679289e-11 1.25672690505037e-09

We can then subset our corresponding trajectory SummarizedExperiment objects to only contain the elements that passed significance above.

trajGSM2 <- trajGSM[corGSM_MM[[1]]$name1, ]
trajMM2 <- trajMM[corGSM_MM[[1]]$name2, ]

To best order these features, we can create a new trajectory where the values of these two trajectories are multiplied. This will allow us to create side-by-side heatmaps that are identically ordered by row.

trajCombined <- trajGSM2
assay(trajCombined) <- t(apply(assay(trajGSM2), 1, scale)) + t(apply(assay(trajMM2), 1, scale))

We can extract the optimal row order from the return of the plotTrajectoryHeatmap() function.

combinedMat <- plotTrajectoryHeatmap(trajCombined, returnMat = TRUE, varCutOff = 0)

## ArchR logging to : ArchRLogs/ArchR-plotTrajectoryHeatmap-15e136ff55add-Date-2020-04-15_Time-15-03-56.log
## If there is an issue, please report to github with logFile!
## Preparing Main Heatmap..
## ArchR logging successful to : ArchRLogs/ArchR-plotTrajectoryHeatmap-15e136ff55add-Date-2020-04-15_Time-15-03-56.log

rowOrder <- match(rownames(combinedMat), rownames(trajGSM2))

With this, we are now ready to create our paired heatmaps. First, we will create the heatmap for the gene score trajectory. We specify the desired row order via the rowOrder parameter.

ht1 <- plotTrajectoryHeatmap(trajGSM2,  pal = paletteContinuous(set = "horizonExtra"),  varCutOff = 0, rowOrder = rowOrder)

## ArchR logging to : ArchRLogs/ArchR-plotTrajectoryHeatmap-15e13aac5f3b-Date-2020-04-15_Time-15-03-57.log
## If there is an issue, please report to github with logFile!
## Preparing Main Heatmap..
## ArchR logging successful to : ArchRLogs/ArchR-plotTrajectoryHeatmap-15e13aac5f3b-Date-2020-04-15_Time-15-03-57.log

Then, we will create the heatmap for the motif trajectory, again specifying the row order via the rowOrder parameter

ht2 <- plotTrajectoryHeatmap(trajMM2, pal = paletteContinuous(set = "solarExtra"), varCutOff = 0, rowOrder = rowOrder)

## ArchR logging to : ArchRLogs/ArchR-plotTrajectoryHeatmap-15e1323caf784-Date-2020-04-15_Time-15-03-57.log
## If there is an issue, please report to github with logFile!
## useSeqnames is NULL or greater than 1 with a Sparse.Assays.Matrix trajectory input.
## 2020-04-15 15:03:58 :
## force=FALSE thus continuing with subsetting useSeqnames = z
## 2020-04-15 15:03:58 :
## Preparing Main Heatmap..
## ArchR logging successful to : ArchRLogs/ArchR-plotTrajectoryHeatmap-15e1323caf784-Date-2020-04-15_Time-15-03-57.log

Plotting these two heatmaps side-by-side, we see that the rows are matched across the two heatmaps. You may notice that the analysis captures both GATA3 and GATA3-AS1 (an anti-sense transcript of GATA3). This is due to how feature matching is performed and the anti-sense transcript entry could be removed manually in post-processing or programatically if desired.

ht1 + ht2

We can repeat this same exact process but use gene expression from the GeneIntegrationMatrix instead of gene scores. Because this is the same analytical workflow, we do not repeat our explanations for each step.

corGIM_MM <- correlateTrajectories(trajGIM, trajMM)

## ArchR logging to : ArchRLogs/ArchR-correlateTrajectories-15e13f7dfd04-Date-2020-04-15_Time-15-04-14.log
## If there is an issue, please report to github with logFile!
## Found 54 Correlated Pairings!
## 2020-04-15 15:04:30 :

corGIM_MM[[1]]

## DataFrame with 54 rows and 12 columns
## idx1 idx2 matchname1 matchname2 name1 name2
##
## 1 295 1601 RUNX3 RUNX3 chr1:RUNX3 z:RUNX3_731
## 2 680 1612 NFIA NFIA chr1:NFIA z:NFIA_742
## 3 1227 1512 MEF2D MEF2D chr1:MEF2D z:MEF2D_642
## 4 1317 930 USF1 USF1 chr1:USF1 z:USF1_60
## 5 1664 1002 ATF3 ATF3 chr1:ATF3 z:ATF3_132
## … … … … … … …
## 50 16981 1696 KLF10 KLF10 chr8:KLF10 z:KLF10_826
## 51 17463 1003 NFIL3 NFIL3 chr9:NFIL3 z:NFIL3_133
## 52 17560 1078 KLF4 KLF4 chr9:KLF4 z:KLF4_208
## 53 17803 1565 RXRA RXRA chr9:RXRA z:RXRA_695
## 54 18090 888 TFE3 TFE3 chrX:TFE3 z:TFE3_18
## Correlation VarAssay1 VarAssay2 TStat
##
## 1 0.506721180509532 0.871189720982743 0.98448275862069 5.81861222342308
## 2 0.899884055575892 0.878017310897264 0.919540229885057 20.4260514493749
## 3 0.787682536829653 0.815493790656416 0.841379310344828 12.656866457322
## 4 0.653303634000368 0.804204075049728 0.875287356321839 8.54234597517229
## 5 0.847460093565636 0.813020805333047 0.936781609195402 15.8034906427774
## … … … … …
## 50 0.57964654735637 0.953335842159024 0.864367816091954 7.04188649257934
## 51 0.849563937747037 0.933928283425622 0.893103448275862 15.9440428922802
## 52 0.820699669277503 0.98876404494382 0.905747126436782 14.2196070232556
## 53 0.899871735441299 0.960485995376593 0.890229885057471 20.424581355724
## 54 0.832322259515336 0.801247244771786 0.886206896551724 14.8650341763864
## Pval FDR
##
## 1 7.47838268446286e-08 2.00835886256963e-07
## 2 4.28849849649256e-37 2.40519926139511e-35
## 3 2.45911764921563e-22 2.0096926995314e-21
## 4 1.7427906036477e-13 7.90509804908143e-13
## 5 1.05531067509617e-28 1.5631789374862e-27
## … … …
## 50 2.63272438071192e-10 8.66605108651006e-10
## 51 5.63648466586524e-29 9.05658892074619e-28
## 52 1.46590261078753e-25 1.5687313439583e-24
## 53 4.31312103836676e-37 2.40519926139511e-35
## 54 7.4030081939147e-27 8.88360983269764e-26

trajGIM2 <- trajGIM[corGIM_MM[[1]]$name1, ]
trajMM2 <- trajMM[corGIM_MM[[1]]$name2, ]

trajCombined <- trajGIM2
assay(trajCombined) <- t(apply(assay(trajGIM2), 1, scale)) + t(apply(assay(trajMM2), 1, scale))

combinedMat <- plotTrajectoryHeatmap(trajCombined, returnMat = TRUE, varCutOff = 0)

## ArchR logging to : ArchRLogs/ArchR-plotTrajectoryHeatmap-15e137bbcc972-Date-2020-04-15_Time-15-04-30.log
## If there is an issue, please report to github with logFile!
## Preparing Main Heatmap..
## ArchR logging successful to : ArchRLogs/ArchR-plotTrajectoryHeatmap-15e137bbcc972-Date-2020-04-15_Time-15-04-30.log

rowOrder <- match(rownames(combinedMat), rownames(trajGIM2))
ht1 <- plotTrajectoryHeatmap(trajGIM2,  pal = paletteContinuous(set = "blueYellow"),  varCutOff = 0, rowOrder = rowOrder)

## ArchR logging to : ArchRLogs/ArchR-plotTrajectoryHeatmap-15e1331e06f29-Date-2020-04-15_Time-15-04-31.log
## If there is an issue, please report to github with logFile!
## Preparing Main Heatmap..
## ArchR logging successful to : ArchRLogs/ArchR-plotTrajectoryHeatmap-15e1331e06f29-Date-2020-04-15_Time-15-04-31.log

ht2 <- plotTrajectoryHeatmap(trajMM2, pal = paletteContinuous(set = "solarExtra"), varCutOff = 0, rowOrder = rowOrder)

## ArchR logging to : ArchRLogs/ArchR-plotTrajectoryHeatmap-15e1345e2da25-Date-2020-04-15_Time-15-04-32.log
## If there is an issue, please report to github with logFile!
## useSeqnames is NULL or greater than 1 with a Sparse.Assays.Matrix trajectory input.
## 2020-04-15 15:04:32 :
## force=FALSE thus continuing with subsetting useSeqnames = z
## 2020-04-15 15:04:32 :
## Preparing Main Heatmap..
## ArchR logging successful to : ArchRLogs/ArchR-plotTrajectoryHeatmap-15e1345e2da25-Date-2020-04-15_Time-15-04-32.log

ht1 + ht2