Identifying functional DNA (and junk) by purifying selection

Functional DNA is best defined as DNA that is currently under purifying selection. In other words, it can't be deleted without affecting the fitness of the individual. This is the "maintenance function" definition and it differs from the "causal role" and "selected effect" definitions [The Function Wars Part IX: Stefan Linquist on Causal Role vs Selected Effect].

It has always been difficult to determine whether a given sequence is under purifying selection so sequence conservation is often used as a proxy. This is perfectly justifiable since the two criteria are strongly correlated. As a general rule, sequences that are currently being maintained by selection are ancient enough to show evidence of conservation. The only exceptions are de novo sequences and sequences that have recently become expendable and these are rare.

However, in the ongoing dispute over junk DNA, sequence conservation is often ignored as a justification for function on the grounds that new functions are not rare but exceedingly numerous. For example, transcripts that are unique to humans would seem to be good candidates for noisy transcription or junk RNA but opponents of junk DNA claim that these are often the product of new genes that have evolved recently in the human lineage. Similarly, putative regulatory sequences are often not conserved but this is dismissed by declaring that regulatory sequences evolve quickly so they must be functional. Similar excuses are used to defend ubiquitous alternative splicing.

Most of these disputes could be settled by showing that these sequences are currently being maintained by purifying selection, or not. The problem, up until now, is that we don't have enough genome sequence data to resolve the controversy. This is about to change with the publications of hundreds of thousands of whole genome sequences (WGS) by the UK Biobank [Whole Genome Sequencing data on 200,000 UK Biobank participants made available for research].

The first substantial analaysis of these sequences was just published in Nature (July 28, 2022) and it contains a lot of information on "Functionally important regions."

Halldorsson, B.V., Eggertsson, H.P., Moore, K.H., Hauswedell, H., Eiriksson, O., Ulfarsson, M.O., Jensson, B.O. et al. (2022) The sequences of 150,119 genomes in the UK biobank. Nature 607:732-740. doi:[doi: 10.1038/s41586-022-04965-x]

Detailed knowledge of how diversity in the sequence of the human genome affects phenotypic diversity depends on a comprehensive and reliable characterization of both sequences and phenotypic variation. Over the past decade, insights into this relationship have been obtained from whole-exome sequencing or whole-genome sequencing of large cohorts with rich phenotypic data1,2. Here we describe the analysis of whole-genome sequencing of 150,119 individuals from the UK Biobank3. This constitutes a set of high-quality variants, including 585,040,410 single-nucleotide polymorphisms, representing 7.0% of all possible human single-nucleotide polymorphisms, and 58,707,036 indels. This large set of variants allows us to characterize selection based on sequence variation within a population through a depletion rank score of windows along the genome. Depletion rank analysis shows that coding exons represent a small fraction of regions in the genome subject to strong sequence conservation. We define three cohorts within the UK Biobank: a large British Irish cohort, a smaller African cohort and a South Asian cohort. A haplotype reference panel is provided that allows reliable imputation of most variants carried by three or more sequenced individuals. We identified 895,055 structural variants and 2,536,688 microsatellites, groups of variants typically excluded from large-scale whole-genome sequencing studies. Using this formidable new resource, we provide several examples of trait associations for rare variants with large effects not found previously through studies based on whole-exome sequencing and/or imputation.

The group sequenced the genomes of 150,119 individuals from the United Kingdom to an average coverage of 32.5X. That's good enough to reliably detect single nucleotide polymorphisms (SNPs) and small deletions and insertions (indels). (Keep in mind that these alleles will almost always be heterozygous so they will only be present in half of the DNA being sequenced.) Some of the genome sequences were unreliable and some participants withdrew consent so the final database had genome sequences from 149,960 individuals.

A total of 585,040,410 SNPs and 58,707,036 short indels were detected in the database. The authors state that the frequency of SNPs is one in 4.8 bp and I assume they arrive at that number by dividing the total number of SNPs into a genome size of 2.8 Gb (2.8/0.0585 = 4.8). Since there are 14.4 (3 x 4.8) possible variants, this represents 7% of all possible SNPs (1/14.4 = 7%).

These variants (alleles) will not be distributed randomly in the genome. Highly conserved regions (e.g. coding sequence) will have far fewer mutations than junk DNA. Other functional regions should also show a deficiency of variants. The paper contains an extensive description of how to go about identifying function by looking at these differences. What the authors are doing, in effect, is looking to see which parts of the genome are under purifying selection.

They do this by examining sequence variants in 500 bp windows and calculating an average expected number of variants based on the base composition and the mutation rate. They then calculate a depletion rank (DR) based on the number of variants that are missing in the windows. The rank ranges from zero (most depletion) to 100 (least depletion). These DR scores are then plotted against the locations of the windows relative to known protein-coding genes. The authors believe that this is a good way of identifying function.

Overall, our results show that DR can be used to help identify genomic regions under constraint across the entire genome and as such provides a valuable resource for identifying non-coding sequence of functional importance.

I'm not convinced that this is the best way to look at the effects of purifying selection but I certainly appreciate the effort.

The chart is very difficult to interpret because there's a large degree of stochastic noise in the data; for example, some regions of junk DNA will, by chance, have fewer mutations than expected (low DP score). Nevertheless, it's obvious that protein-coding genes (orange and red) are over-represented on the left of the chart because they have low DP scores. In other words, they have far fewer variants than expected so they must be subject to purifying selection.

It's more useful to quote the general conclusions reached by the authors.

As expected, coding exons have a low DR (mean DR = 28.4), but a large number of non-coding regions show even lower DR (more depletion), including non-coding regulatory elements. Among the 1% of regions with the lowest DR, 13.0% are coding and 87.0% are non-coding, with an overrepresentation of splice, UTR, gene upstream and downstream regions.

This is pretty much in line with what we expect if 90% of our genome is junk. Assuming that coding regions represent about 1% of the genome then the ratio of that 13% fraction to the other 87% corresponds to about 7% of the rest of the genome. We expect 10% of the genome to be functional and under purifying selection so this split seems reasonable.

We also expect that a significant part of that noncoding DNA will be within the protein-coding genes (splice sites and UTRs) or close by (regulatory sequences).

DR increases with distance from coding exons. After removing coding exons, among the 1% of regions with lowest and highest DR score, we saw a 3.2-fold and 0.4-fold overrepresentation of GWAS variants, respectively (Table 2), suggesting that the DR score could be a useful prior in GWAS analysis. ENCODE candidate cis-regulatory elements are more likely than expected by chance to be found in depleted (low DR) regions. Of note, candidate cis-regulatory elements located in close proximity to transcription start sites, that is, proximal enhancer-like and promoter-like sequences, are more enriched among depleted regions than distal enhancer-like sequences.

This is the most significant part of the analysis. The results indicate that functional regulatory sequences are located close to known genes. This is in line with the analyses of individual genes, which show that regulatory sequences are usually found within a 1000 bp of the transcription start site. It conflicts with the views of ENCODE researchers who claim that up to 20% of the genome consists of regulatory sequences scattered throughout the genome.

I don't think the paper looks at known, well studied, noncoding genes such as ribosomal RNA genes, tRNA genes, snoRNA genes etc. I assume that some of these sequences are included in the "intergenic" fraction on the chart and I assume that their associated regulatory sequences aren't assessed. I think that's because you can't write the appropriate code to look at functional noncoding genes since the annotation of the current reference sequence is unreliable.

Nevertheless, this is a good first step and it looks like a promising way to distinguish functional DNA from junk.

---