Transcriptomics

HPA RNA-seq data overview

In total, 1206 cell lines, 33 human tissues, 192 samples from micro-dissected areas and regions of the human brain, and 18 immune cell types as well as total peripheral blood mononuclear cells (PBMC) have been analyzed by RNA-seq to estimate the transcript abundance of each protein-coding gene.

Immune cell transcriptomics

All blood samples were collected with consent from patients and all samples were anonymized in accordance with approval from the local ethics committee (ref #2015/1552-32) and Swedish rules and legislation. Blood samples were enriched for PBMC and granulocytes, labeled with antibodies and separated into subpopulation by flow sorting.

For a total number of 109 immune cell samples mRNA sequencing was performed on an Illumina NovaSeq 6000 System in four S4 lanes with a read length of 2x150 bases. Transcript abundance estimation was performed using Kallisto v0.48.0. The 18 immune cell types are classified into six different lineages including B-cells, T-cells, NK-cells, monocytes, granulocytes and dendritic cells. More information can be found here.

Normalization of transcriptomics data

For both the HPA and GTEx transcriptomics datasets, the average TPM value of all individual samples for each human tissue or human cell type was used to estimate the gene expression level. To be able to combine the datasets into consensus transcript expression levels, a pipeline was set up to normalize the data for all samples. In brief, all TPM values per sample were scaled to a sum of 1 million TPM (denoted pTPM) to compensate for the non-coding transcripts that had been previously removed. Next, all TPM values of all samples within each data source (HPA + GTEx human tissues, HPA immune cell types, HPA cell lines) were normalized separately using Trimmed mean of M values (TMM) to allow for between-sample comparisons. The resulting normalized transcript expression values, denoted nTPM, were calculated for each gene in every sample. nTPM values below 0.1 are not visualized on the Atlas sections.

For the brain dataset, an additional normalization was performed using linear regression to do the correction for inter-individual variation using the removeBatchEffect in the R package Limma with subject as a batch parameter. To reduce the technical variation between MGI and illumina platforms, 19 reference samples were included and run on both platforms. Intensity normalization based on reference samples was conducted to minimize technical variation between two platforms.

Consensus transcript expression levels for each gene were summarized in 51 human tissues based on transcriptomics data from the two sources HPA and GTEx. The consensus nTPM value for each gene and tissue type represents the maximum nTPM value based on HPA and GTEx. For tissues with multiple sub-tissues (brain regions, immune cells, lymphoid tissues and intestine) the maximum of all sub-tissues is used for the tissue type and the total number of tissue types in the human tissue consensus set is 37.

The FANTOM5 dataset was normalized separately on the sample level using TMM. The normalized Tags Per Million for each gene were calculated based on the average of all individual samples for each human tissue.

Classification of transcriptomics data

The consensus transcriptomics data was used to classify all genes according to their tissue-specific, single cell type-specific, brain region-specific, immune cell-specific or cell line-specific expression into two different schemas: specificity category and distribution category. These are defined based on the total set of all nTPM values in 40 tissues, 155 single cell types, 13 main regions of each mammalian brain,18 immune cell types or 1132 cell lines grouped into 28 cancer types and using a cutoff value of 1 nTPM as a limit for detection across all tissues or cell types.

Explanation of the specificity category

Category	Description
Enriched	nTPM in a particular tissue/region/cell type is at least fourfold higher than in any other tissue/region/cell type
Group enriched	Mean nTPM across a group (of 2-5 tissues, brain regions, or cell lines, or 2-10 immune cell types) is at least fourfold higher than in any other tissue/region/cell line/immune cell type/cell type
Enhanced	nTPM in one or several tissues, brain regions, cell lines or immune cell types is at least fourfold higher than the mean nTPM across all tissues, provided the lowest-expressing enhanced tissue/brain region etc. is at least twofold higher than the average of the next three tissues (a gap ratio of ≥ 2).
Low specificity	nTPM ≥ 1 in at least one tissue/region/cell type but not elevated in any tissue/region/cell type
Not detected	nTPM < 1 in all tissue/region/cell types

Gap ratio

The enhanced criterion compares expression with the mean nTPM across all tissues or cell types. This mean depends on how many tissues or cell types are in the denominator. With many cell types, most of them express a given gene at low levels, which dilutes the mean and makes the fourfold threshold easy to pass. A gene expressed at the same level in k of N cell types passes it whenever k ≤ N/4. Across the 155 single cell types, a gene with broad, gradually declining expression could therefore be called enhanced in dozens of cell types.

The gap ratio prevents this by requiring a clear drop in expression directly after the enhanced set:

gap ratio = nTPM of the lowest-expressing enhanced tissue or cell type / mean nTPM of the next three

A gene is classified as enhanced only if it meets both criteria: at least fourfold higher than the mean, and a gap ratio of ≥ 2. Otherwise it is classified as low specificity.

ACP1 is expressed in many cell types, but spermatogenic cell types stand out at 1001.9 nCPM. That is more than fourfold the mean across all cell types (133.0 nCPM). The next three cell types average 375.2 nCPM, which gives a gap ratio of 2.67. This clear drop after the enhanced cell type means ACP1 is classified as enhanced in spermatogenic cell types.

PACSIN3 shows a gradual decline in expression across many cell types, so the mean (25.7 nCPM) is low and the fourfold criterion is easy to pass. Stromal cells (110.4 nCPM) and basal cells (105.0 nCPM) exceed this threshold. However, the next three cell types follow closely, at 102.6, 81.0 and 72.1 nCPM, giving a gap ratio of 105.0 / 85.2 = 1.23. With no clear break between the candidate cell types and the rest, PACSIN3 is not classified as enhanced.

Elevated genes

The combined category "elevated", containing all genes in the first three categories (tissue/cell line/cell type enriched, group enriched and tissue/cell line/cell type enhanced), has been used for some parts of the analysis. TS/CS-score (Tissue Specificity/Cell Specificity score) is calculated for “elevated” tissues/cell lines. TS/CS-score is calculated as the fold change from the tissue/cell line with highest RNA to the tissue/cell line with second highest RNA.

Explanation of the distribution category

Category	Description
Detected in single	Detected in a single tissue/region/cell type
Detected in some	Detected in more than one but less than one third of tissues/regions/cell types
Detected in many	Detected in at least a third but not all tissues/regions/cell types
Detected in all	Detected in all tissues/regions/cell types
Not detected	nTPM < 1 in all tissues/regions/cell types

Gene clustering of transcriptomics data

The RNA expression data has been used to classify protein-coding genes into expression clusters for tissues, single cell types, immune cells, and cell lines.

Clustering	Number of tissues, cell types or cell lines	Sample aggregation level
Tissue	78	Averaged nTPM expression per tissue type (40 HPA and 38 GTEX tissue types)
Single cell type	1214	Averaged nCPM expression per cell type cluster
Cell lines	1206	nTPM expression of individual cell lines
Immune cells	103	Averaged nTPM expression per immune cell
Brain	192	Averaged nTPM expression per brain region

Data preprocessing

For each dataset, genes with expression level > 1 in at least one of the samples were selected. The data was genewise scaled to Z-scores to account for differences in dynamic ranges between genes across samples. After scaling, the expression data was projected into a lower dimensional space using Principal Component Analysis (PCA), where a number of components were selected to satisfy Kaiser’s rule (eigenvalue ≥ 1) and explaining at least 80% of the total variance. Gene to gene distances were calculated as the Spearman correlation of gene expression across samples, and transformed to Spearman distance (1 - Spearman correlation).

Gene clustering

Based on the distances, a k-nearest neighbors (kNN) graph was computeted based on 20 nearest neighbors, which was subsequently to find clusters of similarly expressed genes via Louvain clustering. To account for the stochasticity in the louvain algorithim, the clustering was performed 100 times. The results were later collapsed into a single consensus clustering. Confidence of the gene-to-cluster assignment was calculated as the fraction of times that the gene was assigned to the cluster.

Cluster annotation

The clustering generated for each of the datasets is manually annotated to assign a specificity and function to each cluster. The annotation is based on overrepresentation analysis towards biological databases, including Gene Ontology, Reactome, PanglaoDB, TRRUST, and KEGG, as well as HPA classifications including subcellular location, protein class, secretion location and classification, and specificity toward tissues, single cell types, immune cells, brain regions, and cell lines. A reliability score is manually set for each cluster indicating the confidence of specificity and function assignment.

Clustering visualization

The clustering results are visualized in a UMAP. Colored polygons were generated to represent the main contiguous masses of genes corresponding to the same cluster. First, for each cluster, the two-dimensional density was estimated in the UMAP, and an area enveloping 95% of the total density was determined. The areas were moderated to include contiguous areas corresponding to at least 5% of the total area in the UMAP space. Finally, contiguous areas were converted to two-dimensional polygons per each cluster.

External immune cell RNA-seq data

In addition to the immune cell type data from blood, generated within the Human Protein Atlas project, data from 15 immune cell types by Schmiedel et al. and 29 immune cell types as well as total PBMC by Monaco et al. have been incorporated into the Blood Atlas.

The Schmiedel dataset is available at the DICE (Database of Immune Cell Expression, Expression quantitative trait loci (eQTLs) and Epigenomics) database, which was established to address how genetic variants associated with risk for human diseases affect gene expression in various cell types. The TPM values per gene for 15 immune cell types were mapped to the corresponding genes in the Ensembl version used in the Human Protein Atlas.

The Monaco dataset contains data for 29 immune cell types within the peripheral blood mononuclear cell (PBMC) fraction of healthy donors using RNA-seq and flow cytometry. Raw data for 29 immune cells as well as total PBMC were analyzed using the same pipeline as for HPA-generated RNA-seq data and also normalized using TMM to allow for between-sample comparisons. Normalized gene expression values are reported as nTPM values.