Mapping to hg19 human genome, quality control and read counting of Ensembl genes was performed by cellranger software with default parameter (v1.3.0) (Figure S3A). Data visualization was performed by cellrangerRkit (v1.1.0). Normalization, dimensionality reduction and clustering of single cells were performed by Seurat (v 1.4.0.14) in Bioconductor (Macosko et al., 2015). Briefly, read count was first normalized per cell and transformed to log2 scale by Setup function. Cells with more than 200 genes and genes detected in more than 3 cells were retained for subsequent analyses. Negative-binomial regression was then performed to variable genes (more than zero dispersion) using batch and UMI as confounding variables. Principal component analysis was conducted using the regressed variable genes by PCAFast function in Seurat package with “pcs.compute=20” parameter. Using top five principal components (singular value > 1,000), we reduced the dimensionality by t-Distributed Stochastic Neighbor Embedding (tSNE). Cells were then grouped by a shared nearest neighbor modularity optimization using FindClusters function with four resolution value. To compare transcriptional similarity of clusters, classification tree was constructed based on average expression of the set of variable
