1. Cell Biology
  2. Developmental Biology
Download icon

Single-cell transcriptome analysis of avian neural crest migration reveals signatures of invasion and molecular transitions

  1. Jason A Morrison
  2. Rebecca McLennan
  3. Lauren A Wolfe
  4. Madelaine M Gogol
  5. Samuel Meier
  6. Mary C McKinney
  7. Jessica M Teddy
  8. Laura Holmes
  9. Craig L Semerad
  10. Andrew C Box
  11. Hua Li
  12. Kathryn E Hall
  13. Anoja G Perera
  14. Paul M Kulesa  Is a corresponding author
  1. Stowers Institute for Medical Research, United States
  2. University of Nebraska Medical Center, United States
  3. University of Kansas School of Medicine, United States
Tools and Resources
  • Cited 7
  • Views 3,129
  • Annotations
Cite this article as: eLife 2017;6:e28415 doi: 10.7554/eLife.28415

Abstract

Neural crest cells migrate throughout the embryo, but how cells move in a directed and collective manner has remained unclear. Here, we perform the first single-cell transcriptome analysis of cranial neural crest cell migration at three progressive stages in chick and identify and establish hierarchical relationships between cell position and time-specific transcriptional signatures. We determine a novel transcriptional signature of the most invasive neural crest Trailblazer cells that is consistent during migration and enriched for approximately 900 genes. Knockdown of several Trailblazer genes shows significant but modest changes to total distance migrated. However, in vivo expression analysis by RNAscope and immunohistochemistry reveals some salt and pepper patterns that include strong individual Trailblazer gene expression in cells within other subregions of the migratory stream. These data provide new insights into the molecular diversity and dynamics within a neural crest cell migratory stream that underlie complex directed and collective cell behaviors.

https://doi.org/10.7554/eLife.28415.001

Introduction

Neural crest migration is crucial to vertebrate organogenesis. Signals within the dorsal neural tube and microenvironments through which neural crest cells travel direct cells in multicellular, discrete streams to precise peripheral targets. Time-lapse analyses have captured the complexity of neural crest cell behaviors that include the transition from the dorsal neural tube to directed migration and movement in collective groups. Cell morphologies and behaviors appear to vary depending on cell position within a migratory stream. Despite these observations at the population level, our understanding of how migrating neural crest cells interpret and respond to microenvironmental signals at the molecular level and coordinate their behaviors remains poorly understood. Thus, there is a tremendous need to interrogate the genomic features of migrating neural crest cells in a spatial and temporal manner with the goal of elucidating the underlying mechanisms of directed and collective cell migration.

The finely-sculpted pattern of discrete neural crest cell migratory streams and their broad differentiation potential have led to fascinating questions of the underlying biology (Kulesa and Gammill, 2010). Neural crest cells undergo an epithelial-to-mesenchymal transition to exit the dorsal neural tube in a rostral-to-caudal timed manner. As cells exit the neural tube and encounter the paraxial mesoderm, signals from within the neural tube and paraxial mesoderm sculpt cells onto stereotypical migratory pathways. By navigating through different microenvironments in a directed manner and moving collectively in discrete streams, subpopulations of neural crest cells are distributed throughout the embryo. In the head, neural crest cells migration directs cells into regions of the face, eye, and branchial arches where cells play an important role in craniofacial patterning. Cranial neural crest cells travel in multicellular streams just underneath the surface ectoderm, presenting the opportunity to visualize and interrogate the in vivo cellular and molecular features of the migration process. Furthermore, the chick embryo is well suited to studying gene expression heterogeneities within neural crest cell migratory streams since cell isolation and single-cell RT-qPCR techniques have been optimized (Morrison et al., 2012, 2015).

Genomic analyses of premigratory neural crest cells or entire stream populations have revealed transcription factors and signaling pathways that may play critical roles during neural crest induction, specification, epithelial-to-mesenchymal transition and early migratory events (Gammill and Bronner-Fraser, 2002; Adams et al., 2008; Simões-Costa et al., 2014; Simões-Costa and Bronner, 2015). This information is typically gained by isolation and bulk RNA-seq analysis of premigratory and migrating neural crest cells at a single time point, yet provides an important foundation to begin functional studies to further elucidate aspects of neural crest biology (Simões-Costa and Bronner, 2015, 2016). However, what remains unclear is our understanding of the dynamic nature of gene expression changes within migrating neural crest cells and how molecular patterns are transduced into directed and collective cell migration behaviors.

To begin to address this, we recently profiled subpopulations of migrating cranial neural crest cells manually dissected and isolated by fluorescence activated cell sorting (FACS) from distinct stream positions at successive developmental stages by RT-qPCR in chick (McLennan et al., 2012, 2015). We discovered regional expression of genes (96 genes selected from the literature) within a typical cranial neural crest cell migratory stream and confirmed the expression of a subset of genes by multiplexed fluorescence in-situ hybridization and quantitative measurement (McLennan et al., 2012, 2015). Intriguingly, by optimizing a single-cell RT-qPCR method (Morrison et al., 2015), we were able to identify a unique molecular signature of 16 out of 96 genes that was stable and consistent in a subset of cranial neural crest cells confined to the invasive front of the migratory stream that we termed ‘trailblazers’ (McLennan et al., 2015). Over-expression of either of two transcription factors (TFAP2A, HAND2) upstream of the trailblazer genes in the follower cell subpopulation led to alterations in stream morphology (McLennan et al., 2015) suggesting a functional role for trailblazer genes.

Here, we leverage our previous RT-qPCR data and exploit advances in whole transcriptome amplification and cell isolation to enable single-cell RNA-seq (scRNA-seq) of neural crest migration. Our goal was to distinguish inherent variations between individual cells within a neural crest cell migratory stream. To do this, we compare gene expression fingerprints from 469 cells collected from the front, lead, and trail subregions of a typical cranial neural crest cell migratory stream at three developmental stages corresponding to initiation of migration from the dorsal neural tube, active migration towards and colonization of the branchial arches. We additionally profiled ~100 migrating neural crest cells collected from explanted cranial neural tube cultures for comparison to in vivo transcriptional signatures. We determine unique transcriptional signatures based on cell position within a stream during migration using bioinformatics tools such as unbiased hierarchical clustering and principle component analyses (PCA). Our bioinformatics analysis provides spatial transcriptomic signatures for neural crest cell subpopulations within the migratory stream and in vitro, independent of manual dissection location. We examine expression of a small subset of the trailblazer transcriptional signature by RNAscope and perform loss-of-function of a subset of trailblazer genes and the effects on neural crest migration. These data provide a comprehensive knowledge base from which to identify and functionally test genes and signaling pathways critical to neural crest migration and organogenesis.

Results

Bulk RNA-seq analysis affirms a rich diversity of gene expression between the invasive front and remainder of the neural crest cell migratory stream

With our previous RT-qPCR analysis that revealed a set of 16 out of 96 genes highly expressed and consistent during migration in a subset of cells restricted to the invasive front, we sought to expand the gene list in an unbiased manner. To do this, we first performed bulk RNA-seq on cranial neural crest cells isolated from the invasive front and remainder of the migratory stream at two successive developmental stages during migration and branchial arch entry (Figure 1A and Figure 1—figure supplement 1; Hamburger and Hamilton, 1951, HHSt13 and HHSt15). At these successive developmental stages, chick cranial neural crest cells of the pre-otic stream adjacent to rhombomere 4 (r4) are mid-migration and colonizing the branchial arch targets, respectively (Figure 1A). We fluorescently labeled premigratory neural crest cells with Gap43-YFP (at HHSt8-9; 4–7 somites), reincubated eggs until HHSt13 or HHSt15 and carefully manually dissected cells from the migratory front (5% of the stream, termed ‘Front’) and remainder of the stream (95% of the stream, termed ‘Stream’) (Figure 1A). We then used fluorescence-activated cell sorting (FACS) to isolate neural crest cells from manually dissected tissues to perform bulk RNA-seq.

Figure 1 with 1 supplement see all
Bulk RNA-seq confirms the unique molecular profile of cells at the invasive front of the neural crest stream.

(A) Schematic representation of method used for harvesting samples from the cranial NC stream. Invasive front (Front) is the ventral-most 5% of NC cells. Stream is the remaining 95% of NC cells. n values as shown. (B) Heatmap of biological replicates of differentially expressed genes in HHSt13 and 15 Front and Stream (n = 477 genes). Black boxes contain gene clusters expressed at significantly higher and lower levels by the Front. (C) Heatmap of mean values of differentially expressed genes in HHSt13 and 15 Front and Stream samples (n = 477 genes). Black boxes contain gene clusters expressed at significantly higher and lower levels by the Front. (D) Venn diagram displaying the numbers of genes enriched at the Front of the NC stream at both HHSt13 and 15 with a union of 23 genes enriched in both comparisons. (E) Venn diagram displaying the numbers of genes reduced at the front of the NC stream at both HHSt13 and 15 with a union of 19 genes reduced in both comparisons. Hamburger and Hamilton (1951); St, stage.

https://doi.org/10.7554/eLife.28415.002

Hierarchical clustering showed three distinct features of the data (Figure 1B). First, we find good correlation of expression levels between the three biological replicates, validating our cells isolations (Figure 1B). Second, there are statistically significant differences in the expression levels of subsets of genes between the different subpopulations isolated, confirming a conserved invasion signature within the front migrating cells (Figure 1B and Figure 1—source data 1). Third, subsets of genes expressed in the front subpopulation were either consistently up- or down-regulated in comparison to the remainder of the stream at both developmental stages (Figure 1B and Figure 1—source data 1).

In order to analyze the bulk RNA-seq data in more general terms, we calculated the mean of the biological replicates and determined the hierarchical clustering based on mean values (Figure 1C). These data reaffirmed the significant differences in the expression levels of subsets of genes between the front and remainder of migrating neural crest cells, and included aspects of the 16 gene scRT-qPCR Trailblazer signature (McLennan et al., 2015) (Figure 1C and Figure 1—source data 1). We find there is a subset of genes enriched in some neural crest cells at the invasive front compared with the remainder of the stream, with consistent expression at both developmental stages and unique to the previous 16 gene scRT-qPCR signature (Figure 1D and Figure 1—source data 1; n = 23 genes). In contrast, we find a subset of genes reduced in migrating front cells compared with the remainder of the stream, with consistent expression at both developmental stages (Figure 1E; m = 19 genes). In comparison to our previous identification of 16 genes out of 100 analyzed by scRT-qPCR (McLennan et al., 2015), bulk RNA-seq analysis validated that ITGB5 and GPC3 (HHSt13), and BAMBI and PKP2 (HHSt15) were enriched in FRONT versus stream samples. EPHA4 (HHSt13) and CDH7 (HHSt15) were reduced. Together, these bulk RNA-seq analyses affirm there is a rich spatio-temporal diversity of gene expression depending on whether a neural crest cell is within the invasive front versus any other position within the stream and reveal genes that are either enhanced or reduced consistently at the invasive front.

Single-cell RNA-seq identifies gene expression variances based upon spatial position within the neural crest cell stream and temporal progression along the migratory pathway

To better characterize unique transcriptional signatures and gene expression heterogeneity during cranial neural crest migration, we isolated and profiled individual cells from different stream positions at three developmental stages (HHSt11,13,15; Figure 2A, Figure 2—source data 1 and Figure 2—figure supplement 1). These three progressive developmental stages were selected based on the different migratory events including recently delaminated from the neural tube (HHSt11), invasion of the paraxial mesoderm (HHSt13) and entry into the second branchial arch (HHSt15) with ~8 hr in between the stages. Since there are few recently emigrated neural crest cells at HHSt11, we could take advantage of single-cell analysis whereas with bulk RNA-seq described above we would not have been able to further dissect into invasive front and stream subpopulations. This unique approach maintained a level of spatial information within our scRNA-seq data set.

Figure 2 with 2 supplements see all
Single-cell RNA-seq shows in vitro and in vivo neural crest have distinct molecular signatures.

(A) Schematic representation of method used for harvesting samples from the cranial NC stream. Front is the ventral-most 5% of NC cells. Lead is the 25% of NC cells immediately following the Front cells. Trail is the remaining dorsal-most 70% of NC cells. Neural crest cells grown in vitro overnight from isolated neural tubes were also isolated for this analysis. (B) Numbers of single cells isolated from different subpopulations. (C) Pearson correlation matrix of the average expression profiles, based upon all differentially expressed genes all subpopulations analyzed (n = 1355 genes). (D) Principal Component Analysis (PCA) of all in vitro and in vivo single cells using all differentially expressed genes (n = 1355 genes). (E) Condensed Gene Ontology (GO) analysis of genes upregulated in in vivo NC cells compared to in vitro NC cells (n = 806 genes). (F) Bar chart displaying NC relevant genes increased in in vivo NC cells compared to in vitro NC cells. HH, Hamburger and Hamilton (1951); St, stage.

https://doi.org/10.7554/eLife.28415.005

Premigratory neural crest cells were fluorescently labeled at HHSt8 with Gap-43-YFP for ease of identification during migration through unlabeled tissue for segregation of neural crest cells from unlabeled tissue by FACS (Figure 2A). We subdivided neural crest cell migratory streams at HHSt13 and HHSt15 into three spatially distinct, non-overlapping subregions by careful manual dissection and single cells were isolated by FACS. This included neural crest cells within the invasive front (5% of the stream, termed ‘Front’), lead subregion of the stream (25% of the stream just proximal to the front, termed ‘Lead’), and remainder of the stream (70% of the stream, termed ‘Trail’) (Figure 2A). At HHSt11, we did not subdivide the stream since there are relatively small numbers of migrating neural crest cells adjacent to rhombomere 4 (r4) and refer to these cells as ‘Front’ (Figure 2A). The total number of single cells analyzed after quality control was 469 (Figure 2B). Hierarchical clustering of single neural crest cell averages revealed statistically enriched and reduced genes between every pairwise combination of subpopulations and concordance between Lead cells at HHSt13 and HHSt15 and Trail cells at HHSt13 and HHSt15 (Figure 2—figure supplement 2). These data show that spatial location within the stream and temporal progression along the migratory pathway are primary determinants of gene expression.

In vitro transcriptional signatures of migrating neural crest cells are distinct from any in vivo profile

Neural tube explant cultures are often used as an assay for neural crest cell migration studies. However, we previously showed that there are gene expression profile differences between neural crest cells grown in vitro and neural crest cells isolated from the invasive front in vivo by RT-qPCR (McLennan et al., 2015). That is, we found significant differences between in vitro and in vivo RT-qPCR gene expression profiles by Euclidean clustering and dissimilarity matrix plots (McLennan et al., 2015) that might explain the lack of collective cell migration in culture.

To better characterize in vitro cranial neural crest cell gene expression profiles, we performed scRNA-seq on neural crest cells derived from explanted neural tubes. To do this, we explanted cranial neural tubes dissected from HHSt8-9 embryos onto fibronectin-coated Mat-tek dishes and allowed cells to invade the dish for 24 hr, comparable to in vivo HHSt15. At t =+ 24 hr, neural tube explants were removed from the dish leaving the migrating cells to be collected and isolated by FACS for single-cell RNA-seq analysis. We identified several key features when we compared in vitro versus in vivo neural crest cell gene expression profiles (Figure 2—source data 2). First, using a Pearson correlation matrix to relate the average expression profile of all migrating neural crest cell subpopulations to show that the average in vitro profile is prominently distinct from all in vivo profiles analyzed (Figure 2C). Specifically, the in vitro neural crest cell gene expression profile correlates poorly with every isolated in vivo subpopulation. The closest correlation is cells most recently exited from the dorsal neural tube (HHSt11) (Figure 2C; Figure 2—figure supplement 2). To definitively ensure that the difference between in vitro and in vivo neural crest cell expression was not the product of averaging, we performed a principal component analysis (PCA) containing all single-cell RNA-seq profiles (Figure 2D). The PCA reveals that all in vitro single neural crest cells cluster together and show minimal overlap with any in vivo single-cell expression profile (Figure 2D).

To gain insight into the gene expression differences between migrating neural crest cells isolated from the embryo and in vitro, we performed a Condensed Gene Ontology (GO) Term analysis. From this analysis, there is enrichment of genes associated with biological processes such as cell signaling, transcription regulation, growth factor response, migration and proliferation by in vivo neural crest cells when compared to in vitro data (Figure 2E). We find that 806 genes are down-regulated in neural crest cells grown in vitro (e.g., NKX2-5, NKX2-6, WNT5B, SNAI1, PAX5, PAX9, HAND1, PHOX2B, EPHA3, EPHA5, NCAM2, FGFs 4, 7, 8, 9, 10, 18, 19, AQPs 1, 3, 5, NEUROD1, NGFR, NRP1, DLXs 1, 3, 5, 6, TWIST1, TGFBR2, DSP, CDHs 5, 10, CXCL12, GATA4, GATA6, GPC3, ISL1, BMP2, and KDR (VEGFR2)) compared with in vivo data (Figure 2F and Figure 2—source data 2). We also find 312 genes up-regulated in migrating neural crest cells in vitro (e.g., WNT9A, SHH, SOX2, PAX7, IL8, MKX, CDH6 and CDH12) (Figure 2—source data 2). However, despite the identification of 312 genes up-regulated in migrating neural crest cells in vitro, we did not find any GO terms associated with these genes. These data strongly suggest that in vivo embryonic neural crest cell microenvironmental signals influence both the spatial and temporal heterogeneities of neural crest cell gene expression and biological processes critical to migration and embryogenesis.

Evidence for unique transcriptional signatures associated with cell subpopulations within the in vivo neural crest cell migratory stream

When we analyzed the heatmap of averaged expression profiles of all subpopulations we found regional and temporal expression differences (Figure 2C; Figure 2—figure supplement 2). This suggested a richness of diversity between each of the subpopulations analyzed. To gain insight into the extent of this diversity, we clustered all single-cell profiles as a heatmap based on subpopulation (Figure 3). First, we confirm that in vitro neural crest cells display up- and down-regulation of genes distinct from the in vivo subpopulations (Figure 3; column of cells marked ‘in vitro’). Second, we identify subsets of genes that either increased (Figure 3, box A) or decreased (Figure 3, box B) with developmental stage and spatial location within the stream. Third, we discover subsets of genes that were uniquely reduced (Figure 3, box A’) or enhanced (Figure 3, box C) in individual subpopulations of migrating neural crest cells. Fourth, we find that some genes show expression changes that correlate with spatial position within the stream at a single developmental time point (Figure 3, box D). Lastly, we identify that heterogeneity of expression among single cells within any subpopulation was common (Figure 3). A Pearson correlation matrix of all single-cell transcriptomes confirmed similarity among some, but not all cells within each subpopulation (Figure 3—figure supplement 1), indicating that our spatial-temporal locations are not completely homogeneous. These results demonstrate that each cranial neural crest cell subpopulation analyzed (Front, Lead, Trail) can be defined by distinct molecular signatures, but also include heterogeneity of expression.

Figure 3 with 1 supplement see all
Single-cell RNA-seq exposes novel subpopulations within the neural crest stream.

Heatmap of all single cells (n = 469 single cells, 1355 genes) analyzed by RNA-seq with single cells clustered by spatiotemporal location. (Box A) Genes increased with developmental time and spatial location. (Box A’) Example of genes uniquely downregulated by an individual subpopulation. (Box B) Genes decreased with developmental time and spatial location. (Box C) Example of genes uniquely upregulated by an individual subpopulation. (Box D) Example of regional differences within one developmental stage.

https://doi.org/10.7554/eLife.28415.010

Hierarchical clustering and PCA analysis identify distinct transcriptional signatures associated with the invasive front of the neural crest migratory stream

Our previous RT-qPCR (McLennan et al., 2015) and bulk RNA-seq (Figure 1) analyses showed that cranial neural crest cells at the invasive front have portions of their expression profiles that are conserved as cells travel through different microenvironments. To determine the transcriptional signatures associated with neural crest cells at the invasive front of the migratory stream, we clustered the invasive front cells in a heatmap based on gene expression, rather than stage of development (Figure 4A). First, hierarchical clustering reveals that the majority of neural crest cells at the invasive front cluster together by developmental stage (Figure 4A; Figure 4—source data 1). Further, a Pearson correlation matrix shows that neural crest cells within the invasive front at HHSt13 and HHSt15 have expression profiles that cluster most similarly with Lead cell subpopulations, which we termed Lead-like (Figure 4B). Second, we discover two distinct cell subpopulations exist at HHSt11 termed Cluster1 and Cluster2 (Figure 4A–B). Third, we find a distinct cluster that contains neural crest cells isolated from the invasive front at all three developmental stages analyzed, termed scRNA-seq Trailblazers (Figure 4A). The single-cell PCA plot shows these five distinct cell subpopulations within the invasive front over time, with very minimal overlap in expression profile between the Trailblazer cluster and any other cell subpopulation cluster (Figure 4C; light blue dots). The identification of a transcriptomic signature associated with Trailblazer cells within the invasive front represents an unbiased and more comprehensive molecular profile than we previously described by scRT-qPCR (McLennan et al., 2015).

Trailblazers are a small subpopulation of cells with a distinct expression signature that exist at the invasive front of the neural crest stream, regardless of developmental stage.

(A) Heatmap with hierarchical clustering of all single cells from the invasive front at HHSt11, 13 and 15 (n = 187 cells), based upon differentially expressed genes from all invasive front comparisons and invasive front compares to stream at HHSt13 and HHSt15 (n = 760 genes). (B) Pearson correlation matrix of the average expression profile of all subpopulations based upon differentially expressed genes from all invasive front comparisons and comparisons between clusters and all other in vivo cells (n = 3555 genes). (C) PCA of all single cells from the invasive front at all three developmental time points, based upon differentially expressed genes from all invasive front comparisons and invasive front compares to stream at HHSt13 and HHSt15 (n = 760 genes).

https://doi.org/10.7554/eLife.28415.012

The scRNA-seq trailblazer signature is uniquely defined and associated with neural crest cells mostly confined to the invasive front of the stream

To characterize the scRNA-seq Trailblazer transcriptional signature, we identified genes that are statistically enriched or reduced (Figure 5A and Figure 4—source data 1). A volcano plot shows that 964 genes are up-regulated (e.g., BAMBI, LUM, MMP2, GATA2, GATA4, GATA5 and GATA6) and 406 genes are down-regulated (e.g., DRAXIN, TFAP2A, TFAP2B, TFAP2C, TFAP2E and SOX10) within the Trailblazers (Figure 5A and Figure 4—source data 1). Nineteen out of 24 genes (79%) enriched in the invasive front from our bulk RNA-seq analysis (Figure 1D) were also enriched in Trailblazers by scRNA-seq (Figure 4—source data 1). Additionally, 7 out of 19 genes (37%) reduced in the invasive front from our bulk RNA-seq analysis (Figure 1E) were also reduced in Trailblazers by scRNA-seq (Figure 4—source data 1). We also find examples of known neural crest-related genes that are up-regulated in Trailblazers (e.g., EDN1, GPC3, NRP1, BMP2) and include genes (AQP1, BAMBI, PKP2) identified in our RT-qPCR analysis (McLennan et al., 2015) (Figure 5B). There are also neural crest-related genes down-regulated in Trailblazers (e.g., MSX1, EDNRB, PAX3, SOX10) (Figure 5B and Figure 4—source data 1).

The Trailblazer signature represents tissue development, adhesion and a lack of differentiation.

(A) Volcano plot of all genes statistically enriched or reduced in Trailblazers (n = 47 cells, 1370 genes). (B) Bar chart showing examples of neural crest-related differentially expressed genes in Trailblazers (red designates previously identified in McLennan et al., 2015). (C) PCA of all single cells analyzed showing minimal overlap of other in vivo (black dots) and in vitro (gray dots) with the Trailblazers (green dots) (n = 1370 genes). (D) Significant signaling pathways predicted to be enriched (red) or reduced (blue) in Trailblazers. VEGF and axonal guidance signaling are highlighted. The statistical significance was set at log (−1.3) [p<0.05]. (E) Comparison to previous RT-qPCR generated trailblazer signature. (F) Genes enriched and detected in >95% of trailblazer cells.

https://doi.org/10.7554/eLife.28415.014

With an ensemble of genes that we use to define Trailblazer cells, we may characterize the presence of these cells within the invasive front. Intriguingly, the number of neural crest cells within the invasive front (2–3% of the entire migratory stream) with the Trailblazer transcriptional signature increased significantly during HHSt11 to HHSt13 and was then consistent from HHSt13 to HHSt15. At HHSt11 we found that approximately 8% of the migrating neural crest cells had a Trailblazer transcriptional signature. By HHSt13, with an estimate of approximately 500 neural crest cells in the typical migratory stream at this stage and axial level, we found that 49% of the cells within the invasive front had a Trailblazer transcriptional signature, or around 12 cells. At HHSt15, approximately 32% (or around 13 cells within the invasive front of a stream of approximately 800 cells) had a Trailblazer transcriptional signature. In an effort to be more precise about the Trailblazer signature, we filtered the 964 genes enriched in Trailblazer cells to include only genes detected in the vast majority (>90%) of trailblazer cells have a more stringent criteria (Figure 5F).

To determine whether the scRNA-seq Trailblazer transcriptional signature is present elsewhere in the neural crest cell migratory stream, independent of the Pearson correlation matrix result (Figure 5C), we generated a PCA plot of all analyzed single-cells based on the Trailblazer transcriptional signature genes (Figure 5C). We find that fewer than 10% of neural crest cells in other subregions of the migratory stream show expression profiles that overlap with any of the Trailblazer cells (Figure 5C; Lead (n = 8 out of 100 cells) and Trail (m = 10 out of 110 cells)). These data suggest that neural crest cells with a Trailblazer transcriptional signature are mostly confined to the invasive front of the migratory stream.

To better understand the correlation between our Trailblazer transcriptional signature and cell invasion we performed signaling pathway enrichment analysis using genes differentially expressed by the Trailblazers compared to all other neural crest cells harvested from the embryo (Figure 5D). We present both the developmentally relevant (Figure 5D) and irrelevant (Figure 5—source data 1) signaling pathways. We find that of the genes differentially expressed by Trailblazers and classified into known pathways, 81% were increased and only 19% were decreased (Figure 5D; red and blue numbers respectively). We find increased genes associated with pathways such as WNT/Beta Catenin, tight junction, Rho-family GTPases, EMT, integrin, gap junction, endothelin-1 and BMP signaling (Figure 5D). Consistent with our previous findings that VEGF/neuropilin-1 chemoattraction guides cranial neural crest cell migration (McLennan and Kulesa, 2007; McLennan et al., 2010; McLennan and Kulesa, 2010; McLennan et al., 2015), we find that increased expression of VEGF receptors 1, 2 and 3 as well as neuropilin-1 (NRP1) in the Trailblazers (Figure 5D). Further, 37 genes associated with axonal guidance signaling are expressed at significantly higher levels by the Trailblazers, including members of the BMP, EPH, MMP and WNT families (Figure 5D). Together, these results highlight the upregulation of genes and gene families within the Trailblazers previously identified as functionally critical to neural crest cell migration.

The scRNA-seq analysis validates some of our original 16 gene RT-qPCR trailblazer signature

When we examined the presence of our 16 gene RT-qPCR Trailblazer signature (McLennan et al., 2015) within the scRNA-seq data, we found four genes (AQP1, BAMBI, CTNNB1 and PKP2) enriched in both analyses (Figure 5E). The HHSt15 Lead-like subpopulation that was included in our RT-qPCR Trailblazer subpopulation (defined by manual dissection and not molecular signature) showed an increase of ITGB5 and NEDD9 (Figure 5E). Of the other 16 genes, EPHB1 was enriched in HHSt11 Cluster1 compared to all other migrating cells profiled from in vivo (Figure 5E) and CDH7 and CXCR4 were enriched in HHSt11 Cluster2 (Figure 5E). Differential expression analysis between HHSt11 Cluster1 and Cluster2 showed that CDH11, CXCR4, NEDD9, and TFAP2A and TFAP2B were enriched in HHSt11 Cluster2 compared to Cluster1, but no gene of the 16 RT-qPCR signature genes were enriched in HHSt11 Cluster1 versus Cluster2 (Figure 4—source data 1).

Expression analysis by RNAscope of individual genes associated with the scRNA-seq Trailblazer transcriptional signature are not necessarily restricted to neural crest cells at the invasive front

To visualize the expression of genes associated with the Trailblazer transcriptional signature, we optimized a novel protocol combining multiplexed fluorescence in situ hybridization (RNAscope; Morrison et al., 2017; Gross-Thebing et al., 2014) with immunohistochemistry (IHC) in chick (Figure 6). We took advantage of image analysis tools to quantify the fluorescence signals within subregions or individual neural crest cells by fluorescence spot counting to measure transcript levels (Figure 6). Together, this allowed us to detect mRNA target expression specifically within HNK1-labeled migrating neural crest cells in 3D within the intact embryo (Figure 6 and Video 1).

Figure 6 with 3 supplements see all
Expression and functional perturbations of selected genes enriched in Trailblazers.

(A–M) Typical cranial neural crest cell migratory streams at HHSt13 adjacent to rhombomere 4 (r4) with HNK-1-positive cells color coded by the number of PKP2 (Plakophilin 2; A), Desmin (B), KAZALD1 (Kazal Type Serine Peptidase Inhibitor Domain1; C), Tescalcin (D), Troponin I1 (E), Annexin1 (F), Glypican 3 (G), BDNF (Brain Derived Neurotrophic Factor; H), Nexilin (I), BAMBI (BMP and Activin Membrane Bound Inhibitor; J), Desmoplakin (K), TFPI2 (Tissue Factor Pathway Inhibitor 2; L), and CDH5 (Cadherin 5; M). RNAscope spots detected per cell. (A’–M’) Single cells outlined with expression were selected from each neural crest cell migratory stream and (A”–M”) the number of RNAscope spots per cell volume for a cell with high expression and adjacent neighboring cells is shown in the bar graphs. (N) Schematic representation of morpholino and electroporation procedure. (O) Morpholino transfected migrating cranial neural crest migration at HHST15 (n = 16 control MO embryos, n = 20 EDN1 MO embryos, n = 20 PKP2 MO embryos, n = 20 MMP2 MO embryos, n = 7 GATA5 MO embryos, n = 5 WNT5B MO embryos, n = 10 GPC3 MO embryos, n = 7 ANXA1 MO embryos, n = 22 LUM MO embryos). (P) Correlation plot of distance migrated and expression in trailblazers. Circle size correlates with statistical significance: PKP2 MO, p=0.017, EDN1 MO, p=0.0004, MMP2 MO, p=0.0032, GATA5 MO, p=0.022, WNT5B MO, p=0.0001 (Q) Box plot of the distance migrated of morpholino transfected neural crest cells as a percentage of the distance from the neural tube to the tip of the branchial arch. HH, Hamburger and Hamilton (1951); St, stage; OV, otic vesicle; MO, morpholino. Bar = 20 um (A–G). Bar = 15 um (H–M). Bar = 50 um (O).

https://doi.org/10.7554/eLife.28415.016
Video 1
RNAscope fluorescence signal of Tescalcin expression

A typical 3D confocal image z-stack of a HHSt13 chick embryo highlighting the neural crest cell migratory stream adjacent to rhombomere four with cells labeled for Tescalcin by RNAscope and neural crest cells marked by HNK1 immunohistochemistry. The HNK1 membrane label was used to mark the boundaries of neural crest cells and mask the Tescalcin channel in 3D. The locations of RNAscope label for Tescalcin were detected and assigned to individual cells using the HNK1 boundaries. The cells were then colored by the number of detected Tescalcin RNAscope label spots with red being the highest and blue the lowest. A typical cell is approximately 20–30 um in diameter.

https://doi.org/10.7554/eLife.28415.021

We selected 13 genes for expression investigation based on the following criteria: (1) genes with enriched expression in Trailblazer cells; (2) high percentage of single Trailblazer cells in which the gene was detected and; (3) high level at which the gene was expressed by the Trailblazer cells. We find that none of the 13 genes tested are expressed in all of, or the majority of cells at the invasive front (Figure 6A–M). That is, some genes appear expressed in a few or some but not in all of the cells (~30%) of the (narrow) front migratory edge (Figure 6A–M). For example, Plakophilin2 (PKP2) and Kazal Type Serine Peptidase Inhibitor Domain 1 (KAZALD1) showed heterogeneity of expression between invasive front neighboring cells (Figure 6A–A”, C–C”; Figure 4—source data 1). Further, some of the 13 genes analyzed by RNAscope display salt and pepper patterns that include strong expression in a few cells in the invasive front and similarly strong signal in some cells in the beginning or middle portions of the stream (not at the front edge of migration). We found Desmin (DES) and Tescalcin (TESC) to be expressed within the invasive front (Figure B,D; green shading), but also in migrating neural crest cells within the remainder of the stream (Figure 6B,D; red sharing). Lastly, some of the 13 genes were not found at the invasive front, and instead appear to be expressed more strongly at various positions within the migratory stream. That is, TroponinI1 Slow Skeletal Type (TNNI1), Glypican3, and TFPI2 showed strong signal in some cells in the middle or back portions of the stream (Figure 6E,G,L).

Knockdown of some trailblazer genes showed changes in neural crest cell distance migrated

To begin to understand the function of Trailblazer genes identified by our scRNA-seq, we performed loss-of-function experiments by morpholino transfection into premigratory neural crest cells (Figure 6N–Q and Figure 6—figure supplements 13). We selected 8 genes from distinct molecular families to test, based on the same criteria used for selecting genes for expression analysis, but also for predicted roles in neural crest migration. These genes included, WNT5B (expressed in 91% of Trailblazer cells), PKP2 (57%), EDN1 (57%), MMP2 (96%), GATA5 (77%), GPC3 (62%), ANXA1 (47%) and LUM (79%) (Figure 4—source data 1; Figure 6N–Q). We found that knockdown of 3 of these genes (WNT5B (p=0.0001), PKP2 (p=0.017), EDN1 (p=0.004)) led to a significant reduction (>13%) in neural crest cell distance migrated (Figure 6N–Q). Specifically, we found that WNT5B showed the most significant change in neural crest cell migration (21% reduction in distance migrated) (Figure 6P,Q). Three of these genes did not significantly reduce migration (Figure 6Q) and for the remaining 2 genes, we could not get sufficient morpholino transfection within the migrating neural crest stream and therefore did not analyze (data not shown). Although there were some genes in which high expression in the larger percentage of Trailblazers (>50%) correlated with the significant reduction in neural crest cell distance migrated (WNT5B, GATA5, PKP2), this was not a general feature, indicating that expression does not necessarily correlate with function.

Splice blocking morpholinos were used so that their effects could be detected by RT-PCR. For example, PKP2 MO and ANXA1 MO both resulted in exclusion of an exon (3 and 4 respectively), which could be seen by using PCR primers that flank the effected exon (Figure 6—figure supplement 2). Other morpholinos were harder to test due to small number of introns/exons (for example LUM), or frame shifts resulting in unknown changes (for example DLL, GATA5). Overall stream migration distance was also calculated using the HNK1 neural crest marker, and although some morpholino transfections caused a statistically significantly shorter stream, they all still migrate over 90% of the distance to the end of the branchial arch (Figure 6—figure supplement 3). This suggests that the morpholino transfections are cell autonomous and that wildtype neural crest cells are able to presume the role of the Trailblazers and therefore migrate into the target site.

Identification of two distinct transcriptional signatures associated with neural crest cells isolated shortly after dorsal neural tube exit

We described above that hierarchical clustering and PCA analyses of scRNA-seq data from neural crest cells isolated from the invasive front revealed two cell subpopulations with distinct gene expression profiles at HHSt11, termed Cluster1 and Cluster2 (Figure 4A–C). Neural crest cells within Cluster1 and Cluster2 are not identified as Trailblazers, but the distinct gene expression profiles hinted at further examination. Specifically, we find that neural crest cells isolated shortly after dorsal neural tube exit (Cluster1) were enriched for genes known to be associated with early neural crest inductive signals and the neural plate border gene regulatory module previously described (Simões-Costa and Bronner, 2015) (Figure 7A and Figure 4—source data 1). These genes include WNT3, WNT3A, WNT7B, WNT8B, ZIC1, ZIC3 MSX1, CDH6, EphA7, and CDH2, some of which are highlighted in a volcano plot (Figure 7A).

Newly emerging neural crest segregate into discrete groups representing a transition from EMT to directed migration.

(A) Volcano plot of all genes statistically enriched or reduced in St11 Cluster1 (n = 20 cells, 381 genes). (B) Volcano plot of all genes statistically enriched or reduced in St11 Cluster2 (n = 55 cells, 1040 genes). (C) GO Term analysis of genes differentially enriched by HHSt11 cluster1. No GO terms were reduced. (D) GO Term analysis of genes differentially reduced by HHSt11 cluster2. No GO terms were enriched. (E) Significant signaling pathways predicted to be enriched (red) or reduced (blue) in St11 Cluster1. Axonal guidance signaling is highlighted. (F) Significant signaling pathways predicted to be enriched (red) or reduced (blue) in St11 Cluster2. Axonal guidance signaling is highlighted. GO terms with a p-value<0.05 and a q-value <0.1 were considered significant. Significant GO terms were considered using the simplify function with an adjusted p-value cutoff of 0.05.

https://doi.org/10.7554/eLife.28415.022

In contrast, other known neural crest border plate genes (e.g., TFAP2B, TFAP2C, TFAP2E and DLX5) and neural crest specification and migration genes (e.g., SNAI1, TWIST1, SOX10, ETS1, RXRG and EBF1, 2 and 3) were reduced in Cluster1 (Figure 7A and Figure 4—source data 1). The transcription factor PAX6, which was not previously identified in the bulk RNA-seq analysis and neural crest gene regulatory network (Simões-Costa and Bronner, 2015) was enriched in HHSt11 Cluster1 cells, and is expressed in 85% of Cluster1 cells (Figure 7A). The HHSt11 Cluster2 subpopulation displays increased expression of genes related to later neural crest cell specification, epithelial-to-mesenchymal transition (EMT) and migration (Figure 7B). These genes include PAX3, TFAP2B, SOX8 and SOX10 (Figure 7B and Figure 4—source data 1). We also uncover additional transcription factors SOX8 (expressed in all Cluster2 cells) and FHL2 (expressed in 70% of Cluster2 cells) that were enriched in the HHSt11 Cluster2 subpopulation (Figure 7B). Furthermore, HHSt11 Cluster2 cells also showed decreased expression of factors thought to describe the neural crest gene regulatory network (e.g., WNT2B, WNT3A, WNT4, WNT5B, WNT7B, WNT9A, BMPR1B, BMP2, BMP10, FGF3, FGF8, FGF10, ZIC1, TFAP2C, DLX1, DLX3, DLX6, SNAI1, ID1 and ID3) (Figure 7B and Figure 4—source data 1). A gene ontology (GO) analysis was compiled using the R package cluster Profiler (3.2.6). GO term analysis reveals distinct molecular characteristics of cells within each cluster (Figure 7C–D). Neural crest cells within Cluster1 display genes enriched in signaling, adhesion and cell guidance, with no depleted GO terms (Figure 7C). Cluster2 identifies with neural crest cells that display genes depleted in signaling and transcription regulation, with no enriched GO terms (Figure 7D). These results more comprehensively describe the molecular characteristics of neural crest cells shortly after dorsal neural tube exit and suggest a dynamic transition from neural tube exit to acquisition of directed migration.

To identify signaling pathways enriched within HHSt11 clusters, we performed pathway enrichment analyses for the lists of genes differentially expressed by HHSt11 Cluster1 and Cluster2 versus all other in vivo NC cells to (Figure 7E,F). For both HHSt11 clusters, we found a mix of up- and down-regulated genes within these pathways, but overall observed increased gene expression associated with signaling pathways in Cluster1 and decreased gene expression associated with signaling pathways in Cluster2 (Figure 7E,F). Specifically, in Cluster1, 59% of differentially expressed genes associated with signaling pathways were up-regulated compared to only 1% in Cluster2. Interestingly, based on the gene expression, pathways such as WNT/Beta Catenin signaling, EMT, PCP pathway and axonal guidance signaling are increased in Cluster1 but then decreased in Cluster2 (Figure 7E,F). When focusing just on axonal guidance signaling, ADAMs, EPHs and WNTs are increased in Cluster1 but then decreased in Cluster2. The presence and absence of signaling pathways are also relevant. For example, Eph/ephrin signaling is significant in Cluster1 but not in Cluster2 (and therefore does not make the Cluster2 list) (Figure 7E,F). Actin cytoskeleton signaling is significant in Cluster2, with 21 out of 22 genes showing decreased expression, but this pathway is not significant in Cluster1 (Figure 7E,F). These results highlight the dramatic differences between Cluster1 and 2 at HHSt11 and suggest that their roles in EMT and migration may be very different.

The transcriptional signatures of lead and trail subpopulations are defined by small numbers of differentially expressed genes

To define the genes and pathways associated with Lead and Trail subpopulations, we analyzed expression differences unique to each Lead and Trail subpopulation relative to all other in vivo cells isolated. Lead and Trail subpopulations were the largest subpopulations isolated (25% (Lead) and 70% (Trail) of the migratory stream, respectively), and exhibit high degrees of heterogeneity (Figure 3 and Figure 3—figure supplement 1; Figure 8—source data 1). First, at HHSt13, we identified 45 enriched and 151 reduced genes in Trailers and 26 enriched and 35 reduced genes in Leaders (Figure 8A,B; Figure 8—source data 1). AQP5, ALOX5AP, and SS5R were enriched and NKX2-6, NRG1 and UNC5B were reduced in Trailers (Figure 8A). PKP1, PKDREJ, and CWH43 were enriched and MSX1 and myosin reduced in Leaders (Figure 8A; Figure 8—source data 1). Strikingly, signaling pathway enriched analysis revealed reduced expression of genes associated with known developmental pathways (Figure 5—source data 1). These included VEGF signaling as well as genes associated with actin cytoskeleton and integrin signaling (Figure 8C).

Lead and trailing subpopulations exhibit high degrees of heterogeneity.

(A) Volcano plot of all genes statistically enriched or reduced in St13 Trail (n = 60 cells, 196 genes). (B) Volcano plot of all genes statistically enriched or reduced in St13 Lead (n = 49 cells, 61 genes). (C) Significant signaling pathways predicted to be enriched (red) or reduced (blue) in St13 Trail and Lead. (D) Volcano plot of all genes statistically enriched or reduced in St15 Trail (47 cells, n = 306 genes). (E) Volcano plot of all genes statistically enriched or reduced in St15 Lead (n = 46 cells, 215 genes). (F) Significant signaling pathways predicted to be enriched (red) or reduced (blue) in St15 Trail and Lead.

https://doi.org/10.7554/eLife.28415.023

Second, at HHSt15, we identified 174 enriched and 132 reduced genes in Trailers (Figure 8D; Figure 8—source data 1). Enriched genes in HHSt15 Trailers included GATA5, HOXA3 and NOG; reduced genes included FGF16, NEFM and NEUROD4 (Figure 8D; Figure 8—source data 1). Genes associated with developmental pathways such as Rho GTPases, RHOA, CXCR4, BMP and Actin cytoskeleton signaling were reduced in HHSt15 Trailers (Figure 8F; Figure 5—source data 1). HHSt15 Leaders had 77 enriched genes, including FGFs, NKX2-6 and SHH, as well as 138 reduced genes, including CDH7, PAX7 and SOX10 (Figure 8E; Figure 8—source data 1). From these gene lists, signaling pathway enrichment analysis showed reduced activity in axonal guidance, ephrin receptor and netrin signaling pathways (Figure 8F; Figure 5—source data 1).

Discussion

Cell migration in the embryo positions cells into precise target tissues at appropriate axial levels to ensure proper patterning and organ development. Directed motion and collective movement in groups are characteristics of embryonic cell migration and differ from processes where precursor cells are born and in the same location undergo shape changes and sorting to give rise to tissue structures. To study cell migration during embryonic development, we used the highly invasive neural crest as a model. We focused on cranial neural crest cell migration and exploited the strengths of the avian embryo to perform manual cell isolation at progressive developmental stages, bulk and single-cell RNA-seq to determine gene expression profiles within cells of a typical migratory stream over time. Our goal was to identify novel transcriptional signatures and molecular transitions associated with subpopulations of cells within a cranial neural crest cell migratory stream and determine how these transcriptional signatures changed as cells traveled through different microenvironments. This information has the potential to lead us to a better understanding of the molecular features underlying observed complex neural crest cell behaviors and mechanistically help to explain directed and collective cell migration processes.

We found that four features contribute to the precise directed and collective migration of the cranial neural crest to precise head targets. First, the majority of neural crest cells within the invasive front change their gene expression in a consistent manner during migration; however, a small subpopulation of Trailblazer cells narrowly confined to the invasive front have a conserved transcriptional signature throughout migration. Second, cranial neural crest cells newly exited from the dorsal neural tube cluster into two subpopulations with profiles distinct from the Trailblazer signature and indicative of a rapid switch from epithelial-to-mesenchymal transition to directed migration. Third, Lead neural crest cells (near the front of the stream, but not narrowly confined to the invasive front) cluster more closely with cells within the invasive front that do not display a Trailblazer transcriptional signature. Fourth, Trailing neural crest cells display a gene expression profile distinct from Lead and Trailblazer cells and cluster more closely with newly exited cells and in vitro migrating cells, the latter of which have a molecular signature distinct from any other in vivo profile.

The significant molecular distinction between cranial neural crest cells isolated from the invasive front and compared to other cells throughout the stream during migration (HHSt13 and HHSt15) suggests cells at the invasive front rapidly adapt to changes in the local microenvironments through which cells travel. Heatmaps of both the bulk RNA-seq 3X replicates (Figure 1B) and the mean values (Figure 2—figure supplement 2) showed that there are distinct clusters of genes expressed at significantly higher or lower levels in the invasive front when compared to cells within the remainder of the stream. The number of genes up-regulated in cells within the invasive front at HHSt13 (n = 61) and HHSt15 (m = 55) was consistent with an overlap that included 23 genes (Figure 1). In contrast, greater than 50% (255 out of 477) of the differentially expressed genes on the heatmaps (Figure 1B) at HHSt13 but less than 8% (36 out of 477) at HHSt15 were down-regulated in cells within the invasive front, with an overlap of 19 genes (Figure 1B, blue bands). The nearly 4-fold difference in the number of genes down-regulated versus up-regulated in cells within the invasive front at HHSt13 confirmed rapid changes take place as neural crest cells encounter the microenvironment of the branchial arch entrance. Yet, the 42 genes (23 up-regulated and 19 down-regulated) out of the 477 differentially expressed genes represented a conserved molecular signature at the two developmental stages and level of consistent gene expression during this microenvironmental transition.

In vivo expression analysis using RNAscope and immunolabeling identified heterogeneous expression within neighboring Trailblazers, but could not confirm gene expression unique to cells within the invasive front suggesting an inconsistency between the Trailblazer transcriptional signature and the in vivo mRNA expression patterns. While we expected the RNAscope expression analysis to identify more Trailblazer cells than were detected we found that none of the 13 genes tested were expressed in all of, or the majority of the cells at the invasive front. Instead, RNAscope analysis showed some genes either appeared in a few Trailblazer cells or displayed salt and pepper patterns that included strong expression in a few cells in the middle or back portion of the migratory stream (Figure 6). Further, some genes were not found expressed at the front edge at all (Figure 6). This could be due to differences in sensitivity between scRNAseq and RNAscope technologies. That is, the chemistry of scRNAseq may detect more transcripts than the newly established RNAscope recently optimized for avian embryos (Morrison et al., 2017). Further, the low number of cells within the invasive front that we detected to have a Trailblazer signature (for example, at HHSt13 there are approximately 10–15 Trailblazer cells out of 500 total cells in the migratory stream) make it difficult to visualize such a low number of expressing cells. In addition, the 3D nature of the cranial neural crest cell migratory stream and depth into the embryo make it difficult to more accurately assess fluorescence transcript levels in trailblazer cells in comparison to follower cells located more dorsal and accessible to signal detection by confocal microscopy. Future studies aided by technological advances that increase resolution may better link scRNA-seq data with mRNA expression detection.

Several genes discovered within the scRNA-seq Trailblazer transcriptional signature and previously identified in our RT-qPCR analysis have been implicated in human cancer cell invasion. These genes included AQP1, BAMBI, CTTNB1, and PKP2 (Figure 5; McLennan et al., 2015). BAMBI, a transforming growth factor beta (TGFB) inhibitor is up-regulated in metastatic carcinomas including colorectal and osteosarcoma (Fritzmann et al., 2009; Zhou et al., 2013). Similarly, the expression of PKP2, a desmosomal adhesions protein has been implicated in glioma and bladder cancer metastasis (Takahashi et al., 2012; Zhang et al., 2017) and high expression of Tescalcin, a gene thought to play a role in cell differentiation and implicated in invasive colorectal cancer (Kang et al., 2016). Furthermore, WNT5B, a gene identified within our scRNA-seq Trailblazers has been implicated in cancer aggressiveness in a number of cancer cell lines and tissues (Harada et al., 2017). When we knocked down PKP2 and WNT5B in premigratory neural crest cells by morpholino we found significant reduction in cell distance migrated (Figure 6N–Q; WNT5B (21% reduction, p=0.0001), PKP2 (13% reduction, p=0.017). Thus, genes identified and associated with some of the most invasive cranial neural crest cells may have a critical functional role in migration and may be part of a broader molecular signature in cancer aggressiveness.

Pathway analysis of the Trailblazer transcriptional signature identified genes involved in VEGF signaling, strengthening support for our previous discovery that VEGF chemoattraction is critical for neural crest migration (McLennan and Kulesa, 2007; McLennan and Kulesa, 2010; McLennan et al., 2010). We found that all three VEGF receptors as well as neuropilin-1 (NRP1) were expressed at significantly higher levels in the Trailblazer neural crest cells (Figure 5D). When we previously knocked down neuropilin-1 in chick cranial neural crest cells, cells were able to migrate away from the neural tube but failed to invade the second and rostral portion of the third branchial arches (McLennan and Kulesa, 2007; McLennan and Kulesa, 2010; McLennan et al., 2010). These data provide strong evidence to support a critical role of the Trailblazer neural crest cells to guide invasion into the branchial arches.

Axonal guidance signaling in general has provided insights that have led to discoveries of mechanisms of neural crest migration in the head (Osborne et al., 2005; Gammill et al., 2007), eye (Lwigale and Bronner-Fraser, 2009) and gut (Anderson et al., 2007). Furthermore, there is exciting emerging evidence that genes implicated in nervous system development have functional similarities with vascular patterning, including endothelial tip cell invasion (Wälchli et al., 2015). When we performed pathway analysis, we found a significant number of axonal growth cone signaling genes (n = 26) identified within the Trailblazer transcriptional signature (Figure 5D) and early emigrating neural crest cells (Figure 8) indicative of shared molecular features between the two invasion phenomena. However, our investigation of the neural crest literature showed that only 7 out of the 26 genes we found enriched in the Trailblazer cells have been examined in neural crest cell migration, including BMP2, BMP4, MMP2, RAC2, SEMA3C, WNT1, and WNT11 (Anderson, 2010; Correia et al., 2007; de Melker et al., 2004; Goldstein et al., 2005; Matthews et al., 2008; Tang et al., 2016; Toyofuku et al., 2008). Thus, determination of the function of the other 19 out of 26 axonal growth signaling genes in neural crest migration may lead to exciting insights and further our knowledge of common developmental principles of neurovascular patterning (Figure 5D).

The scRNA-seq analysis validated some but not all of the genes identified by our 16 gene RT-qPCR Trailblazer signature (McLennan et al., 2015), which is understandable given the key variances between the two analyses. Specifically, the RT-qPCR Trailblazer signature was defined by cells narrowly confined within the invasive front in HHSt13 and 15 embryos. Our scRNA-seq data defined Trailblazers as a molecularly distinct subpopulation within the invasive front at HHSt11, 13 and 15. Second, RT-qPCR and scRNA-seq are distinct methodologies. Notably, RT-qPCR is the gold standard of sensitivity for detecting lowly abundant transcripts, whereas scRNA-seq is capable of uncovering many more genes within a single analysis. Thus, our new scRNA-seq analysis added information from HHSt11 cells to refine and expand the Trailblazer signature based upon unbiased clustering of a comprehensive transcriptome analysis.

Newly emerging neural crest cells clustered as two distinct subpopulations with transcriptional signatures unique from that of the Trailblazer cells suggesting a rapid molecular transition to directed migration shortly after neural tube exit. At HHSt11, there are few neural crest cells that have delaminated from the dorsal neural tube at the axial level of r4 such that we grouped together all cells that were isolated (Figure 2). Hierarchical clustering identified the Trailblazer signature and two other transcriptional signatures that we associated with Cluster1 and Cluster2 (Figure 4). Genes within the transcriptional signature of Cluster1 closely resembled those related to the premigratory and delamination events of the neural crest (Figure 7). In contrast, genes within the transcriptional signature of Cluster2 more closely resembled those related to directed neural crest cell migration (Figure 7). Interestingly, we found very little overlap between the Cluster1 and Cluster2 transcriptional signatures that was indicative of a rapid switch in the gene expression profile rather than continuous progression of gene expression changes (Figure 7). Further analysis of mRNA and protein expression of multiple genes using RNAscope as shown above or in other emerging strategies (Choi et al., 2016; Roellig et al., 2017) will help to shed light on the spatial positions of cells with these unique transcriptional signatures identified by scRNA-seq.

The significant differences in gene expression profiles comparing in vivo transcriptional signatures with in vitro data was not surprising given the absence of the neural crest microenvironment and lack of formation of discrete migratory streams in culture. PCA analysis highlighted the poor correlation in gene expression patterns (Figure 2C,D) with the closest subpopulation to the in vitro neural crest cells being a subpopulation of the newly exited neural crest cells at HHSt11 (Figure 2C). This made logical sense since both in vitro and newly exited in vivo neural crest cells are within the range of influence from signals within the neural tube (Figure 2C). Furthermore, no GO terms were upregulated in in vitro neural crest cells and instead there was a down-regulation of terms such as development, cell signaling, migration and proliferation (Figure 2E). Furthermore, our previous RT-qPCR study showed that migrating chick cranial neural crest cells did not display distinct Lead and Trailer gene expression profiles in culture (McLennan et al., 2015). However, exposure to VEGF in vitro resulted in the upregulation of a small subset of genes associated with an in vivo lead cell signature (McLennan et al., 2015). Together, this evidence supported our model that the in vivo embryonic neural crest cell microenvironment critically influences transcriptional signatures of migrating neural crest cells.

Advances in cell isolation and single-cell RNA-seq have proven invaluable for elucidating cell-type specific transcriptional signatures, often defined by an ensemble of genes rather than high expression of any single gene (Dueck et al., 2016Tasic et al., 2016), during single time-point embryonic and adult neurogenic events (Chu et al., 2016; DeLaughter et al., 2016; Nelson et al., 2016; Scholz et al., 2016; Tirosh et al., 2016a, 2016b; Dulken et al., 2017; Litzenburger et al., 2017). However, the identification of distinct transcriptional signatures associated with dynamic cell populations during cell migratory events has remained challenging and unclear. By using single-cell transcriptome analysis applied to embryonic cranial neural crest cell migration at three progressive developmental stages, we provide the first comprehensive analysis of neural crest migration. Our data identify and establish hierarchical relationships between cell position- and time-specific transcriptional profiles, including a Trailblazer transcriptional signature characterized by a large number of differentially expressed genes within a small subpopulation of cells at the invasive front. Functional in vivo knockdown of a subset of individual Trailblazer genes showed significant but modest changes in neural crest cell distance migrated. However, expression analysis by RNAscope multiplex fluorescence in situ hybridization could not definitively confirm expression restriction to the invasive front. The present study provides new insights into the molecular patterns that underlie directed, collective cell migration of the neural crest and reveals individual targets and signaling pathways for future mechanistic studies and to other neural crest migratory streams throughout the body and other cell invasion phenomena in cancer, wound healing, and the immune system response.

Materials and methods

In Ovo electroporation

Fertilized, white leghorn chicken eggs (NCBI Taxonomy ID:9031; Centurion Poultry Inc., Lexington, GA, USA) were incubated in a humidified chamber at 38C to the developmental desired stage. For NC cell isolations, pre-migratory NC were labeled by electroporation of a Gap43-YFP plasmid DNA into the dorsal neural tube of HHSt9 embryos as previously described (McLennan and Kulesa, 2007) and re-incubated (8, 16 and 24 hr for HHSt11, 13 and 15, respectively). After re-incubation, embryos were harvested into chilled 0.1% DEPC PBS and screened for health and transfection efficiency before tissue isolations were performed.

For functional analysis, H2B mCherry (2.5 μg/μl) and fluorescein-tagged morpholinos (0.5 mM) (Gene Tools, Philomath, OR, USA) were injected into the neural tube and electroporated at HHSt8 (Figure 6N). The embryos were re-incubated to HHSt15. After re-incubation, embryos were harvested and fixed in 4% paraformaldehyde at room temperature for 2 hr or at 4C overnight. Immunohistochemistry was then performed on the embryos for HNK1 (neural crest cell marker) using a standard protocol (McLennan and Kulesa, 2007) before being imaged on a LSM 800, using a Plan-Apochromat 10x/0.45 M27 objective (Zeiss, Oberkochen, Germany) (Figure 6—figure supplement 3).

For splice-blocking morpholino validation by RT-PCR, the chicken LMH cell line was transfected with either control or experimental splice blocking morpholino using Endo-porter (Gene Tools, Philomath, OR). 48 hr after transfection, cell lysate was harvested and total RNA isolated (RTN70 Sigma, St. Louis, MO). RNA was converted to cDNA (4368814, Thermo Fisher Scientific, Waltham, MA) for RT-PCR (M7122, Promega, Madison, WI, USA) using primers flanking the region affected by each morpholino to demonstrate perturbation of translation (Figure 6—source data 1). PCR products were analyzed on a LabChip GX Touch (PerkinElmer, Waltham, MA) (Figure 6—figure supplement 2).

Tissue isolations

All manual dissections of Gap43-YFP transfected neural crest cells were performed using a fluorescent dissecting microscopy. For bulk RNA-seq samples, the NC streams adjacent to rhombomere 4 (r4) were isolated by manual dissection from HHSt13 and 15 embryos and divided into three biological replicates. The streams were subsequently cut into two spatially distinct, non-overlapping portions: the most invasive front (5% of the stream, termed ‘Front’) and the remaining proximal-most trailing portion of the stream (95% of the stream, termed ‘Stream’). Developmental stage- and stream position-matched pools of tissue were dissociated as previously described (Morrison et al., 2015).

For scRNA-seq, the NC r4 stream was isolated by manual dissection at HHSt11, 13 and 15. HHSt11 r4 streams containing few migrating NC were pooled. For HHSt13 and 15 embryos, the r4 NC stream was further subdivided by manual dissection into three spatially distinct, non-overlapping portions: the most invasive front (5% of the stream, termed ‘Front’), the leading portion of the stream (25% of stream proximal of the most invasive front, termed ‘Lead’) and the proximal-most trailing portion of the stream (70% of the stream, termed ‘Trail’). Developmental stage- and stream position-matched tissues were pooled and dissociated before cytometric isolation as previously described (Morrison et al., 2015). Optimizing tissue dissociations and limiting processing time from embryo to lysis are critical for maintaining cell health and generating high quality scRNAseq data (Morrison et al., 2015).

In Vitro neural crest cell cultures

Cranial neural tubes (NT; r3-r5) were removed from embryos at HHSt9, plated as previously described (Krull et al., 1997). NT were incubated for 24 hr to allow NC to emerge from the NT and migrate onto the substrate. NT were then gently removed with a tungsten needle, so that only migratory NC remained. Media was gently removed and the adherent cells rinsed once with 0.1% DEPC PBS before trypsinization (25200–056, Thermo Fisher Scientific, Waltham, MA) for 5 min at 37C and trypsin inactivation with fetal bovine serum. The cells were then pelleted and resuspended in chilled 0.1% DEPC PBS.

Cytometry

Cells were isolated by FACS, which included forward scatter, side scatter, pulse width, live/dead stain and YFP gates as previously described (Morrison et al., 2015). Cells were sorted directly into 6 ul of Clontech lysis solution containing 0.05% RNAse inhibitor. Following lysis for 5 min at room temperature, lysates were immediately frozen on dry ice and stored at −80C.

Bulk RNA-seq

Bulk RNA-seq lysates were thawed on ice. cDNA synthesis and library preparation were performed with SMART-seq v4 Ultra Low Input RNA-seq (634892, Takara, Kusatsu, Shiga, Japan) and Nextera XT DNA dual indexing library preparation kits as recommended by the manufacturer (FC-131–1096, Illumina, San Diego, CA). All 12 libraries (4 samples with 3 biological replicates each) were pooled and run across 2 lanes of an HiSeq 2500 generating 50 bp single reads (Illumina). Quality control analyses, including reads per samples, coverage by normalized position along transcript and correlation among all samples and biological replicates are shown in Figure 1—figure supplement 1.

Single-cell RNA-seq

Single-cell lysates were thawed on ice. cDNA synthesis and library preparation were performed with SMART-seq v4 Ultra Low Input RNA-seq (634892, Takara, Kusatsu, Shiga, Japan) and Nextera XT DNA dual indexing library preparation kits as recommended by the manufacturer (FC-131–1096, Illumina). Despite more cost-effective scRNAseq strategies, SMARTseq v4 was employed in order to detect the greatest number of genes per single cell, which is critical for distinguishing highly similar cell subpopulations (Ziegenhain et al., 2017; Dueck et al., 2016). All single-cell cDNAs and libraries were first quantified on a Qubit. Representative single-cell cDNAs and libraries were also quantified on a Bioanalyzer 2100 High Sensitivity DNA Assay (Agilent, Santa Clara, CA). A maximum of 96 single-cell libraries were pooled per lane of a HiSeq 2500 generating 50 bp single reads (Illumina). Quality control analyses for scRNA-seq are shown in Figure 2—figure supplement 1.

Computational analysis

Fastq files were mapped to the chicken genome galGal4 from UCSC using Tophat (2.0.13) with options -x 1 g 1. Ensembl 80 annotations were used to define gene coordinates. Quality of the samples was assessed using FastQC (0.10.1). The R environment (RRID:SCR_001905) was used for the statistical analysis of the data. All annotated code is located in Supplementary file 1. RPKM values generated and differentially expressed genes were found using edgeR. Genes were considered differentially expressed if they had an Benjamini-Hochberg adjusted p-value less than 0.05 and a log2 fold change greater than 2 or less than −2. Hierarchical clustering was performed using Pearson correlations. Clustering was achieved using base principal component functions in R. Gene ontology (GO) analysis was completed using ClusterProfiler (3.2.6). GO terms with a p value less than 0.05 and a q value below 0.1 were considered significant. Significant GO terms were condensed using the simplify function with an adjusted p value cutoff of 0.05. Pathway enrichment analyses were performed using Qiagen’s Ingenuity Pathway Analysis (RRID:SCR_008653; v31813283, Hilden, Germany).

RNAscope ISH, Immunohistochemistry and Image Analysis

RNAscope multiplex fluorescent ISH (Advanced Cell Diagnostics, Newark, CA) combined with IHC in the whole avian embryo was adapted from Gross-Thebing et al.’s work in the zebrafish (Morrison et al., 2017; Gross-Thebing et al., 2014). Major optimizations for success in the avian included fixation time, avoiding air drying the embryos, volumes and incubation times. IHC to label migratory NC in the whole embryo was performed as described previously (McLennan and Kulesa, 2007). Images were collected on a LSM 800 confocal microscope. Images were analyzed with the Cell module in Imaris (RRID:SCR_007370; Bitplane, Belfast, Northern Ireland). Using the HNK1 as a neural crest membrane marker (RRID:AB_10013722), cells were segmented and then spots of RNAscope signal detected and assigned to each cell. Smoothing and background subtraction aided in segmentation of the cells. Manual modifications to the cell borders was completed in some cases where Imaris failed to segment properly using the Surface module and marching cubes detection. Cell surfaces were colored to represent the number of detected transcripts.

For morpholino transfected fixed embryos at HHSt15, a sharpened tungsten needle was used to remove the trunk portion of the embryo and cut the cranial embryo in half down the midline. Embryos were then mounted into glass slides as previously described (Teddy and Kulesa, 2004) and imaged by confocal microscopy on a LSM 800 (Zeiss). The percentage of distance the morpholino transfected cells migrated was calculated by measuring the total distance of the ba2 and distance the morpholino transfected cells migrated.

References

  1. 1
  2. 2
  3. 3
  4. 4
  5. 5
  6. 6
  7. 7
  8. 8
  9. 9
  10. 10
  11. 11
  12. 12
  13. 13
  14. 14
  15. 15
  16. 16
  17. 17
  18. 18
  19. 19
  20. 20
  21. 21
  22. 22
  23. 23
  24. 24
  25. 25
  26. 26
  27. 27
  28. 28
  29. 29
  30. 30
  31. 31
  32. 32
  33. 33
  34. 34
  35. 35
  36. 36
  37. 37
  38. 38
  39. 39
  40. 40
  41. 41
  42. 42
  43. 43
  44. 44
  45. 45
  46. 46
  47. 47
  48. 48
  49. 49
  50. 50

Decision letter

  1. Marianne Bronner
    Reviewing Editor; California Institute of Technology, United States

In the interests of transparency, eLife includes the editorial decision letter and accompanying author responses. A lightly edited version of the letter sent to the authors after peer review is shown, indicating the most substantive concerns; minor comments are not usually included.

[Editors’ note: a previous version of this study was rejected after peer review, but the authors submitted for reconsideration. The first decision letter after peer review is shown below.]

Thank you for choosing to send your work, "Single-cell transcriptome analysis of neural crest migration reveals signatures of invasion and molecular transitions", for consideration at eLife. Your initial submission has been assessed by a Senior Editor, who is also a member of our Board of Reviewing Editors, and two reviewers. Although the work is of interest, we regret to inform you that the findings at this stage are too preliminary for further consideration at eLife.

Specifically, the reviewers think you should confirm expression of some of trailblazer genes in embryos and perform further in vivo tests for a few more candidate genes. In addition, the revised manuscript should compare the new data with your results previously published in 2015, to establish if single cell RNA-seq approach agrees with or changes our perspective. It would be important to emphasize what is new. We would be very interested in receiving a new version of this paper after addressing all of the reviewers' comments and every effort would be made to return the paper to the same reviewers. I have included the full reviews which I hope will be helpful in revising the manuscript.

Reviewer #1:

This study, titled "Single-cell transcriptome analysis of neural crest migration reveals signatures of invasion and molecular transitions" submitted by Jason Morrison et al., is a follow up on their group study in 2015 (McLennan et al., 2015 Development 142:2014-25), uses transcriptomic analyses to dissect the molecular signatures of the known and novel subpopulations within the avian neural crest migratory streams. The study identifies molecular signatures for the 'trailblazer', a subpopulation of cells at the invasive front of the migratory stream, and two distinct subpopulations of neural crest cells exiting the dorsal neural tube. The study partially confirms previous findings in terms of the trailblazer molecular signatures but also raises questions regarding the experimental designs and analyses that should be performed to engage in an unbiased manner the transcriptomic heterogeneity of the migratory neural crest cells (see comments below). The reviewer found a scarcity of functional assays to explore the biological meaning of cellular heterogeneity in neural crest (only one supplemental figure described in the Discussion section, which perhaps should be brought to the main portion of the manuscript). It seems fundamental to verify the RNA seq data sets by independent assays, e.g., PCR or in situ, on key genes/pathways that had been identified in the bioinformatics analyses, as well as with the inclusion of known markers for each relevant population. The reviewer also suggests removal of some figure panels that appear to be irrelevant to the study to streamline the flow of the manuscript (see comments below).

Results section 1

• Why were these two stages being chosen? Stage 13 appears to be very close to stage 15 in terms of gene expression. Bioinformatics analyses also did not appear to justify there was a stage specific difference.

• At what stage did the authors label the premigratory neural crest cells?

• There was no description whether the 16 genes identified in McLennan et al., 2015 were found in the bulk RNA seq and what were the new players. For example the analyses should have been performed on genes such as, Hand2/GPC3 as invasive front signature, PDGFRL/CXCR1 as stream signature, and EPHA4 and FOXD3 as trailing signature identified in McLennan et al. 2015. These need to be first established in the bulk RNA seq analysis.

Results section 2

• "we isolated and profiled individual cells from different stream positions at three developmental stages". “We subdivided neural crest cell migratory streams at HHSt13 and HHSt15 into three spatially distinct, non-overlapping subregions by careful manual dissection and single cells were isolated by FACS". The authors need to justify the rationales behind the selection of stages (two versus three) and dissection methods (front/trail vs. front/lead/trail) for the bulk and single cell RNA seq. Were there other previous experiments and studies that the authors or others had performed that would help justify the differences?

• "Hierarchical clustering of single neural crest cell averages….." The initial analysis should have been performed without averaging the single cell transcriptome values according to their migration positions. The current experimental design, by clustering single cell samples, could have unavoidably biased the results. The conclusion reached also did not add any new insight into the study. Figure 2C could be removed. The aim of the experiment was to identify new genes/signatures that can define the different subpopulations and the authors could have used the previous 16 genes (McLennan et al., 2015) or other known targets from the bulk RNA seq in Figure 1 or other studies as a guide to distinguish the front cells.

Results section 3

• "in vitro transcriptional signatures are distinct from any in vivo profile". This section appears not relevant to the rest of the study. This section does not give any new insight into the cellular heterogeneity of migrating neural crest cells.

Rest of the Results sections

• “Our previous RT-qPCR (McLennan et al., 2015) and bulk RNA-seq presented above (Figure 1) show that cranial neural crest cells at the invasive front of a stream have aspects of their transcriptome that are conserved during migration.” This rationale does not fit into the experimental design for Figure 5. If the invasive front were to have a conserved transcriptome during migration, then all cells should have been clustered for analysis. The experiment presented in Figure 5 instead appears to try to identify heterogeneity of the front cells dissected from the embryos.

Figure 5: Where did the front cells for stage 11 come from? It was stated in the first paragraph of the subsection “Single-cell RNA-seq identifies gene expression variances based upon spatial position within the neural crest cell stream and temporal progression along the migratory pathway” that there was no subdivision of stage 11 cells.

• Use different colors for Trail and Trailblazers in Figure 6C.

• Subsection “The scRNA-seq Trailblazer signature is uniquely defined and associated with neural crest cells confined to the invasive front of the stream”, last paragraph: It is interesting that the authors were able to identify previously known pathways, such as the VEGF signaling, that were associated with Trailblazer identity. However, as the purpose was to discover new pathways, an unbiased search should be performed instead of a focus on developmental relevant signaling pathways. As noted, a significantly higher number of pathways were found to be upregulated and many of the p-values were relatively high (Figure 6D) suggesting that there might exist other "developmentally irrelevant" pathways that were associated with trailblazer signatures.

• Please highlight the p-value in Figure 6D that represent significance difference.

Figure 6D can be better presented with the p-values align according to their values from top to bottom.

• What were the criteria to select the gene list for presentation in Figure 6B? Were these genes most highly enriched or depleted?

• Highlight in the text and figure legend for Figure 6B that those genes highlighted in red were previously found in McLennan et al., 2015.

• The conclusion "These results highlight the upregulation of many genes within the trailblazer subpopulation, supporting the hypothesis that trailblazers are a dramatically active subpopulation of cells at the most invasive front during neural crest migration" needs to be substantiated (and if possible reworded to omit superlatives.

• " We also found that CDH7, CXCR4 and EPHB1 were upregulated in the HHSt11 single cells (Figure 6E)." What were these single cells? Were there any bias in stage 11 clusters 1 and 2 identified earlier?

• Subsection “The scRNA-seq analysis validates some of our original 16 gene RT-qPCR Trailblazer signature”: the authors said there was "no significant enrichment" for a number of genes found in previous studies meaning the expression of these genes were upregulated but not statistically significant or were downregulated instead? Following up on this, if these genes were not upregulated in the trailblazer population, were they found to be upregulated in any of the subpopulations of the migrating neural crest cells?

• Subsection “Identification of two distinct transcriptional signatures associated with neural crest cells isolated shortly after dorsal neural tube exit”, second paragraph: How was the GO search performed for clusters 1 and 2? Why were there no depleted and enriched terms for clusters 1 and 2 respectively?

• Please correct typos. “To identify signaling pathways enriched within HHSt11 clusters, we performed pathway enrichment analyses "for" the lists of genes differentially expressed by HHSt11 Cluster1 and Cluster2 versus all other in vivo NC cells to[…]”.

• "There may be even more dramatic differences between the pathways associated with each cluster in the form of significance." It is not clear what was the message the authors were trying to convey in this sentence.

• Subsection “Identification of two distinct transcriptional signatures associated with neural crest cells isolated shortly after dorsal neural tube exit”, last paragraph: It is intriguing to find that WNT/PCP signaling were respectively up and downregulated in clusters 1 and 2. It might shed important functional roles of these distinct subpopulations of stage 11 neural crest cells by altering WNT at this stage and to see if and how skewed production of these subpopulations lead to phenotypes in neural crest cell migration and formation of the Trailblazer population at later stages.

• "This conclusion is most likely the result of increased heterogeneity of these subpopulations due to the increased size when compared to the other subpopulations examined." The opposite appears to be true to common sense: reduced sample size will likely lead to increased heterogeneity. The findings presented here suggest that trail and lead subpopulations were highly heterogeneous and no defined signatures can be identified. Increasing the sample size (number of cells) in future experiments on the other hand might help to identify small subpopulations.

• The result description for Figure 6—figure supplement 1 should be elaborated and should be placed within the Results section instead of in the DDiscussion section.

Reviewer #2:

In the manuscript "Single-cell transcriptome analysis of neural crest migration reveals signatures of invasion and molecular transitions", Morrison and colleagues employ transcriptional profiling of individual neural crest cells to investigate cell heterogeneity during migration. The authors begin by performing bulk RNA-seq experiments to identify transcripts enriched and depleted in the invasive front. Next, they survey the transcriptome of single neural crest cells to establish the molecular signatures of distinct subpopulations within the R4 stream (also showing that there are substantial molecular differences between cells collected directly form embryos or cultured in vitro). The study is centered on the characterization of the transcriptional signature of trailblazers (cells that lead the neural crest migratory stream) and the molecules that set this subpopulation apart.

I found the manuscript to be well structured, given the challenge of condensing and organizing transcriptomic data from hundreds of cells. This is partly due to the decision (a good one, in my view) of focusing specifically on the differences between leading vs. trailing cells. Indeed, the manuscript represents a substantial addition to the work the group has done on this topic – while also resulting in a useful resource for future studies. Nevertheless, I believe it would greatly benefit from additional in vivo analysis. While reading the paper, I found myself trying to visualize how the distinct populations characterized through transcriptomics were distributed within the R4 stream. Where would the Trailblazer cells be and how many of them would be present? Mapping the expression of the novel trailblazer genes back into the embryo would be an important validation of the analysis, and yield in spatial information that could clarify the molecular heterogeneity observed.

1) How do the authors reconcile the results from the bulk to the single cell experiments? I would expect to find the genes that are part of the Trailblazer molecular signature to be enriched in the invasive front, but that is only the case for a handful of genes. Is this a matter of sensitivity / sequencing depth?

2) The use of single-cell transcriptomics allowed the authors to refine their previous findings on the molecular signature of the trailblazer cells. Yet I found the discussion on the significance of this signature to be somewhat vague. Are Trailblazer cells characterized by the co-expression of these genes? Are there any Trailblazer-specific individual factors? How heterogeneous is the Trailblazer population? In short – if one needs to identify Trailblazer cells, what are the factors they should be focusing on?

3) While I appreciate the authors' use of microdissection to incorporate a degree of spatial resolution in the analysis, the work would greatly benefit from additional in vivo experiments. The authors should confirm co-localization of at least some of the transcripts that are part of the trailblazer molecular signature in embryonic tissue. Is it possible to map individual trailblazer cells utilizing these markers? How large is the population of trailblazers in the stages analyzed, and is this consistent with the results of the single cell transcriptome data?

4) The morpholino experiments on Figure 6 come across as an afterthought. Have these morpholinos been published (I do not see the morpholino sequences in the method sections). How are these genes expressed? Proper controls for efficiency and specificity should also be included.

5) Surprisingly, Trailblazing cells have weaker expression of genes that are thought to be central to the neural crest regulatory program such as Tfap2b, Sox10, and Msx1. I was wondering how the authors would interpret these results, given the data showing that these genes are thought to be necessary for neural crest migration, and differentiation. Would the transcript levels present in these cells be sufficient for these processes to occur?

[Editors’ note: what now follows is the decision letter after the authors submitted for further consideration.]

Thank you for submitting your article "Single-cell transcriptome analysis of neural crest migration reveals signatures of invasion and molecular transitions" for consideration by eLife. Your article has been reviewed by two peer reviewers, and the evaluation has been overseen by Marianne Bronner as the Senior and Reviewing Editor. The reviewers have opted to remain anonymous.

The reviewers have discussed the reviews with one another and the Reviewing Editor has drafted this decision to help you prepare a revised submission. Although the essential changes are summarized below, the full reviews are included for your reference.

Summary:

In the present manuscript, Morrison and colleagues conduct a single-cell transcriptome analysis of the R4 neural crest to investigate cell heterogeneity during collective cell migration. The study builds upon previous work from this lab that has shown that a subset of neural crest cells (trailblazers) guides migration and is molecularly distinct from the rest of the migratory stream. The authors (1) perform bulk RNA-seq experiments to identify transcripts enriched and depleted in the invasive migratory front; (2) survey the transcriptome of single neural crest cells to establish the molecular signatures of distinct subpopulations within the stream; (3) characterize a transcriptional signature of trailblazer cells; (4) conduct perturbation assays to investigate how trailblazer genes affect neural crest migration; (5) analyze cell heterogeneity within leading and following cells.

Essential revisions:

1) Comparison of the bulk seq results presented in Figure 1D with the other datasets.

2) More robust message from the seq analysis from which one could envision a "logic" for the synexpression of certain genes/pathways in trailblazers. Figure 6 does not show restricted expression in Trailblazers – are these the best examples?

3) Provide RNAscope quantification to validate their Trailblazer identity.

4) Perform more RNAscope analysis on their possible highest confidence targets.

5) Provide 3D reconstruction of streams from RNAscope?

6) The mixed expression in leading edge and in more caudal regions of the stream needs to be pointed in Results and explained in the Discussion.

7) Morpholino experiments need to include controls that show loss of the target protein.

8) To address the modest effects from the morpholinos, compare MO targeted vs. MO untargeted HNK1 positive cells in streams. Perhaps some untargeted cells can continue to perform as Trailblazers, and the MO-targeted cells then become followers… explaining the somewhat limited effects. This could be easily addressed by following the label on the MOs and showing the difference against HNK1. One would expect that MO. negative HNK1+ cells would be in front of MO+ cells. Also this might enable a much better quantification of the effects.

9) It is crucial to match RNAscope and MO-analysis. MO analysis on all their RNAscope samples and/or RNAscope on all their MO samples. If they are doing new RNAscope and or MOs, this "matching RNAscope and MO analysis" argument must be addressed.

10) Quantification of RNAscope data to confirm positive cells are indeed in the leading edge of the stream.

11) Please add a Discussion paragraph addressing differences between their previous qPCR study and this work.

Reviewer #1:

This study has technical advances of potential high relevance for the NC and developmental biology communities. The research presented focuses on single cell RNA sequencing to identify profiles of migrating NC cells. Specifically directed to pinpoint differences between single cells along the migratory streams that have been proposed to occupy significant positions as front runners or trailblazers, immediate followers and more distant trackers. This type of work could identify the specific contribution of certain genes to the proposed actions needed to perform as a Trailblazer, the cells that open the road for followers. Perhaps it could identify distinct required functions enabling cells to simply follow more advanced cells. While this could certainly be a significant message for developmental biologists at large, the work presented does not reach these targets.

Here the authors have done few changes and while more functional analysis seem essential, they instead provide minimal evidence of either restricted expression, or requirement for a specific function during migration.

In this new version, the authors added one set of image analysis for mRNA levels within a immunodetected HNK1 positive migratory stream. While this particular technique could be of great use, it really does not provide clear confirmation of the status assigned to the genes tested since their expression is actually apparent at multiple distinct positions in the stream and not restricted or graded according to the expected Trailblazer status designated through the mRNA seq analysis. The true value should have been the identification of key genes expressed in a pattern associated with specific position according to their migration "status", as Front, Lead and Trail, or Trailblazers or leading edge, vs. less advanced cells. This validation is missing. While the authors describe a quantification of mRNA spots on HN1 labeled r4-strem cells. The data provided is not that clear. In fact, the images provided suggest that the genes identified and selected through the single cell RNAseq do not provide a restricted expression to the expected positions.

As critical as the expression profile or perhaps even more, it would have been the identification of the functions provided from such genes. To this end the authors transferred some of the supplementary data from their original version, with identical images, to the main Results section of the revised manuscript. And tise set of KD results shows modest effects on the migration of neural crest cells at best. The evidence does not suggest a robust correlation towards a specific migratory status or position. While there are multiple genes in the table provided with apparent higher correlation to the Trailblazer, they only tried Wnt5B from those.

Similarly, while they identify two interesting groups or sets from the stage 11 embryos, which seem pretty interesting, again, no validation via expression is provided, nor a functional assessment.

Finally, while they have made an effort to compare their new data, to their 2015 data from QPCR, a clear message broadly assessing the differences, and commonalities, is still missing. Even more, no clear outcome of the analysis is evident. It is difficult to pinpoint an eloquent conclusion or message.

Reviewer #2:

As expected, the study is rich in transcriptomic data – the authors sequenced hundreds of cells of chick embryos from different stages. The results are properly organized and presented, and despite the abundance of heat maps and bar graphs the story flows well. The study is not without limitations – single cell transcriptomes are inherently noisy and fully distinguishing technical noise from true biological heterogeneity is challenging, if not impossible. But taken together, I believe the results represent a significant contribution. Indeed, this study perhaps makes the most compelling case for the existence/importance of Trailblazing cells in the migratory neural crest.

Of particular note is how the scRNA-seq results relate to the in vivo development of the neural crest. By comparing native cells to neural crest cultured in vitro, the authors show that significant transcriptional changes emerge with cell culture. Furthermore, the authors also attempt to bring some of the scRNA-seq findings back to the embryo, by conducting a number of in vivo experiments. They map the expression of some of the trailblazing genes by RNAscope/FRUIT, and subsequently they also show that a subset of those is required for proper neural crest migration. The first experiment produced specially compelling data – one is able to visualize individual Trailblazer cells within the migratory stream. The latter is potentially interesting but in my view needs additional controls.

[Editors' note: further revisions were requested prior to acceptance, as described below.]

Thank you for submitting your article "Single-cell transcriptome analysis of neural crest migration reveals signatures of invasion and molecular transitions" for consideration by eLife. Your article has been reviewed by two peer reviewers, and the evaluation has been overseen by Marianne Bronner as the Senior and Reviewing Editor. The reviewers have opted to remain anonymous. The consensus opinion is that your paper could be acceptable as a Tools and Resources (TR) paper with some additional validation as discussed below. Our sense is that in order to qualify as a Research Article (RA), you would need to perform significant additional functional tests that may take quite a bit of time to complete.

Summary:

The consensus opinion is that the manuscript has been improved by the additional of functional data but it remains critical to show that a "Trailblazer" signature is actually unique to the leading cells in order to support the authors claims.

Essential revisions:

The full reviews of the reviewers are included below for the authors' reference. Both reviewers feel that further validation of the Trailblazer signature is necessary. Thus, it is critical that you do one of the following, with the first suggestion best in line for a TR paper:

1) Identify 5-10 genes that are either restricted at the front, or have graded expression with highest intensity towards the leading edge; and/or

2) Perform additional functional experiments to see if any of the identified genes is functionally important with a major phenotype.

Please also note that you should include better controls for the morpholino experiments.

I hope you find the reviewers comments helpful in revising the paper.

Reviewer #1:

This revised version has added changes in a short time, and the authors have put effort in amending the previous version. However, the main objection to this study has been the biological evidence for the significance of the transcriptomic signatures suggested for the trailblazers, the foremost advanced migratory crest cells leading the way for follower neural crest cells. Unfortunately, the evidence presented seems insufficient.

Key to the significance is Figure 6, looking at expression and function. Expression profiles were monitored through quantifiable in situ hybridization (RNAscope) and co-immunostaining with HNK1 to reveal the full migratory stream. The results are intriguing, because while there are clearly some cells near the front with high expression, there is no apparent restriction to the leading cells. Instead each gene presented shows a) expression in some (not all, not even the majority) cells near the edge being very strong, and b) multiple patterns, with some displaying strong expression at the opposite end of the stream (back delayed followers), or leading edge or just more anterior cells within the stream with intermediate levels of expression. These results offer a complex pattern not easily ascribed to a restricted subset of cells in a particular level. While the quantification provided may suggest that there is a pattern, their analysis is limited to a restricted number of cells, nearby their ideal high-expressing target. We previously asked, are these the best candidates? This being the third version of the manuscript, may suggest so.

Then they assess the function of some target genes in the migration. They do not address exactly the same monitored for expression (Figure 6A) which would have been a nice demonstration had they observed a clear phenotype even if the expression was not optimal. They suggest to have used the same criteria carefully laid out for the expression analysis, but do not restrict themselves to the top Trailblazer signature and instead they added targets based on migration-likelihood as criteria. If these genes are expressed in a pattern to support the Trailblazer claim, that would be perfect, but no additional expression evidence is presented. But the most critical aspect is simply that the Morpholino treatment of any of these genes, do not provide a strong phenotype. The average migration effect seen by HNK1 in treated streams only displays a reduction of 10%! Control Mo displays a range around 96% migration, their best candidate Wnt5b, has an average of 91% and the next two best candidates GATA5 and LUM show about 94% migration. This does not seem to be a robust phenotype. They suggest that this is due to migration from wt, cells that did not received the morpholinos, but no evidence supports this claim. The results presented in the main Figure 6G-I do not provide a robust effect. Few of their electroporations leads to a reduced expression compared to the control range, their top candidate Wnt5b being perhaps 15% reduced over the control. It would have been good to show the images of that specific KD, so even if it is limited at least people could see it. The images for panel G provide a confusing picture. IN one hand, EDN1 and PKP2 MO treated cells appear near the edge of the stream…

Reviewer #2:

In the revised version of manuscript, Morrison and colleagues conduct an extensive single-cell transcriptome analysis to investigate cell heterogeneity and Trailblazer cell identity in neural crest cells from the R4 migratory stream. The results reveal a transcriptional signature that remains essentially unchanged in cells isolated from embryos in distinct stages of development (HH11, HH13 and HH15) (Figure 4). The authors postulate that this signature would define the identity of the Trailblazer cells, which are located at the migratory front and have a crucial role in the establishment of the migratory stream. The unbiased identification of a subpopulation of neural crest cells through scRNA-seq is an important result; however, as highlighted in previews comments, it is crucial that such findings be validated in vivo.

This reviewer finds that the new version of the manuscript does not unequivocally show that Trailblazer cells express the genes that are part of their 'transcriptional signature'. In Figure 6, the authors employ RNAscope to visualize and quantify the expression of the genes mostly enriched in the Trailblazer signature. The results are difficult to interpret, as the cells that have the highest expression of the genes are in many cases far away from the leading edge (Desmin). In another example, the target gene seems to be downregulated in the migratory front (Tescalcin). I found the figure confusing, as the authors opt to quantify transcript levels in the cells that have the highest expression of the transcript, rather than focusing on the cells at the leading edge (the request from the reviewers was RNAscope quantification in Trailblazer cells). Thus, while I agree that this data shows neural crest heterogeneity, I don't see how it supports the hypothesis of a Trailblazer transcriptional signature.

I acknowledge the authors' comment on how the scRNA-seq delineates a complex trailblazing signature, composed of hundreds of genes, but I don't see how this this should hinder in vivo validation of the transcriptome data. Figure 4 clearly shows a group of genes that are enriched in the trailblazing cells and depleted in other populations. If the analysis and the hypothesis are correct, one would expect that expression mapping would show higher levels and co-expression of these genes in Trailblazers. It is this reviewer's opinion that in vivo validation of RNA-seq – especially in the case of noisy single cell data – should be emphasized. Additionally, the lack of consistency in the genes that are highlighted in the bioinformatic analysis, the genes that are analyzed with RNAscope and the genes that are chosen for loss of function analysis raises questions about the dataset. I have no problem with the authors focusing in a particular subset of genes for the in vivo experiments, but a certain level of consistency would make the manuscript easier to follow.

As previously stated, the scRNA-seq analysis performed by the authors represents an important resource for community, yielding an interesting hypothesis. However, the fact that the link between the scRNA-seq results and the in vivo analysis is tenuous creates a disconnect between the two parts of the manuscript. This should be resolved for improvement on the consistency and reliability of the study.

[Editors' note: further revisions were requested prior to acceptance, as described below.]

Thank you for submitting your article "Single-cell transcriptome analysis of neural crest migration reveals signatures of invasion and molecular transitions" for consideration by eLife. Your article has been reviewed by two peer reviewers, and the evaluation has been overseen by a Marianne Bronner as the Senior and Reviewing Editor. The reviewers have opted to remain anonymous.

The reviewers have discussed the reviews with one another and the Reviewing Editor has drafted this decision to help you prepare a revised submission. I have attached the full reviews of both reviewers which are highly concordant. While they feel that the data are valuable as a Tools and Resource paper, both reviewers felt that the way the manuscript was presented was confusing and that it would be important to rewrite the manuscript in a manner that doesn't over interpret the results. I think that you would be able to do this with textual revisions alone but would urge you to spend particular care on addressing the points raised in the full reviews, as both reviewers commented that they felt you had not addressed their earlier concerns and were trying too hard to put your conclusions in a positive light rather than addressing the caveats of the single cell RNA-seq analysis and the difficulties of linking this to in vivo expression analysis and function.

Given the history of this manuscript, this must be your last opportunity to produce an acceptable version.

Essential revisions:

1) The authors should aim to provide a clear, and accurate description of their findings. These in vivo findings need to be spelled clearly: These genes were identified by scRNAseq and bioinformatics (through a stringent protocol with harsh significance in mind), as specific to the expected Trailblazer at the narrow front edge of the migratory front, however, they do not seem to appear restricted in that fashion in vivo, and insteada) None of the 13 genes tested so far (which include some of the top candidates) are expressed in all of, or the majority of the cells at the front edge.b) some may appear in a few or some but not in all of the cells (~30% or 1/3?) of the (narrow) front migratory edge.c) some of the genes may display salt and pepper patterns which may include strong expression in a few cells in the narrow front edge and similarly strong signal in some cells in the beginning or middle portions of the stream (not at the front edge of migration).d) Some of the genes where not found at the front edge at all, and instead appear to be expressed more strongly at various positions within the migratory stream.

2) The authors should provide a revised version of the so called "complex profile/signature of expression for a Trailblazer" If the genes are not specifically expressed in all of them, or the majority, and if these genes may be expressed in the stream, then what is the so called complex signature of the trail blazer? (Not to mention its relevance…) The authors could offer some useful insights regarding the limitations or caveats of the scRNAseq at the current time/or in these complex settings. Perhaps cumulative RNA vs. freshly made RNA.

3) The novel RNAscope data shows heterogeneous expression of genes that are part of the "Trailblazer signature". While this reviewer finds that a thorough transcriptomic characterization of single RNA cells may be a useful resource for the community, the results of the in vivo analysis require further clarification.

4) Upon presenting the expanded RNAscope analysis, the authors conclude that "the combination of fluorescence in situ hybridization with IHC and tissue clearing allowed the confirmation of a small subset of the Trailblazer transcriptional signature". This statement seems to imply that heterogeneous expression is sufficient to validate transcripts as being part of the Trailblazer signature. According the Introduction of the manuscript such genes should be expressed by cells "narrowly confined by the migratory front". This discrepancy needs to be addressed or the reader will not be able to link the transcriptomic data with the in vivo expression analysis.

Reviewer #1:

The Authors previously used cell sorting and q-RT-PCR to identify a unique molecular signature of 16 out of 96 genes that was stable and consistent in cranial neural crest cells narrowly confined to the most invasive front of the migratory stream that we termed 'Trailblazers' (McLennan et al., 2015). In this manuscript, they have further scrutinized the expression profile of individual cells to address molecular changes and singularities that could help us better understand the behavior of migrating cells, specifically those at the fore front of the migration. The appealing idea being that these front runners might be endowed with unique capabilities to determine direction, remove ECM, perhaps even push other cells to the side, open a path, and lay the road and or instructions for follower cells.

There is a lot of value in this hard and fine work. They had accomplished a technically challenging task addressing the expression profiles of migrating neural crest cells. They provide significant information regarding the differences between in vitro cultured and in vivo migrating neural crest cells. They confirm some of the genes they previously identified as enriched in the "Trailblazer" cells using a single cell RNA seq approach. Yet the outcome of this huge endeavor in its present format, is misleading. As presented by the authors, one would assume that they have identified a unique "expression signature" of Trailblazer cells, which migrate at the front edge of the migratory stream in chick embryos, which they define as "narrowly confined by the migratory front". The implication being that a few genes would be expressed in a specific pattern associated with the migratory front edge. Does this mean a strong expression of a few genes consistent at the front edge? This argument is made throughout the paper in multiple places. But they do not provide the evidence necessary to grant such claim, and should be corrected throughout.

Key to these arguments would be a clear, measurable, objective definition of the front edge. In terms of its distance from the dorsal neural tube or departure point, truly the edge, near the edge, 2, 3 5 cells behind the apparent edge? And then the width or range of the front edge, all the cells visible at the edge? Some angle or curvature to restrict it? And what is the volume of the front edge? The use of ill-defined, subjective terms provides ample room for misinterpretation and confusion.

If we were to think of our arm and hand as a migratory stream departing from our trunk, what could be considered as the Trailblazers, the nails? The tip of the thumb? Or the index fingers only, the distal, middle or proximal portions of some or all of the phalanges? Would it be enough to look at the back of the hand/dorsal portion of the phalanges, or the palm portion/ ventral, or the middle? etc…

While performing the transcriptomic portion of the work, a lot of effort was placed to define and collect cells from specific regions of the stream to be analyzed. Yet upon validation of expression in embryos, a much laxer approach is taken and confusion ensues.

The heading of one of the Results sections reads: "Trailblazer cells are identified by an ensemble of genes", and after presenting the results in cryptic and limited format, concludes: "Thus, the combination of fluorescence in situ hybridization with IHC and tissue clearing allowed the confirmation of a small subset of the Trailblazer transcriptional signature." Yet the expression of the genes here presented is not restricted or biased to the front edge.

The authors perform single cell transcriptomics and come up with a selected list of genes thought to decipher the molecular signature of the "Trailblazer" cell. They try to corroborate their thorough and delicate transcriptomics and bioinformatics work via two assessments. The first one consisting of expression analysis, and the second on Morpholino mediated knockdown. The knockdown provided minimal effects at best and thus could not confirm the tenants of a "relevant" trailblazer expression. Perhaps combinatorial approaches could reveal function, but at this point no more evidence is provided.

To demonstrate expression profile associated with Trailblazer character, the authors performed a nice combination of immunofluorescence to identify all the migratory NC via HNK1, combined with quantifiable in situ hybridization analysis, "RNAscope". The criteria provided to select genes to prove and validate the trailblazer character states:

1) genes with enriched expression in Trailblazer cells;

2) high percentage of single Trailblazer cells in which the gene was detected and;

3) high level at which the gene was expressed by the Trailblazer cells.

(The genes tested are PKP2, DESMIN, KAZALD1, TESCALCIN, TROPONIN, ANEXIN, BAMBI, GLYPICAN3, DESMOPLAKIN, BDNF, TFPI2, NEXILIN, CDH5).

In no case, a unique pattern is presented that strictly labels with the highest signal the expected trailblazer cells at the front edge of the migratory stream defined as "narrowly confined by the migratory front". Instead, in the best cases for their argument (PKP2, and KAZALD1), expression seem to appear strong or strongest in a few of the visible trailblazer cells. These genes do not appear consistently in their strongest expression profile at the front edge, in all the trailblazers. Some Trailblazer cells clearly do not express these genes and or express them weaker than follower cells. And in multiple cases stream or follower cells display a wide range of expression levels varying from the strongest to no expression. And one has to wonder what is the depth of the image presented? Does this represent the whole edge, or would one imagine 3 similar projections to cover the volume of the front edge? If so, what is the number of cells at the front edge that strongly express the interrogated marker? Further, when they analyze the expression relative to neighbors, the effort is limited to very few cells, when the same arguments made for distance to the edge, width and depth should apply.

This is even more clear when analyzing the other 11 genes they present. Desmin, and Troponin1, display 1 and 2 cells respectively, with the strongest expression (red) far behind the expected Trailblazers (are not at the very front) and instead are behind multiple cells with 2/3 or 1/2 as much expression. And unfortunately, the evidence provided for their analysis of the additional other 8 genes, does not support their case either, and actually are not a better example of what one would have expected for a trailblazer profile of expression.

Importantly, for some undisclosed reason, they do not include the same targets on expression and function (only 3 of the 9 targets assessed via Morpholinos are also tested amongst the 13 targets used for expression analysis with RNSAscope).

Given the lack of evidence of relevant role for the genes tested, and that no signature is confirmed by expression in the embryo, the current version would be a misleading tool or resource.

The authors should aim to provide a clear, and accurate description of their findings. These in vivo findings need to be spelled clearly:

These genes were identified by scRNAseq and bioinformatics (through a stringent protocol with harsh significance in mind), as specific to the expected Trailblazer at the narrow front edge of the migratory front, however, they do not seem to appear restricted in that fashion in vivo, and insteada) None of the 13 genes tested so far (which include some of the top candidates) are expressed in all of, or the majority of the cells at the front edge.b) some may appear in a few or some but not in all of the cells (~30% or 1/3?) of the (narrow) front migratory edge.c) some of the genes may display salt and pepper patterns which may include strong expression in a few cells in the narrow front edge and similarly strong signal in some cells in the beginning or middle portions of the stream (not at the front edge of migration).d) some of the genes where not found at the front edge at all, and instead appear to be expressed more strongly at various positions within the migratory stream.

The authors should provide a revised version of the so called "complex profile/signature of expression for a Trailblazer" If the genes are not specifically expressed in all of them, or the majority, and if these genes may be expressed in the stream, then what is the so called complex signature of the trail blazer? (Not to mention its relevance…) The authors could offer some useful insights regarding the limitations or caveats of the scRNAseq at the current time/or in these complex settings. Perhaps cumulative RNA vs freshly made RNA…

Reviewer #2:

This revised version of the manuscript by Morrisson and colleagues on scRNA-seq analysis of neural crest cells includes an expanded RNAscope analysis of putative Trailblazer genes. Reviews of the previous versions of the manuscript emphasize the need for further in vivo validation of the transcriptomic datasets. The novel RNAscope data shows heterogeneous expression of genes that are part of the "Trailbrazer signature". While this reviewer finds that a thorough transcriptomic characterization of single RNA cells may be a useful resource for the community, the results of the in vivo analysis require further clarification.

Upon presenting the expanded RNAscope analysis, the authors conclude that "the combination of fluorescence in situ hybridization with IHC and tissue clearing allowed the confirmation of a small subset of the Trailblazer transcriptional signature". This statement seems to imply that heterogeneous expression is sufficient to validate transcripts as being part of the Trailblazer signature. According to the Introduction of the manuscript such genes should be expressed by cells "narrowly confined by the migratory front". This discrepancy needs to be addressed or the reader will not be able to link the transcriptomic data with the in vivo expression analysis.

https://doi.org/10.7554/eLife.28415.028

Author response

[Editors’ note: the author responses to the first round of peer review follow.]

Reviewer #1:

[…] Results section 1

• Why were these two stages being chosen? Stage 13 appears to be very close to stage 15 in terms of gene expression. Bioinformatics analyses also did not appear to justify there was a stage specific difference.

We agree and have clarified this in the text. Stage 11, 13 and 15 are separated developmentally by approximately 8 hours. Stage 11 coincides with the initiation of migration from the neural tube, Stage 13 is active migration and Stage 15 is colonization of the target tissue. Only after doing these experiments did we discover the similarity between Stage 13 and 15.

• At what stage did the authors label the premigratory neural crest cells?

We have clarified this in the text.

• There was no description whether the 16 genes identified in McLennan et al., 2015 were found in the bulk RNA seq and what were the new players. For example the analyses should have been performed on genes such as, Hand2/GPC3 as invasive front signature, PDGFRL/CXCR1 as stream signature, and EPHA4 and FOXD3 as trailing signature identified in McLennan et al. 2015. These need to be first established in the bulk RNA seq analysis.

We have clarified this in the text and in a figure. The complete differential expression analysis of the bulk subpopulations (both validating our previous findings and establishing new signatures) may be found in Supplementary file 1.

Results section 2

• "we isolated and profiled individual cells from different stream positions at three developmental stages". “We subdivided neural crest cell migratory streams at HHSt13 and HHSt15 into three spatially distinct, non-overlapping subregions by careful manual dissection and single cells were isolated by FACS". The authors need to justify the rationales behind the selection of stages (two versus three) and dissection methods (front/trail vs. front/lead/trail) for the bulk and single cell RNA seq. Were there other previous experiments and studies that the authors or others had performed that would help justify the differences?

We have expanded our rationale in the text.

• "Hierarchical clustering of single neural crest cell averages….." The initial analysis should have been performed without averaging the single cell transcriptome values according to their migration positions. The current experimental design, by clustering single cell samples, could have unavoidably biased the results. The conclusion reached also did not add any new insight into the study. Figure 2C could be removed. The aim of the experiment was to identify new genes/signatures that can define the different subpopulations and the authors could have used the previous 16 genes (McLennan et al., 2015) or other known targets from the bulk RNA seq in Figure 1 or other studies as a guide to distinguish the front cells.

We have moved Figure 2C (bulk RNA-seq averages) to the supplementary material. Although either the RT-qPCR or bulk RNA-seq signatures could have been used to identify new genes/signatures to define cell subpopulations, these data were constrained by the manual cuts to subdivide the stream. Instead, we chose an unbiased approach to identify novel cell subpopulations independent of the position of the manual dissections to subdivide the stream and have clarified this in the text.

Results section 3

• "in vitro transcriptional signatures are distinct from any in vivo profile". This section appears not relevant to the rest of the study. This section does not give any new insight into the cellular heterogeneity of migrating neural crest cells.

We appreciate the reviewer’s suggestion. Since discrete neural crest cell migratory streams are not re-capitulated in vitro, but in vitro assays are readily used for cell migration studies, we wished to confirm our hypothesis and previous RT-qPCR data (McLennan et al., 2015) that in vitro transcriptional signatures are very distinct from our in vivo results and is useful information. Thus, further supporting the model that in vivo neural crest microenvironmental signals play an important role in establishing in vivo transcriptional signatures. We have moved the in vitro data to later in the manuscript and clarified this in the text.

Rest of the Results sections

• “Our previous RT-qPCR (McLennan et al., 2015) and bulk RNA-seq presented above (Figure 1) show that cranial neural crest cells at the invasive front of a stream have aspects of their transcriptome that are conserved during migration.” This rationale does not fit into the experimental design for Figure 5. If the invasive front were to have a conserved transcriptome during migration, then all cells should have been clustered for analysis. The experiment presented in Figure 5 instead appears to try to identify heterogeneity of the front cells dissected from the embryos.

We have clarified this in the text.

Figure 5: Where did the front cells for stage 11 come from? It was stated in the first paragraph of the subsection “Single-cell RNA-seq identifies gene expression variances based upon spatial position within the neural crest cell stream and temporal progression along the migratory pathway” that there was no subdivision of stage 11 cells.

We refer to HHSt11 cells as FRONT, since the small number of migrating neural crest cells adjacent to r4 does not allow for subdivision of cells into subgroups. We have clarified this in the text.

• Use different colors for Trail and Trailblazers in Figure 6C.

Done.

• Subsection “The scRNA-seq Trailblazer signature is uniquely defined and associated with neural crest cells confined to the invasive front of the stream”, last paragraph: It is interesting that the authors were able to identify previously known pathways, such as the VEGF signaling, that were associated with Trailblazer identity. However, as the purpose was to discover new pathways, an unbiased search should be performed instead of a focus on developmental relevant signaling pathways. As noted, a significantly higher number of pathways were found to be upregulated and many of the p-values were relatively high (Figure 6D) suggesting that there might exist other "developmentally irrelevant" pathways that were associated with trailblazer signatures.

This is an excellent suggestion and we have added lists of all statistically significant pathways for all single cell RNA-seq subpopulations discussed in Supplementary file 6.

• Please highlight the p-value in Figure 6D that represent significance difference.

Done.

Figure 6D can be better presented with the p-values align according to their values from top to bottom.

We appreciate the reviewer’s suggestion but wished to highlight the links to the VEGF and axonal guidance signaling pathways and thus plotted by alphabetical order. P values can still be easily seen in the graph.

• What were the criteria to select the gene list for presentation in Figure 6B? Were these genes most highly enriched or depleted?

We selected genes that we considered of general interest to the neural crest and embryonic cell migration fields. The most enriched and depleted genes may be found in the differential gene lists contained in Supplementary file 4.

• Highlight in the text and figure legend for Figure 6B that those genes highlighted in red were previously found in McLennan et al., 2015.

Done.

• The conclusion "These results highlight the upregulation of many genes within the trailblazer subpopulation, supporting the hypothesis that trailblazers are a dramatically active subpopulation of cells at the most invasive front during neural crest migration" needs to be substantiated (and if possible reworded to omit superlatives.

We agree and have reworded the text.

• " We also found that CDH7, CXCR4 and EPHB1 were upregulated in the HHSt11 single cells (Figure 6E)." What were these single cells? Were there any bias in stage 11 clusters 1 and 2 identified earlier?

We have included text within the manuscript to describe these results more accurately and have added the differential expression analysis between HHSt11 Cluster 1 and 2 to Supplementary file 4.

• Subsection “The scRNA-seq analysis validates some of our original 16 gene RT-qPCR Trailblazer signature”: the authors said there was "no significant enrichment" for a number of genes found in previous studies meaning the expression of these genes were upregulated but not statistically significant or were downregulated instead? Following up on this, if these genes were not upregulated in the trailblazer population, were they found to be upregulated in any of the subpopulations of the migrating neural crest cells?

We agree and have clarified this in the text.

• Subsection “Identification of two distinct transcriptional signatures associated with neural crest cells isolated shortly after dorsal neural tube exit”, second paragraph: How was the GO search performed for clusters 1 and 2? Why were there no depleted and enriched terms for clusters 1 and 2 respectively?

The Gene Ontology (GO) analysis was completed using the R package ClusterProfiler (3.2.6). GO terms with a p value less than 0.05 and a q value below 0.1 were considered significant. Significant GO terms were condensed using the simplify function with an adjusted p value cutoff of 0.05. While the results are surprising, there simply were no depleted or enriched terms because no pathways contained enough genes from our differential expression analysis to be statistically significant. We have clarified this in the text.

• Please correct typos. “To identify signaling pathways enriched within HHSt11 clusters, we performed pathway enrichment analyses "for" the lists of genes differentially expressed by HHSt11 Cluster1 and Cluster2 versus all other in vivo NC cells to[…]”.

Done.

• "There may be even more dramatic differences between the pathways associated with each cluster in the form of significance." It is not clear what was the message the authors were trying to convey in this sentence.

We have clarified this in the text.

• Subsection “Identification of two distinct transcriptional signatures associated with neural crest cells isolated shortly after dorsal neural tube exit”, last paragraph: It is intriguing to find that WNT/PCP signaling were respectively up and downregulated in clusters 1 and 2. It might shed important functional roles of these distinct subpopulations of stage 11 neural crest cells by altering WNT at this stage and to see if and how skewed production of these subpopulations lead to phenotypes in neural crest cell migration and formation of the Trailblazer population at later stages.

We agree that the role of WNT/PCP signaling deserves future investigation. We have added a list of genes differentially expressed between HHSt11 Cluster 1 and 2 to Supplementary file 4 and highlighted genes associated with WNT/Β-Catenin and PCP signaling pathways.

• "This conclusion is most likely the result of increased heterogeneity of these subpopulations due to the increased size when compared to the other subpopulations examined." The opposite appears to be true to common sense: reduced sample size will likely lead to increased heterogeneity. The findings presented here suggest that trail and lead subpopulations were highly heterogeneous and no defined signatures can be identified. Increasing the sample size (number of cells) in future experiments on the other hand might help to identify small subpopulations.

We have clarified this in the text. The measurement of a population’s heterogeneity is based upon both the sample size (as the reviewer notes) as well as the actual heterogeneity of the population. While similar numbers of cells were analyzed from each subpopulation (50-60 cells/each; Figure 2B), the total numbers of Leader and Trailer cells within the stream is much larger than the number of Front cells.Additionally, McLennan 2015 showed regional diversity throughout the stream, thus we would expect a larger subpopulation of the stream (Trailers, 70% of the stream) would have more heterogeneity than a smaller subpopulation (Invasive Front, 5% of the stream). Less statistically significant differential expression is likely to occur when larger subpopulations with more heterogeneity are averaged.

• The result description for Figure 6—figure supplement 1 should be elaborated and should be placed within the Results section instead of in the Discussion section.

Done.

Reviewer #2: […] I found the manuscript to be well structured, given the challenge of condensing and organizing transcriptomic data from hundreds of cells. This is partly due to the decision (a good one, in my view) of focusing specifically on the differences between leading vs. trailing cells. Indeed, the manuscript represents a substantial addition to the work the group has done on this topic – while also resulting in a useful resource for future studies. Nevertheless, I believe it would greatly benefit from additional in vivo analysis. While reading the paper, I found myself trying to visualize how the distinct populations characterized through transcriptomics were distributed within the R4 stream. Where would the Trailblazer cells be and how many of them would be present? Mapping the expression of the novel trailblazer genes back into the embryo would be an important validation of the analysis, and yield in spatial information that could clarify the molecular heterogeneity observed.

1) How do the authors reconcile the results from the bulk to the single cell experiments? I would expect to find the genes that are part of the Trailblazer molecular signature to be enriched in the invasive front, but that is only the case for a handful of genes. Is this a matter of sensitivity / sequencing depth?

We agree and appreciate the suggestion to clarify this in the text. Briefly, single cell RNA-seq identified Trailblazers as a subset of the Front populations (8, 50 and 32% at HHSt11, 13 and 15, respectively). Therefore, when entire Front populations are averaged, as in the bulk analyses, some of the Trailblazer genes may or may not be statistically significant.

2) The use of single-cell transcriptomics allowed the authors to refine their previous findings on the molecular signature of the trailblazer cells. Yet I found the discussion on the significance of this signature to be somewhat vague. Are Trailblazer cells characterized by the co-expression of these genes? Are there any Trailblazer-specific individual factors? How heterogeneous is the Trailblazer population? In short – if one needs to identify Trailblazer cells, what are the factors they should be focusing on?

We have clarified this in the text.Briefly, our scRNAseq analysis did not reveal any individual gene that is able to distinguish the Trailblazers from the remaining migratory stream. Rather, we assert that Trailblazers are characterized by a large number of differentially expressed genes (added in Supplementary file 4).

3) While I appreciate the authors' use of microdissection to incorporate a degree of spatial resolution in the analysis, the work would greatly benefit from additional in vivo experiments. The authors should confirm co-localization of at least some of the transcripts that are part of the trailblazer molecular signature in embryonic tissue. Is it possible to map individual trailblazer cells utilizing these markers? How large is the population of trailblazers in the stages analyzed, and is this consistent with the results of the single cell transcriptome data?

We agree and have selected a subset of genes to confirm expression based on the criteria: 1) genes with enriched expression in Trailblazers, 2) the percentage of single Trailblazer cells in which the gene is detected, and 3) the level at which the gene is expressed by the Trailblazers.This is presented in a revised Figure 7 and expanded in the Materials and methods section.

4) The morpholino experiments on Figure 6 come across as an afterthought. Have these morpholinos been published (I do not see the morpholino sequences in the method sections). How are these genes expressed? Proper controls for efficiency and specificity should also be included.

We agree and have included the additional data. We have expanded our investigation of Trailblazer function via morpholinos and have included the data in Figure 7 (Morpholino sequences are included in Figure 6—source data 1).

5) Surprisingly, Trailblazing cells have weaker expression of genes that are thought to be central to the neural crest regulatory program such as Tfap2b, Sox10, and Msx1. I was wondering how the authors would interpret these results, given the data showing that these genes are thought to be necessary for neural crest migration, and differentiation. Would the transcript levels present in these cells be sufficient for these processes to occur?

Sox10 for example, is highest in the most recently delaminated neural crest cells and its expression decreases over developmental time. That is, by HH St15, the front of the stream is no longer Sox10 positive. We have seen this in qPCR, RNAseq and in situ approaches over the years. We see a similar expression profile for FoxD3. We have not done in depth expression analysis of other known neural crest genes such as TFAP2B or MSX1, but they are mainly known to be involved in the specification of neural crest, not necessarily the active invasion of unexplored microenvironments. In this study, both TFAP2B and MSX1 are expressed higher at Stage 11 that at the later stages, consistent with their roles in the neural crest cell gene regulatory network (Supplementary file 2). Furthermore, this study did not compare neural crest cells to other cell types and therefore genes that are reasonably homogeneously expressed by neural crest cells will not be highlighted.

[Editors' note: the author responses to the re-review follow.]

Essential revisions:

1) Comparison of the bulk seq results presented in Figure 1D with the other datasets.

We have included the requested comparisons. Supplementary file 1 compares differentially expressed genes between the Front and Stream by bulk RNAseq and highlights the presence of genes from the 16 gene scRT-qPCR Trailblazer signature. Supplementary file 4 highlights the presence of the enriched and reduced genes from Figure 1D-E in the lists of genes differentially expressed by Trailblazers from scRNA-seq. We have included text in the Results section to describe the comparisons.

2) More robust message from the seq analysis from which one could envision a "logic" for the synexpression of certain genes/pathways in trailblazers. Figure 6 does not show restricted expression in Trailblazers – are these the best examples?

Done.

3) Provide RNAscope quantification to validate their Trailblazer identity.

Done (revised Figure 6).

4) Perform more RNAscope analysis on their possible highest confidence targets.

We agree with the reviewer that the validation of our Trailblazer transcriptional signature would enhance the value of our data by identifying genes and gene families for further functional analysis. To begin to address this, we confirmed that the Trailblazer transcriptional signature validated many of the genes we identified by RT-qPCR analysis (Supplementary files 1 and 4) to be expressed high and consistent over time (McLennan et al., 2015). Second, the Trailblazer transcriptional signature is defined as an ensemble of differentially expressed genes (900 enhanced, 400 reduced) from our PCA analysis and highlights the heterogeneity of expression within cells at the invasive front. To address the heterogeneity of expression, we confirmed that neighboring cells within the invasive front and stream may have dramatic differences in mRNA expression of an individual gene by using multiplexed fluorescence in-situ hybridization (RNAscope) and quantitative single cell analysis (revised Figure 6).

We completely understand the reviewer’s request to visualize the expression of an individual gene associated with the Trailblazer transcriptional signature, but suggest this may or may not uniquely identify a Trailblazer cell. For example, if Gene A is enhanced within the Trailblazer transcriptional signature ensemble, its expression may also be high within follower cells, but since the follower cell lacks the full complement of differentially expressed genes as part of the ensemble, it would not uniquely identify a Trailblazer cell. To clarify this, we have: (i) added quantitation and a 3D movie of our RNAscope data to strengthen the message of the heterogeneity of individual gene expression of Trailblazer transcriptional signature (revised Figure 6; Video 1); (ii) added the analysis of 3 other genes of the Trailblazer transcriptional signature and morpholino knockdown experiments (revised Figure 6 and Figure 6—figure supplement 2) and; (iii) edited the text to elucidate the ensemble aspect of the Trailblazer transcriptional signature.

5) Provide 3D reconstruction of streams from RNAscope?

Done.

6) The mixed expression in leading edge and in more caudal regions of the stream needs to be pointed in Results and explained in the Discussion.

Done.

7) Morpholino experiments need to include controls that show loss of the target protein.

We have included supplementary material that shows alterations in target mRNA in chick cells subjected to splice blocking morpholinos and added text to the Results and Methods sections (revised Figure 6figure supplement 2), as we do not have antibodies to these genes to visualize protein levels.

8) To address the modest effects from the morpholinos, compare MO targeted vs. MO untargeted HNK1 positive cells in streams. Perhaps some untargeted cells can continue to perform as Trailblazers, and the MO-targeted cells then become followers… explaining the somewhat limited effects. This could be easily addressed by following the label on the MOs and showing the difference against HNK1. One would expect that MO. negative HNK1+ cells would be in front of MO+ cells. Also this might enable a much better quantification of the effects.

Done (revised Figure 6—figure supplement 3).

9) It is crucial to match RNAscope and MO-analysis. MO analysis on all their RNAscope samples and/or RNAscope on all their MO samples. If they are doing new RNAscope and or MOs, this "matching RNAscope and MO analysis" argument must be addressed.

For expression analysis by RNAscope, we selected genes for expression confirmation based on the criteria: (1) genes with enriched expression in Trailblazer cells; (2) high percentage of single Trailblazer cells in which the gene was detected and; (3) high level at which the gene was expressed by the Trailblazer cells. For knock down experiments, we selected 10 genes from distinct molecular families to test, based on the same criteria used for selecting genes for expression analysis, but also for predicted roles in neural crest migration. Thus, different rationales produced to divergent lists of genes. We have now added the expression analysis of 3 genes from the knock down list (revised Figure 6).

10) Quantification of RNAscope data to confirm positive cells are indeed in the leading edge of the stream.

See #3 and #4.

11) Please add a Discussion paragraph addressing differences between their previous qPCR study and this work.

Done.

[Editors' note: further revisions were requested prior to acceptance, as described below.]

Essential revisions:

The full reviews of the reviewers are included below for the authors' reference. Both reviewers feel that further validation of the Trailblazer signature is necessary. Thus, it is critical that you do one of the following, with the first suggestion best in line for a TR paper:

1) Identify 5-10 genes that are either restricted at the front, or have graded expression with highest intensity towards the leading edge; and/or

2) Perform additional functional experiments to see if any of the identified genes is functionally important with a major phenotype.

Please also note that you should include better controls for the morpholino experiments.

We have addressed the final Essential revision based on the reviewers’ request to add in situ data for several Trailblazer genes. We now include expression analysis of 6 more Trailblazer genes (for a total of 13) expressed in but not exclusive in these cells. We are very grateful for the reviewer’s helpful suggestions that have strengthened the manuscript.

[Editors' note: further revisions were requested prior to acceptance, as described below.]

Essential revisions:

1) The authors should aim to provide a clear, and accurate description of their findings. These in vivo findings need to be spelled clearly: These genes were identified by scRNAseq and bioinformatics (through a stringent protocol with harsh significance in mind), as specific to the expected Trailblazer at the narrow front edge of the migratory front, however, they do not seem to appear restricted in that fashion in vivo, and insteada) None of the 13 genes tested so far (which include some of the top candidates) are expressed in all of, or the majority of the cells at the front edge.b) some may appear in a few or some but not in all of the cells (~30% or 1/3?) of the (narrow) front migratory edge.c) some of the genes may display salt and pepper patterns which may include strong expression in a few cells in the narrow front edge and similarly strong signal in some cells in the beginning or middle portions of the stream (not at the front edge of migration).d) Some of the genes where not found at the front edge at all, and instead appear to be expressed more strongly at various positions within the migratory stream.

We have revised the Results and Discussion sections to provide a clear and accurate description of our findings, including a more accurate description of the gene expression analysis and inconsistencies between the RNAseq and in vivo expression analysis data exactly as the reviewer’s describe above.

2) The authors should provide a revised version of the so called "complex profile/signature of expression for a Trailblazer" If the genes are not specifically expressed in all of them, or the majority, and if these genes may be expressed in the stream, then what is the so called complex signature of the trail blazer? (Not to mention its relevance…) The authors could offer some useful insights regarding the limitations or caveats of the scRNAseq at the current time/or in these complex settings. Perhaps cumulative RNA vs. freshly made RNA.

We have clarified this in the text and also provide a list of genes expressed in greater than 95% of trailblazer neural crest cells (Figure 5F) and have clarified the modest changes in neural crest cell distance migrated after knockdown of several trailblazer genes, including the addition of Wnt5B MO and p-values for each of the morpholino experiments (Results section, Figure 6, and Figure 6 legend).

3) The novel RNAscope data shows heterogeneous expression of genes that are part of the "Trailbrazer signature". While this reviewer finds that a thorough transcriptomic characterization of single RNA cells may be a useful resource for the community, the results of the in vivo analysis require further clarification.

We appreciate the reviewer’s helpful suggestions and have revised the text to be clear.

4) Upon presenting the expanded RNAscope analysis, the authors conclude that "the combination of fluorescence in situ hybridization with IHC and tissue clearing allowed the confirmation of a small subset of the Trailblazer transcriptional signature". This statement seems to imply that heterogeneous expression is sufficient to validate transcripts as being part of the Trailblazer signature. According to the Introduction of the manuscript such genes should be expressed by cells "narrowly confined by the migratory front". This discrepancy needs to be addressed or the reader will not be able to link the transcriptomic data with the in vivo expression analysis.

We completely agree with the reviewers and have revised both the Abstract, Introduction and Results to clarify our results and establish a better link between the transcriptomic data and in vivo expression analysis.

https://doi.org/10.7554/eLife.28415.029

Article and author information

Author details

  1. Jason A Morrison

    Stowers Institute for Medical Research, Kansas City, United States
    Contribution
    Conceptualization, Formal analysis, Supervision, Funding acquisition, Validation, Investigation, Visualization, Methodology, Writing—original draft, Project administration, Writing—review and editing
    Contributed equally with
    Rebecca McLennan
    Competing interests
    No competing interests declared
    ORCID icon "This ORCID iD identifies the author of this article:" 0000-0002-1654-8076
  2. Rebecca McLennan

    Stowers Institute for Medical Research, Kansas City, United States
    Contribution
    Conceptualization, Data curation, Formal analysis, Supervision, Validation, Investigation, Visualization, Methodology, Writing—original draft, Project administration, Writing—review and editing
    Contributed equally with
    Jason A Morrison
    Competing interests
    No competing interests declared
  3. Lauren A Wolfe

    Stowers Institute for Medical Research, Kansas City, United States
    Present address
    Department of Biostatistics, Bioinformatics and Epidemiology, Fred Hutchinson Cancer Research Center, Seattle, United States
    Contribution
    Conceptualization, Data curation, Software, Formal analysis, Supervision, Validation, Investigation, Visualization, Methodology, Writing—original draft, Project administration, Writing—review and editing
    Competing interests
    No competing interests declared
  4. Madelaine M Gogol

    Stowers Institute for Medical Research, Kansas City, United States
    Contribution
    Software, Formal analysis, Validation, Methodology, Writing—original draft
    Competing interests
    No competing interests declared
    ORCID icon "This ORCID iD identifies the author of this article:" 0000-0002-8738-0995
  5. Samuel Meier

    Stowers Institute for Medical Research, Kansas City, United States
    Present address
    Broad Institute, Cambridge, United States
    Contribution
    Software, Formal analysis, Validation, Investigation, Methodology
    Competing interests
    No competing interests declared
  6. Mary C McKinney

    Stowers Institute for Medical Research, Kansas City, United States
    Contribution
    Software, Formal analysis, Validation, Investigation, Visualization, Methodology
    Competing interests
    No competing interests declared
  7. Jessica M Teddy

    Stowers Institute for Medical Research, Kansas City, United States
    Contribution
    Software, Formal analysis, Validation, Investigation, Visualization, Methodology
    Competing interests
    No competing interests declared
  8. Laura Holmes

    Stowers Institute for Medical Research, Kansas City, United States
    Contribution
    Data curation, Software, Validation, Investigation, Methodology
    Competing interests
    No competing interests declared
  9. Craig L Semerad

    University of Nebraska Medical Center, Omaha, United States
    Contribution
    Data curation, Software, Validation, Investigation, Methodology
    Competing interests
    No competing interests declared
  10. Andrew C Box

    Stowers Institute for Medical Research, Kansas City, United States
    Contribution
    Data curation, Software, Validation, Investigation, Methodology
    Competing interests
    No competing interests declared
  11. Hua Li

    Stowers Institute for Medical Research, Kansas City, United States
    Contribution
    Data curation, Software, Investigation, Visualization, Methodology
    Competing interests
    No competing interests declared
  12. Kathryn E Hall

    Stowers Institute for Medical Research, Kansas City, United States
    Contribution
    Data curation, Software, Investigation, Visualization, Methodology
    Competing interests
    No competing interests declared
  13. Anoja G Perera

    Stowers Institute for Medical Research, Kansas City, United States
    Contribution
    Data curation, Supervision, Methodology
    Competing interests
    No competing interests declared
  14. Paul M Kulesa

    1. Stowers Institute for Medical Research, Kansas City, United States
    2. Department of Anatomy and Cell Biology, University of Kansas School of Medicine, Kansas City, United States
    Contribution
    Conceptualization, Data curation, Supervision, Funding acquisition, Methodology, Writing—original draft, Project administration, Writing—review and editing
    For correspondence
    pmk@stowers.org
    Competing interests
    No competing interests declared
    ORCID icon "This ORCID iD identifies the author of this article:" 0000-0001-6354-9904

Funding

Stowers Institute for Medical Research

  • Paul M Kulesa

National Institute of Neurological Disorders and Stroke (R21NS092001)

  • Paul M Kulesa

National Institute of Child Health and Human Development (R03HD089190)

  • Paul M Kulesa

The funders had no role in data collection and interpretation, or the decision to submit the work for publication.

Acknowledgements

This research was supported by the generosity of the Stowers Institute for Medical Research and partially supported by the National Institute of Neurological Disorders and Stroke as well as the National Institute of Child Health and Human Development of the National Institutes of Health under Award Numbers R21NS092001 and R03HD089190, respectively. We also thank members of the Flow Cytometry, Molecular Biology and Bioinformatics core facilities at the Stowers Institute. Original data underlying this manuscript can be accessed from the Stowers Original Data Repository at http://www.stowers.org/research/publications/LIBPB-1158.

Ethics

Animal experimentation: All experiments were performed according to institutional (IBC-2003-23-pmk) and federal ethical standards.

Reviewing Editor

  1. Marianne Bronner, California Institute of Technology, United States

Publication history

  1. Received: May 5, 2017
  2. Accepted: December 2, 2017
  3. Accepted Manuscript published: December 4, 2017 (version 1)
  4. Version of Record published: December 13, 2017 (version 2)

Copyright

© 2017, Morrison et al.

This article is distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use and redistribution provided that the original author and source are credited.

Metrics

  • 3,129
    Page views
  • 468
    Downloads
  • 7
    Citations

Article citation count generated by polling the highest count across the following sources: Crossref, Scopus, PubMed Central.

Download links

A two-part list of links to download the article, or parts of the article, in various formats.

Downloads (link to download the article as PDF)

Download citations (links to download the citations from this article in formats compatible with various reference manager tools)

Open citations (links to open the citations from this article in various online reference manager services)