The Evolutionary Resource - Methods summary

Data Overview

Under the evolutionary resource, we present primary transcriptomic data spanning the whole body of four mammalian species: Macaque (Macaca fascicularis), Mouse (C57BL/6, Mus musculus), Rat (Sprague Dawley, Rattus norvegicus), and Pig (Chinese Bama minipigs, Sus scrofa). For each species, the sample collection covered four healthy adult organisms: 2 females and 2 males, resulting in 3-4 replicates per sampled tissue.

Samples were collected in accordance with approval from the local ethics committee: Mouse & Rat (Stockholms djurförsöksetiska nämnd N631-19), Pig (Experimental animal permission number: SYXK2017-0123. Ethical permission number: 44410500000078 and BGI-IRB18135), Macaque (Animal Care and Use Committee of Huazhen Bioscience: HZ2019027, Institutional Review Board on Ethics Committee of Beijing Genomics Institute: BGI-IRB 19125-T2 and BGI-IRB 21136)

Data Processing

For macaque, mouse, rat, and pig samples, total RNA was extracted from snap-frozen tissue samples using TRIzol reagent (Invitrogen) and evaluated for integrity. Samples went through ribosomal RNA depletion followed by library construction using MGIEasy RNA library kits. Paired-end sequencing ($2 \times 100$ bp) was performed across DNBSEQ platforms.

The human data used under the evolutionary resource is based on the Tissue Resource transcriptomic data. Details regarding the human data can be found here.

Transcript abundance was quantified using Kallisto (v. 0.50.1) against Ensembl v109 annotation of the species' cDNA references (Macaca_fascicularis_6.0, GRCm39, GRCr8, and Sscrofa11.1, respectively for macaque, mouse, rat, and pig). The resulting quantitative transcript TPM (Transcript per Million) data are filtered for protein-coding transcripts, collapsed into gene-level counts, and rescaled to one million to obtain the protein coding TPM counts (pTPM). The pTPM values were then TMM normalized, resulting in the final normalized TPM values (nTPM).

As protein coding genes, we consider those genes flagged as protein coding ("protein_coding", "protein_coding_LoF"), immunoglobulin genes ("IG_V_gene", "IG_D_gene", "IG_C_gene", "IG_J_gene"), and T-cell receptor genes ("TR_V_gene", "TR_J_gene", "TR_D_gene", "TR_C_gene") by Ensembl. TMM normalization (Trimmed means of M values) is based on NOISeq’s (v.2.46.0) implementation, specifying refColumn to median, logratioTrim to 0.3, sumTrim to 0.3, and doWeighting to FALSE.

Tissue Hierarchical Organisation

To facilitate effective data presentation and comparison, we make use of a three-tier tissue organisation system. Tissues are grouped into Grouped Tissues, which are organized into Organ Systems. A Tissue expression profile is obtained by mean pooling its constituent samples, while a Grouped tissue profile is obtained by max pooling the constituent Tissue profiles.


Organization of Brain Tissues

Brain tissues follow a dedicated regional hierarchy wherein each individual brain tissue maps to one of 13 brain regions (listed below). Regional profiles are derived by max pooling the constituent Tissue profiles.



Intra-Species Specificity Annotation

For every species, based on the complete repertoire of available grouped tissues, each gene was assigned a tissue specificity category and a tau score as defined below. Similarly, based on the available brain regions, a brain region specificity was assigned to each gene, together with a brain regional Tau score. Additionally, an analogous regional specificity category and brain tau score were calculated based on the expression profiles of the 13 defined brain regions.

Intra-Species Specificity Category

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


Genes assigned to any of the first three categories (enriched, group enriched and enhanced), are collectively also referred to as "elevated".

Tau Score

A tau (τ) specificity score was calculated for each protein-coding gene following the formula proposed by Yanai et al. Formally, its calculation is presented below. The calculation uses the loge(nTPM+1) transformed expression levels, where n represents the total number of tissues, xi represents the loge(nTPM+1) expression level in tissuei, and x^i represents the expression in tissue i normalized by the maximum expression across all n tissues.

τ = ∑ i = 1 n 1 - x ^ i n - 1 ; x ^ i = x i max 1 ≤ j ≤ n x j

Orthologs and Orthogroups

Pairwise non-human-to-human orthologs defined by Ensembl Compara (v109) were filtered to retain all one-to-one relationships, alongside select one-to-many pairs that could be resolved to map to a unique human partner based on Ensembl's high-confidence flag (is_high_confidence = 1). These pairwise links were subsequently clustered around each human gene to construct multi-species orthogroups.

For the cross-species annotations explained below, only one-to-one orthogroups were considered (n = 18084), that is, orthogroups which are made up of a single gene per species).

Phylogenetic Gene Trees

Evolutionary gene trees were constructed for all cross-species orthogroups in three steps:

Multiple Sequence Alignment: The canonical amino-acid sequences defined by Ensembl were chosen for each member gene. These were then aligned using MAFFT (v7.526) with the --auto and --anysymbol parameters.

Gene Tree Inference: For orthogroups with ≥ 3 genes, maximum-likelihood trees were inferred using FastTree (v2.2.0) under the Le-Gascuel amino acid substitution model with a 20-category CAT approximation and local Shimodaira-Hasegawa (SH)-like branch support. Two-gene families were assigned canonical two-taxon trees.

Species-Tree Reconciliation: Raw gene trees are reconciled against the species reference tree (((Macaque:1,Human:1)1:1,(Mouse:1,Rat:1)1:1)1:1,Pig:1); using Treerecs (v1.2). Trees were parsimoniously re-rooted (-r), and branches with support values < 0.8 were collapsed and resolved against the species tree topology (-t 0.8). Raw FastTree outputs were retained as fallbacks.

Cross-Species Comparison

Similar to the tissue and regional specificity annotation, we performed two iterations of the Cross-Species comparison explained below.One based upon the tissue view, and a second one specific at brain regional level. We explain below what this entailed.

Defining Comparable Tissues

To compare the overlap and conservation of gene expression and overall tissue specificty, we picked a selection of tissues available across all 5 species. At the tissue level comparison we take the human RNA data from the HPA's Tissue resource. Based on the profile of these 31 tissues, we computed the specificity categories as mentioned above for each species.

For the brain regional comparison, we settled on taking all 13 regions. Here, we make use of the human brain RNA data from the Brain resource. To be noted, however the human brain dataset is lacks olfactory bulb, while the human and mouse dataset is missing a dedicated white matter profile.

Pairwise Tissue Specificity Overlap

We compare between species whether a gene has preserved its specificity status (Tissue Elevated, Low Tissue Specificity, or Not Detected). For elevated genes particularly, we assign a numerical score reflecting the ratio of tissues towards which the specificity was retained.

We describe the categorical overlap through a numerical pairwise specificity overlap score (PSO), a Jaccard similarity index computed as the number of intersecting elevated tissues divided by the total number of available elevated tissues between both human and target species:

PSO ( A ref , A tgt ) = | A ref ∩ A tgt | | A ref ∪ A tgt |

where Aref and Atgt denote the sets of elevated tissues identified in the reference (human) and target species respectively.

Each comparable human gene is assigned one of seven categories:

Category Description PSO
Complete overlap Identical set of elevated tissues overlap, preserved low tissue specificity, or preserved not detected in human and target species. 100%
Partial overlap At least 40% of elevated tissues are shared but below complete overlap (e.g. ≥2 intersecting tissues out of 5 tissues). 40-99%
Minimal overlap The specificity is retained at a minimal level (e.g. 1 intersecting tissue out of 3 tissues). 1-39%
Tissue elevation shift Elevated in both species but in different tissues. 0%
Human specific tissue elevated Elevation only in human but low tissue specificity or not detected in target species 0%
Non-Human tissue elevated Elevation only in the target species but low tissue specificity or not detected in human. 0%
Inconsistent detection Not detected in human and low tissue specificity in target or vice versa 0%

A seperate fallback category, "Unasessed (Missing tissue)", is assigned in case a gene is enrielevated in a tissue that was not sampled in one of the compared species. More practically, this category is relevant for the brain specific comparison, where genes elevated in olfactory bulb or white matter of a target species, will deliver "Unasessed", given we are lacking these tissues in the human dataset.


Tissue Specificity Conservation Category

We summarise the multiple pairwise categories into a single unified specificity conservation category. Each gene is assigned one based on the pairwise specificity overlap score.

Category Description
Complete conservation Min(PSO) = 100%; either identical set of elevated tissues, low tissue specificity, or not detected across all species.
Partial conservation Either median(PSO) or mean (PSO) > 37.5%, but fails complete conservation.
Minimal conservation Max(PSO) > 0% but fails complete or partial conservation.
Tissue elevation shift Max(PSO) = 0%, while human is elevated and at least one species shows shifted elevation.
Human specific tissue elevated Max(PSO) = 0%, while human is elevated, but no other species shows elevation.
Non-Human tissue elevated Max(PSO) = 0%, while ≥1 target species is elevated, but is not elevated in human.
Human specific detected Max(PSO) = 0%, while human shows low specificity, while all target species show not detected.
Human not detected Max(PSO) = 0%, while human is not detected and targets show low tissue specificity

Genes are assigned as "Unassessed", in case all pairwise overlaps resulted as "Unasessed (Missing tissue)". This occus when a gene is consistently elevated in a tissue that was not sampled in human. In practice, this is relevant for the brain comparison, where genes elevated in olfactory bulb or white matter across all species will return as "Unassessed".

We performed a random permutation benchmark to define the partial conservation PSO threshold of 37.5.

Cross-Species Cosine Similarity

To complement the discrete specificity categories, we introduce a continuous numerical metric to assess a gene's expression profile conservation across species. For each ortholog pair, we conducted pairwise comparisons between the human reference and the target mammalian species across the coinciding tissues.

First, a pairwise cosine similarity score (Scos) is computed between expression vectors u (human) and v (target species):

Scos ( u , v ) = ∑ t ∈ T u t v t ∑ t ∈ T u t 2 ⋅ ∑ t ∈ T v t 2

where ut and vt represent the expression levels in tissue t, and T is the set of shared tissues between both species.

Secondly, the pairwise human to target scores were then averaged to return one human focused cross-speces cosine similairty score per gene. This metric ranges from 0 (orthogonal non conserved profiles) to 1 (identical relative expression distribution across tissues). T denotes the set of coinciding tissues.