Atlas of RNA sequencing profiles for normal human tissues

Comprehensive analysis of molecular pathology requires a collection of reference samples representing normal tissues from healthy donors. For the available limited collections of normal tissues from postmortal donors, there is a problem of data incompatibility, as different datasets generated using different experimental platforms often cannot be merged in a single panel. Here, we constructed and deposited the gene expression database of normal human tissues based on uniformly screened original sequencing data. In total, 142 solid tissue samples representing 20 organs were taken from post-mortal human healthy donors of different age killed in road accidents no later than 36 hours after death. Blood samples were taken from 17 healthy volunteers. We then compared them with the 758 transcriptomic profiles taken from the other databases. We found that overall 463 biosamples showed tissue-specific rather than platform- or database-specific clustering and could be aggregated in a single database termed Oncobox Atlas of Normal Tissue Expression (ANTE). Our data will be useful to all those working with the analysis of human gene expression.


Background & Summary
High throughput gene expression (transcriptomic) analyses can be used in every aspect of biomedicine 1-3 , including fundamental research [4][5][6] and molecular diagnostics 7 . Growing amount of transcriptomic data is deposited in the special public repositories like Gene Expression Omnibus, GEO 8 and Array-Express 9 which have already accumulated over two million individual transcriptomic profiles obtained in over 100,000 series of experiments 10 . The data cover a wide spectrum of specific human physiological conditions, including most of known diseases and developmental features 8,9,11,12 .
The enclosed transcriptional profiles were obtained using different experimental platforms of either microarray hybridization or next-generation sequencing (NGS) of mRNA. However, the gene expression may be poorly comparable across the different platforms because of use of different equipment and reagents [13][14][15][16][17] . The NGS methods have become the standard in the field of gene expression studies because of higher reproducibility and lower platform bias 18 .
Comprehensive studies of molecular pathology or tissue specific patterns require a collection of reference samples 19,20 . Ideally, they should represent normal tissues from healthy donors, profiled in a single series of experiments using the same equipment and reagents. Nowadays, there is a shortage of such gene expression data. The largest published dataset, GTEx 21 (11,688 samples), lacks publicly available data on the donors' age, and this puts limits on data analysis and makes impossible age-matched studies, which are crucial for many practical applications and aging research. For that reason, we didn't include GTEx data in this study.
The other relevant databases of significant size do include the information on the donors' age: TCGA 22  www.nature.com/scientificdata www.nature.com/scientificdata/ lack one or several of the previously mentioned features. For example, in The Cancer Genome Atlas (TCGA) project database, specimens of histologically normal tissue adjacent to surgically removed tumors 24 are considered normal. However, these tissues may not be completely normal due to numerous effects tumors may have on the neighboring cells, including biased growth factors and cytokine balances 25 , pathological inflammation 26 , and altered vascularization 27 . The ENCODE polyA RNA-seq and ENCODE total RNA-seq datasets were generated  www.nature.com/scientificdata www.nature.com/scientificdata/ based on the normal tissues subjected to NGS using different library preparation methods. They have only 1-4 samples profiled per tissue type (both male and female donors included) and most of the cases can not form a statistically significant reference group.  www.nature.com/scientificdata www.nature.com/scientificdata/ During our study we have designed and assembled the gene expression database of normal human tissues based on uniformly obtained original sequencing data using Illumina HiSeq-3000 engine. A total of 142 solid tissue samples representing 20 organs were taken from post-mortem human healthy donors, who had died in road accidents, no later than 36 hours after death. Blood samples were taken from 17 healthy volunteers. The materials, that had been being collected since 2012, were stored until gene expression profiles were obtained in one series of experiments using the same reagents and protocols. The gene expression profiles were then submitted to Gene Expression Omnibus under accession id GSE120795 28 . We also compared our data for consistency with the transcriptomic profiles taken from databases The Cancer Genome Atlas (TCGA), ENCODE polyA RNA-seq, and ENCODE-total  www.nature.com/scientificdata www.nature.com/scientificdata/ RNA-seq. We found that 463 biosamples showed tissue-specific rather than platform-or database-specific clustering. These have beeen aggregated in a single database named Oncobox Atlas of Normal Tissue Expression (ANTE) (can be found on Figshare 29 in the file "ANTE overview"), including 11 sex-matched statistically significant tissue groups. Our data will be useful to all those working with the analysis of human gene expression.

Methods
Biosamples. In the period from March 2012 to May 2018, we obtained 142 solid tissue samples from 20 human organs at the Department of Pathology at the Faculty of Medicine, Moscow State University, Russia, from autopsies taken from 23 non-related adult healthy donors killed in road accidents, no later than 36 hours after death. We also collected normal blood samples from 17 healthy volunteers. For each biosample an informed written consent to participate in the study was obtained from the patient's legal representative. The consent procedure and the design of the study were approved by the ethical committees of the Faculty of Medicine, Moscow State University, and of the National Research Center for Hematology. The biosamples were evaluated by a pathologist to confirm the tissue origin of every specimen prior to analyses. www.nature.com/scientificdata www.nature.com/scientificdata/ Overall, blood and solid tissues were taken from 8 male and 9 female and 14 male and 9 female donors, respectively.  Fig. 7 The dendrogram of normal samples from Oncobox (no prefix in sample names), TCGA database ("TCGA_" prefix in samples names) and ENCODE ("ENCODE_" prefix in samples names). For the TCGA data 10 random samples per tissue type were selected for visualization, in cases when more norms were available. Euclidian distance between the samples was measured using gene expression data. The dendrogram was built using R ward.D2 method. The color markers indicate the tissue type. The lower scales indicate the number of uniquely mapped reads.  www.nature.com/scientificdata www.nature.com/scientificdata/ Library preparation. RNA Integrity Number (RIN) was measured using Agilent 2100 bioanalyzer. Agilent RNA 6000 Nano or Qubit RNA Assay Kits were used to measure RNA concentration. KAPA RNA Hyper with RiboErase (KAPA Biosystem) Kit was used for further depletion of ribosomal RNA and library preparation. Different adaptors were used for multiplexing samples in one sequencing run. Library concentrations and quality were measured using Qubit ds DNA HS Assay kit (Life Technologies) and Agilent Tapestation (Agilent). RNA sequencing was performed using Illumina HiSeq 3000 equipment for single end sequencing, 50 bp read length, for approximately 30 million raw reads per sample. Data quality check was conducted using Illumina SAV. De-multiplexing was performed using Illumina Bcl2fastq2 v 2.17 software.
Processing of RNA sequencing data. RNA sequencing FASTQ files were processed with STAR aligner 30 in 'GeneCounts' mode with the Ensembl human transcriptome annotation (Build version GRCh38 and transcript annotation GRCh38.89). Ensembl gene IDs were converted to HGNC gene symbols using Complete HGNC dataset (https://www.genenames.org, database version of July 13, 2017. In total, expression levels were established for 36596 annotated genes with corresponding HGNC identifiers. Additional QC metrics for obtained data were generated using NCBI MAGIC software 18,31,32 . All metrics and detailed protocol for each sample can be found on Figshare 29 in the file "QC and additional meta-information".
Data clustering. '1' was added to all raw gene counts prior to cluster analyses, to avoid zero expression values, as described by Dillies et al. 33 , the gene expression data were merged into single datasets and quantile normalized 34 . Hierarchical clustering was performed using R ward.D2 method. The dendrogram was visualized using custom R script.

Data Records
Gene expression profiles were deposited at the Gene Expression Omnibus database (GEO) under accession number GSE120795 28 . The data is provided as a matrix of raw counts as generated by STAR. The raw data can be downloaded from the Sequence Read Archive (SRA 35 ). The mapping statistics for corresponding dataset can be found on Figshare 29 in the file "QC and additional meta information". The available RNA sequencing profiles for normal and cancer tissues matched normal samples were extracted from the websites of projects TCGA and ENCODE (portal.gdc.cancer.gov and www.encodeproject.org, respectively). Combined meta information for compatible gene expression profiles from different databases can be found on Figshare 29 in the file "ANTE overview".

technical Validation
RNa sequencing data quality control and consistency tests. To assess if the obtained gene expression profiles are in correlation with the biological nature of the biosamples tested and identify samples that might be of low quality, we performed cluster analysis of all the RNA sequencing data we had obtained (Fig. 1). The color on the dendrogram indicates technical type of biosamples: FFPE or tissue in RNAlater for solid tissue, and RNA in ethanol for blood samples.
The FFPE and RNAlater samples formed several mixed clusters on the dendrogram, which suggests that the clustering was independent of the sample preparation technique for solid tissues. The blood samples formed a clearly distinct cluster on the left side of the dendrogram, which is in line with the histological status of this tissue (Fig. 1), which stand apart among other samples.
However, we then observed a mixed cluster next to the one formed by the blood samples, which included both blood and solid tissue specimens (Fig. 1, red), and such clustering did not correspond to physiological origin of biosamples. We noticed that that cluster was formed exclusively by the samples with relatively low (less than 2.5 million) number of uniquely mapped sequencing reads (Fig. 1, scale) and hypothesized that this may represent a deviation which arose due to insufficiency of data. We analyzed contents of the samples with respect to the number of uniquely mapped reads (Fig. 2) and found that the samples with lowest number of reads indeed formed a distinct cluster with the upper threshold of ~2.5 million reads uniquely mapped to HGNC genes (Fig. 2). We, therefore, used this threshold, which enabled us to separate effectively solid tumors from hematological ones, as an indicator for quality control (QC) of the sequenced gene expression profiles. 132 samples out of a total of 159 passed the QC threshold (Table 1). After the QC filter was applied RNA sequencing profiles became clustered in the hierarchical way following the histological origin of the tissues (Fig. 3).
The established threshold effectively marked samples with low quality values of other QC metrics, e.g. proportion of genomic counts, high rate of mismatches, number of reads spanning splice junction, high percentage of ribosomal counts. Full list of mapping statistics generated by STAR aligner and additional QC metrics generated by NCBI MAGIC is available on Figshare 29 in the file "QC and additional meta information".
We then looked for a correlation between whether a biosample passes the QC and its internal characteristics. The following parameters were investigated: (i) RIN, that shows the level of RNA degradation (lower RIN points to more degraded RNA). All samples with high RIN (RIN > 4) passed the QC (Fig. 4). However, low RIN turned out not to be an informative marker of the insufficient number of reads, and most of the samples with 1 < RIN < 2 passed the QC as well. (ii) RNA concentration. We found no correlation between RNA concentration and number of uniquely mapped reads (Fig. 5). Therefore, low RNA concentration was not a significant indicator for quality assessment. In most of the cases, low RNA concentration still allowed the sample to pass the QC.
Only the samples that had more reads than the QC threshold were analyzed in further comparisons.
www.nature.com/scientificdata www.nature.com/scientificdata/ Finally, we assessed the reproducibility of the RNA sequencing output data for the characterization of normal human tissues. We performed RNA sequencing of four different fragments of each of the following tissue specimens: two human liver (ID 16_5) and esophagus (ID E_3). The tissue fragments were blinded and sent for sequencing separately. In each case, one of quadruplicates was sequenced and mapped in a separate batch. In both cases we observed high pairwise correlation coefficients between quantile normalized gene expression values (Spearman's rho ≥ 0.94 in all cases) (Fig. 6a,b). We concluded that the obtained gene expression profiles are highly reproducible among the replicates for materials taken from different portions of the same biosample.
Compatibility with other gene expression datasets. After that we checked if our gene expression data are compatible with previously published transcriptomic datasets. The largest published datasets of RNA sequencing data for normal tissues were TCGA database 22 of 'normal' tissues adjacent to tumors (625 samples, see the file "TCGA gene expression data" on Figshare 29 ), ENCODE 23 database, poly(A) priming library preparation (38 samples, see the file "ENCODE polyA RNA-seq gene expression data" on Figshare 29 ), and ENCODE 23 database, random priming library preparation (see the file "ENCODE total RNA-seq gene expression data" on Figshare 29 ). All the above mentioned gene expression profiles, both obtained in our experiments and published earlier, were merged and quantile normalized 34 . We then performed cluster analysis of these data and found that both ENCODE databases showed clustering similar to that of the current experimental database, with clear tissue specific, rather than database specific, clustering patterns (Fig. 7). Contrariwise, most of the samples from the TCGA database clustered in a database specific way (Fig. 7), probably due to cancer-associated pathological changes in the respective tissue specimens. However, the TCGA samples of lung, liver, kidney, thyroid gland, adrenal gland, and brain clustered according to the tissue they represent along with the samples from other datasets (Fig. 7). We then assembled all the biosamples in all four datasets, that clustered in the tissue-specific way, in a single database named Oncobox Atlas of Normal Tissue Expression (ANTE). The database overview is presented in the Table 2, the detailed information can be found on Figshare 29 in the file "ANTE overview". In total, the ANTE database contains 463 transcriptomic profiles for 11 tissue types. All the types of tissue are represented by seven samples or more, with the mean representation of 42 samples per tissue.

Usage Notes
We have assembled the gene expression database of normal human tissues termed Oncobox Atlas of Normal Tissue Expression (ANTE). The database includes 67 original experimental and 396 previously published gene expression profiles for 11 normal human tissues. The experimental profiles that are published here for the first time can be considered a golden standard data from the histopathological point of view as they had been obtained for the post-mortem human healthy donors, killed in road accidents, no later than 36 hours after death. Blood samples were taken from 17 healthy volunteers. The ANTE database provides unparalleled source of normal human tissues for transcriptomic research. The data in the database are in a machine-readable format and annotated by the available patients' data which include tissue type, sex, and age of a donor. Besides comparisons of individual gene expressions, e.g. of pathological v normal human tissues, this collection of molecular data could be beneficial for a more complex level of data analysis such as molecular pathway activation scoring 36,37 and bioinformatic ranking of drugs 38 . Additionally, the database can be used to revise lists of tissue specific gene expression biomarkers and housekeeping genes.