in this study, we built a large - scale integrated cell atlas of human neural organoids, the hnoca, by integrating 1. 8 million cells spanning 36 scrna - seq datasets generated by 15 different laboratories worldwide using 26 different differentiation protocols as well as diverse scrna - seq technologies. the resulting atlas revealed the high complexity of neuronal, glial and non - neural cell types that can develop in neural organoids grown under existing protocol conditions. mapping the hnoca data to various human developing brain cell reference atlases allowed comprehensive evaluation of neural organoid protocols to generate cell types of different brain regions. we found that organoids in the first 3 months of culture best match to first - trimester primary data, whereas organoids around 3 months of culture and older best match second - trimester primary cell states. we did not observe significant neuronal maturation and diversification signatures matching older developmental stages, suggesting a limitation of neuronal maturation in current neural organoid protocols. we performed de analysis between organoid neuron types and their primary counterparts to evaluate transcriptomic fidelity, and identified metabolic changes related to the glyco lysis pathway as a main factor that disting uishes organoid and primary cell states, consistent with previous reports. despite the negative effects of metabolic stress on overall transcriptomic fidelity, the molecular identity of regional cell types is maintained as evidenced by transcription factor coexpression patterns that are highly consistent with primary counterparts. we showcased the mapping of query data, a recently published single - cell transcriptomic neural organoid morphogen screen, to the hnoca and the primary reference, which enabled a refined cell type annotation, as well as a compositional comparison with existing neural organoid datasets. our powerful framework will facilitate quantitative and comparative analysis of scrna - seq data of human neural organoids, and for the benchmarking of new neural organoid protocols. consistent with earlier reports, we find that unguided protocols generate neural cells with high brain regional variability, which is useful when studying broader fate determination during neurodevelopment. guided protocols resulted in a strong enrichment of the targeted brain regions. we also note that some guided protocols, particularly those targeting midbrain, show relatively low specificity and generate neural cells from the nearby brain regions. this issue may be due to a differential response of neural stem cells in the organoid to the same morphogen cue, or to the lack of a full understanding of the timing, concentration and combinations of morphogens required to precisely define cells of the deeper regions in the central nervous system. the integrated hnoca is also an excellent resource for analysis of disease - modelling neural organoid data. it facilitates cell type annotation and provides a large control cohort of single - cell transcriptomes for comparison. for example, we observed discrepancy of cell type and regional composition between control and disease model samples in many studies. at the same time, the hnoca provides the opportunity to identify disease-specific molecular features against a multi - line multi - protocol large - scale control cohort. we demonstrate how the hnoca can be extended and updated by projecting extra single - cell transcriptomic data of neural organoids to the atlas. further, we have developed a computational toolkit, hnoca - tools, which will enable other researchers to recapitulate the analytic framework applied in our study. together, we imagine that the hnoca will be kept up to date and continue to reflect the landscape of human neural cell states generated in organoids in vitro, serving as a living resource for the neural organoid community that enables the assessment of organoid fidelity, the characterization of perturbed and diseased states and the development of new protocols. we included 33 human neural organoid data from a total of 25 publications plus three unpublished datasets in our atlas ( supplementary table 1 ). we curate d all neural organoid datasets used in this study through the sfaira framework ( github dev branch, 18 april 2023 ). for this, we obtained scrna - seq count matrices and associated metadata from the location provided in the data availability section for every included publication or directly from the authors in case of unpublished data. we harmonized metadata according to the sfaira standards ( https: / / sfaira. readthedocs. io / en / latest / adding _ datasets. html ) and manually curate d an extra metadata column organoid _ age _ days, which described the number of days the organoid had been in culture before collection. we next removed any non - applicable subsets of the published datasets: diseased samples or samples expressing disease-associated mutations ( refs. ), primary fetal data ( refs. ), data collected before neural induction ( refs. ) and share - seq data ( ref. ). we harmonized all remaining datasets to a common feature space using any genes of the biotype ‘ protein _ coding ’ or ‘ lncrna ’ from ensembl release 104 while filling any genes missing in a dataset with zero counts. on average, 50 % of the full gene space ( 36, 842 genes ) was reported in each of the constituent datasets. we then concatenated all remaining datasets to create a single anndata object. all processing and analyses were carried out using scanpy ( v. 1. 9. 3 ) unless indicated otherwise. for quality control and filtering of hnoca, we removed any cells with fewer than 200 genes expressed. we next removed outlier cells in terms of two quality control metrics: the number of expressed genes and percentage mitochondrial counts. to define outlier cells on the basis of each quality control metric, z - transformation is first applied to values across all cells. cells with any z - transformed metric less than −1. 96 or greater than 1. 96 are defined as outliers. for any dataset collected using the v. 3 chemistry by 10x genomics, which contains more than 500 cells after the filtering, we fitted a gaussian distribution to the histogram denoting the number of expressed genes per cell. if a bimodal distribution was detected, we removed any cell with fewer genes expressed than defined by the valley between the two maxima of the distribution. we then normalized the raw read counts for all smart - seq2 data by dividing it by the maximum gene length for each gene obtained from biomart. we next multiplied these normalized read counts by the median gene length across all genes in the datasets and treated those length - normalized counts equivalently to raw counts from the datasets obtained with the help of unique molecular identifiers in our downstream analyses. as a next step we generated a log - normalized expression matrix by first dividing the counts for each cell by the total counts in that cell and multiplying by a factor of 1, 000, 000 before taking the natural logarithm of each count + 1. we computed 3, 000 highly variable features in a batch - aware manner using the scanpy highly _ variable _ genes function ( flavor = ‘ seurat _ v3 ’, batch _ key = ‘ bio _ sample ’ ). here, bio _ sample represents biological samples as provided in the original metadata of the datasets. on average, 72 % of the 3, 000 highly variable genes were reported in each of the constituent hnoca datasets. we used these 3, 000 features to compute a 50 - dimensional representation of the data using principal component analysis ( pca ), which in turn we used to compute a k - nearest - neighbour ( knn ) graph ( n _ neighbors = 30, metric = ‘ cosine ’ ).