Pseudotemporal Analysis of AER Cell Development

Pseudotemporal ordering of AER cells, forelimb or hindlimb was done with Monocle 2⁵⁷ . Briefly, differentially expressed genes across five development stages were identified with the differentialGeneTest function of Monocle 2⁵⁷ . The top 500 genes with the lowest q value were used to construct the pseudotime trajectory using Monocle 2⁵⁷ , with UMI count per cell as a covariate in the tree construction. Each cell was assigned a pseudotime value based on its position along the trajectory. Smoothed gene marker expression change along pseudotime were generated by plot_genes_in_pseudotim function in Monocle 2⁵⁷ . Cells in the trajectory were grouped in the same method as a previous study⁶⁴ . Briefly, cells were grouped first at similar positions in pseudotime by k-means clustering along the pseudotime axis (k = 10). These clusters were subdivided into groups containing at least 50 and no more than 100 cells. We then aggregated the transcriptome profiles of cells within each group. The gene expression along pseudotime was calculated in the same approach as a previous study⁶⁴ . Briefly, genes passing significant test (FDR of 5%) across different treatment conditions were selected and a natural spline was used to fit the gene expression along pseudotime, with mean_number_genes included as a covariate. The gene expression for each gene was subtracted by the lowest expression and then divided by the highest expression. Genes with max expression within the early 20% of pseudotime were labeled as repressed genes. Genes with max expression in the last 20% of pseudotime were labeled as activated genes. Other genes were labeled as transient genes. Enriched reactome terms (Reactome_2016) and transcription factors (ChEA_2016) were identified using EnrichR/v1.0 package⁶⁵ .