Genome-wide functional screen of 3'UTR variants uncovers causal variants for human disease and evolution.

3' untranslated region (3'UTR) variants are strongly associated with human traits and diseases, yet few have been causally identified. We developed the massively parallel reporter assay for 3'UTRs (MPRAu) to sensitively assay 12,173 3'UTR variants. We applied MPRAu to six human cell lines, focusing on genetic variants associated with genome-wide association studies (GWAS) and human evolutionary adaptation. MPRAu expands our understanding of 3'UTR function, suggesting that simple sequences predominately explain 3'UTR regulatory activity. We adapt MPRAu to uncover diverse molecular mechanisms at base pair resolution, including an adenylate-uridylate (AU)-rich element of LEPR linked to potential metabolic evolutionary adaptations in East Asians. We nominate hundreds of 3'UTR causal variants with genetically fine-mapped phenotype associations. Using endogenous allelic replacements, we characterize one variant that disrupts a miRNA site regulating the viral defense gene TRIM14 and one that alters PILRB abundance, nominating a causal variant underlying transcriptional changes in age-related macular degeneration.

The impact of rare variation on gene expression across tissues.

Rare genetic variants are abundant in humans and are expected to contribute to individual disease risk. While genetic association studies have successfully identified common genetic variants associated with susceptibility, these studies are not practical for identifying rare variants. Efforts to distinguish pathogenic variants from benign rare variants have leveraged the genetic code to identify deleterious protein-coding alleles, but no analogous code exists for non-coding variants. Therefore, ascertaining which rare variants have phenotypic effects remains a major challenge. Rare non-coding variants have been associated with extreme gene expression in studies using single tissues, but their effects across tissues are unknown. Here we identify gene expression outliers, or individuals showing extreme expression levels for a particular gene, across 44 human tissues by using combined analyses of whole genomes and multi-tissue RNA-sequencing data from the Genotype-Tissue Expression (GTEx) project v6p release. We find that 58% of underexpression and 28% of overexpression outliers have nearby conserved rare variants compared to 8% of non-outliers. Additionally, we developed RIVER (RNA-informed variant effect on regulation), a Bayesian statistical model that incorporates expression data to predict a regulatory effect for rare variants with higher accuracy than models using genomic annotations alone. Overall, we demonstrate that rare variants contribute to large gene expression changes across tissues and provide an integrative method for interpretation of rare variants in individual genomes.

Allele-specific expression reveals interactions between genetic variation and environment.

Identifying interactions between genetics and the environment (GxE) remains challenging. We have developed EAGLE, a hierarchical Bayesian model for identifying GxE interactions based on associations between environmental variables and allele-specific expression. Combining whole-blood RNA-seq with extensive environmental annotations collected from 922 human individuals, we identified 35 GxE interactions, compared with only four using standard GxE interaction testing. EAGLE provides new opportunities for researchers to identify GxE interactions using functional genomic data.

An Efficient Multiple-Testing Adjustment for eQTL Studies that Accounts for Linkage Disequilibrium between Variants.

Methods for multiple-testing correction in local expression quantitative trait locus (cis-eQTL) studies are a trade-off between statistical power and computational efficiency. Bonferroni correction, though computationally trivial, is overly conservative and fails to account for linkage disequilibrium between variants. Permutation-based methods are more powerful, though computationally far more intensive. We present an alternative correction method called eigenMT, which runs over 500 times faster than permutations and has adjusted p values that closely approximate empirical ones. To achieve this speed while also maintaining the accuracy of permutation-based methods, we estimate the effective number of independent variants tested for association with a particular gene, termed Meff, by using the eigenvalue decomposition of the genotype correlation matrix. We employ a regularized estimator of the correlation matrix to ensure Meff is robust and yields adjusted p values that closely approximate p values from permutations. Finally, using a common genotype matrix, we show that eigenMT can be applied with even greater efficiency to studies across tissues or conditions. Our method provides a simpler, more efficient approach to multiple-testing correction than existing methods and fits within existing pipelines for eQTL discovery.

Impact of the X Chromosome and sex on regulatory variation.

The X Chromosome, with its unique mode of inheritance, contributes to differences between the sexes at a molecular level, including sex-specific gene expression and sex-specific impact of genetic variation. Improving our understanding of these differences offers to elucidate the molecular mechanisms underlying sex-specific traits and diseases. However, to date, most studies have either ignored the X Chromosome or had insufficient power to test for the sex-specific impact of genetic variation. By analyzing whole blood transcriptomes of 922 individuals, we have conducted the first large-scale, genome-wide analysis of the impact of both sex and genetic variation on patterns of gene expression, including comparison between the X Chromosome and autosomes. We identified a depletion of expression quantitative trait loci (eQTL) on the X Chromosome, especially among genes under high selective constraint. In contrast, we discovered an enrichment of sex-specific regulatory variants on the X Chromosome. To resolve the molecular mechanisms underlying such effects, we generated chromatin accessibility data through ATAC-sequencing to connect sex-specific chromatin accessibility to sex-specific patterns of expression and regulatory variation. As sex-specific regulatory variants discovered in our study can inform sex differences in heritable disease prevalence, we integrated our data with genome-wide association study data for multiple immune traits identifying several traits with significant sex biases in genetic susceptibilities. Together, our study provides genome-wide insight into how genetic variation, the X Chromosome, and sex shape human gene regulation and disease.

FIRE: functional inference of genetic variants that regulate gene expression.

MOTIVATION: Interpreting genetic variation in noncoding regions of the genome is an important challenge for personal genome analysis. One mechanism by which noncoding single nucleotide variants (SNVs) influence downstream phenotypes is through the regulation of gene expression. Methods to predict whether or not individual SNVs are likely to regulate gene expression would aid interpretation of variants of unknown significance identified in whole-genome sequencing studies. RESULTS: We developed FIRE (Functional Inference of Regulators of Expression), a tool to score both noncoding and coding SNVs based on their potential to regulate the expression levels of nearby genes. FIRE consists of 23 random forests trained to recognize SNVs in cis-expression quantitative trait loci (cis-eQTLs) using a set of 92 genomic annotations as predictive features. FIRE scores discriminate cis-eQTL SNVs from non-eQTL SNVs in the training set with a cross-validated area under the receiver operating characteristic curve (AUC) of 0.807, and discriminate cis-eQTL SNVs shared across six populations of different ancestry from non-eQTL SNVs with an AUC of 0.939. FIRE scores are also predictive of cis-eQTL SNVs across a variety of tissue types. AVAILABILITY AND IMPLEMENTATION: FIRE scores for genome-wide SNVs in hg19/GRCh37 are available for download at https://sites.google.com/site/fireregulatoryvariation/. CONTACT: nilah@stanford.edu. SUPPLEMENTARY INFORMATION: Supplementary data are available at Bioinformatics online.

Incorporation of Biological Knowledge Into the Study of Gene-Environment Interactions.

A growing knowledge base of genetic and environmental information has greatly enabled the study of disease risk factors. However, the computational complexity and statistical burden of testing all variants by all environments has required novel study designs and hypothesis-driven approaches. We discuss how incorporating biological knowledge from model organisms, functional genomics, and integrative approaches can empower the discovery of novel gene-environment interactions and discuss specific methodological considerations with each approach. We consider specific examples where the application of these approaches has uncovered effects of gene-environment interactions relevant to drug response and immunity, and we highlight how such improvements enable a greater understanding of the pathogenesis of disease and the realization of precision medicine.

The origin, evolution, and functional impact of short insertion-deletion variants identified in 179 human genomes.

Short insertions and deletions (indels) are the second most abundant form of human genetic variation, but our understanding of their origins and functional effects lags behind that of other types of variants. Using population-scale sequencing, we have identified a high-quality set of 1.6 million indels from 179 individuals representing three diverse human populations. We show that rates of indel mutagenesis are highly heterogeneous, with 43%-48% of indels occurring in 4.03% of the genome, whereas in the remaining 96% their prevalence is 16 times lower than SNPs. Polymerase slippage can explain upwards of three-fourths of all indels, with the remainder being mostly simple deletions in complex sequence. However, insertions do occur and are significantly associated with pseudo-palindromic sequence features compatible with the fork stalling and template switching (FoSTeS) mechanism more commonly associated with large structural variations. We introduce a quantitative model of polymerase slippage, which enables us to identify indel-hypermutagenic protein-coding genes, some of which are associated with recurrent mutations leading to disease. Accounting for mutational rate heterogeneity due to sequence context, we find that indels across functional sequence are generally subject to stronger purifying selection than SNPs. We find that indel length modulates selection strength, and that indels affecting multiple functionally constrained nucleotides undergo stronger purifying selection. We further find that indels are enriched in associations with gene expression and find evidence for a contribution of nonsense-mediated decay. Finally, we show that indels can be integrated in existing genome-wide association studies (GWAS); although we do not find direct evidence that potentially causal protein-coding indels are enriched with associations to known disease-associated SNPs, our findings suggest that the causal variant underlying some of these associations may be indels.

Phenotypic instability between the near isogenic substrains BALB/cJ and BALB/cByJ.

Closely related substrains of inbred mice often show phenotypic differences that are presumed to be caused by recent mutations. The substrains BALB/cJ and BALB/cByJ, which were separated in 1935, have been reported to show numerous highly significant behavioral and morphological differences. In an effort to identify some of the causal mutations, we phenotyped BALB/cJ and BALB/cByJ mice as well as their F1, F2, and N2 progeny for behavioral and morphological phenotypes. We also generated whole-genome sequence data for both inbred strains (~3.5× coverage) with the intention of identifying polymorphic markers to be used for linkage analysis. We observed significant differences in body weight, the weight of the heart, liver, spleen and brain, and corpus callosum length between the two substrains. We also observed that BALB/cJ animals showed greater anxiety-like behavior in the open field test, less depression-like behavior in the tail suspension test, and reduced aggression compared to BALB/cByJ mice. Some but not all of these physiological and behavioral results were inconsistent with prior publications. These inconsistencies led us to suspect that the differences were due to, or modified by, non-genetic factors. Thus, we did not perform linkage analysis. We provide a comprehensive summary of the prior literature about phenotypic differences between these substrains as well as our current findings. We conclude that many differences between these strains are unstable and therefore ill-suited to linkage analysis; the source of this instability is unclear. We discuss the broader implications of these observations for the design of future studies.

Spatiotemporal patterns of Aedes aegypti populations in Cairns, Australia: assessing drivers of dengue transmission.

OBJECTIVES: To identify the meteorological drivers of dengue vector density and determine high- and low-risk transmission zones for dengue prevention and control in Cairns, Australia. METHODS: Weekly adult female Ae. aegypti data were obtained from 79 double sticky ovitraps (SOs) located in Cairns for the period September 2007-May 2012. Maximum temperature, total rainfall and average relative humidity data were obtained from the Australian Bureau of Meteorology for the study period. Time series-distributed lag nonlinear models were used to assess the relationship between meteorological variables and vector density. Spatial autocorrelation was assessed via semivariography, and ordinary kriging was undertaken to predict vector density in Cairns. RESULTS: Ae. aegypti density was associated with temperature and rainfall. However, these relationships differed between short (0-6 weeks) and long (0-30 weeks) lag periods. Semivariograms showed that vector distributions were spatially autocorrelated in September 2007-May 2008 and January 2009-May 2009, and vector density maps identified high transmission zones in the most populated parts of Cairns city, as well as Machans Beach. CONCLUSION: Spatiotemporal patterns of Ae. aegypti in Cairns are complex, showing spatial autocorrelation and associations with temperature and rainfall. Sticky ovitraps should be placed no more than 1.2 km apart to ensure entomological coverage and efficient use of resources. Vector density maps provide evidence for the targeting of prevention and control activities. Further research is needed to explore the possibility of developing an early warning system of dengue based on meteorological and environmental factors.

Peripubertal vitamin D(3) deficiency delays puberty and disrupts the estrous cycle in adult female mice.

The mechanism(s) by which vitamin D(3) regulates female reproduction is minimally understood. We tested the hypothesis that peripubertal vitamin D(3) deficiency disrupts hypothalamic-pituitary-ovarian physiology. To test this hypothesis, we used wild-type mice and Cyp27b1 (the rate-limiting enzyme in the synthesis of 1,25-dihydroxyvitamin D(3)) null mice to study the effect of vitamin D(3) deficiency on puberty and reproductive physiology. At the time of weaning, mice were randomized to a vitamin D(3)-replete or -deficient diet supplemented with calcium. We assessed the age of vaginal opening and first estrus (puberty markers), gonadotropin levels, ovarian histology, ovarian responsiveness to exogenous gonadotropins, and estrous cyclicity. Peripubertal vitamin D(3) deficiency significantly delayed vaginal opening without affecting the number of GnRH-immunopositive neurons or estradiol-negative feedback on gonadotropin levels during diestrus. Young adult females maintained on a vitamin D(3)-deficient diet after puberty had arrested follicular development and prolonged estrous cycles characterized by extended periods of diestrus. Ovaries of vitamin D(3)-deficient Cyp27b1 null mice responded to exogenous gonadotropins and deposited significantly more oocytes into the oviducts than mice maintained on a vitamin D(3)-replete diet. Estrous cycles were restored when vitamin D(3)-deficient Cyp27b1 null young adult females were transferred to a vitamin D(3)-replete diet. This study is the first to demonstrate that peripubertal vitamin D(3) sufficiency is important for an appropriately timed pubertal transition and maintenance of normal female reproductive physiology. These data suggest vitamin D(3) is a key regulator of neuroendocrine and ovarian physiology.

Prediction by statistical overlap theory of fraction of baseline occupied by chromatographic peaks.

The average minimum resolution required for separating adjacent single-component peaks (SCPs) in one-dimensional chromatograms is an important metric in statistical overlap theory (SOT). However, its value changes with changing chromatographic conditions in non-intuitive ways, when SOT predicts the average number of peaks (maxima). A more stable and easily understood value of resolution is obtained on making a different prediction. A general equation is derived for the sum of all separated and superposed widths of SCPs in a chromatogram. The equation is a function of the saturation α, a metric of chromatographic crowdedness, and is expressed in dimensionless form by dividing by the duration of the chromatogram. This dimensionless function, f(α), is also the cumulative distribution function of the probability of separating adjacent SCPs. Simulations based on the clustering of line segments representing SCPs verify expressions for f(α) calculated from five functions for the distribution of intervals between adjacent SCPs. Synthetic chromatograms are computed with different saturations, distributions of intervals, and distribution of SCP amplitudes. The chromatograms are analyzed by calculating the sum of the widths of peaks at different relative responses, dividing the sum by the duration of the chromatograms, and graphing the reduced sum against relative response. For small values of relative response, the reduced sum approaches the fraction of baseline that is occupied by chromatographic peaks. This fraction can be identified with f(α), if the saturation α is defined with the average minimum resolution equaling 1.5. The identification is general and is independent of the saturation, the interval distribution, or the amplitude distribution. This constant value of resolution corresponds to baseline resolution, which simplifies the interpretation of SOT.

Estimation of migration-time and mobility distributions in organelle capillary electrophoresis with statistical-overlap theory.

The separation of organelles by capillary electrophoresis (CE) produces large numbers of narrow peaks, which commonly are assumed to originate from single particles. In this paper, we show this is not always true. Here, we use established methods to partition simulated and real organelle CEs into regions of constant peak density and then use statistical-overlap theory to calculate the number of peaks (single particles) in each region. The only required measurements are the number of observed peaks (maxima) and peak standard deviation in the regions and the durations of the regions. Theory is developed for the precision of the estimated peak number and the threshold saturation above which the calculation is not advisable due to fluctuation of peak numbers. Theory shows that the relative precision is good when the saturation lies between 0.2 and 1.0 and is optimal when the saturation is slightly greater than 0.5. It also shows the threshold saturation depends on the peak standard deviation divided by the region's duration. The accuracy and precision of peak numbers estimated in different regions of organelle CEs are verified by computer simulations having both constant and nonuniform peak densities. The estimates are accurate to 6%. The estimated peak numbers in different regions are used to calculate migration-time and electrophoretic-mobility distributions. These distributions are less biased by peak overlap than ones determined by counting maxima and provide more correct measures of the organelle properties. The procedure is applied to a mitochondrial CE, in which over 20% of peaks are hidden by peak overlap.

Estimation of low-level components lost through chromatographic separations with finite detection limits.

The search for biomarkers allowing the assessment of disease by early diagnosis is facilitated by liquid chromatography. However, it is not clear how many components are lost due to being present in concentrations below the detection limit and/or being obscured by chromatographic peak overlap. First, we extend the study of missing components undertaken by Enke and Nagels, who employed the log-normal probability density function (pdf) for the distribution of signal intensities (and concentrations) of three mixtures. The Weibull and exponential pdfs, which have a higher probability of small-concentration components than the log-normal pdf, are also investigated. Results show that assessments of the loss of low-intensity signals by curve fitting are ambiguous. Next, we simulate synthetic chromatograms to compare the loss of peaks from superposition (overlap) with neighboring peaks to the loss arising from lying below the limit of detection (LOD) imposed by a finite signal-to-noise ratio (SNR). The simulations are made using amplitude pdfs based on the Enke-Nagels data as functions of relative column efficiency, i.e., saturation, and SNR. Results show that at the highest efficiencies, the lowest-amplitude peaks are lost below the LOD. However, at small and medium efficiencies, peak overlap is the dominant loss mechanism, suggesting that low-level components will not be found easily in liquid chromatography with single channel detectors regardless of SNR. A simple treatment shows that a multichannel detector, e.g., a mass spectrometer, is necessary to expose more low-level components.

Transcriptomic signatures across human tissues identify functional rare genetic variation.

Rare genetic variants are abundant across the human genome, and identifying their function and phenotypic impact is a major challenge. Measuring aberrant gene expression has aided in identifying functional, large-effect rare variants (RVs). Here, we expanded detection of genetically driven transcriptome abnormalities by analyzing gene expression, allele-specific expression, and alternative splicing from multitissue RNA-sequencing data, and demonstrate that each signal informs unique classes of RVs. We developed Watershed, a probabilistic model that integrates multiple genomic and transcriptomic signals to predict variant function, validated these predictions in additional cohorts and through experimental assays, and used them to assess RVs in the UK Biobank, the Million Veterans Program, and the Jackson Heart Study. Our results link thousands of RVs to diverse molecular effects and provide evidence to associate RVs affecting the transcriptome with human traits.

Effective saturation: a more informative metric for comparing peak separation in one- and two-dimensional separations.

A theoretical comparison is made of the numbers of observed peaks in one-dimensional (1D) and two-dimensional (2D) separations having the same peak capacity, as calculated from the traditional metric of resolution. The shortcoming of the average minimum resolution of statistical overlap theory (SOT) for this comparison is described. A new metric called the "effective saturation" is introduced to ameliorate the shortcoming. Unlike the "saturation", which is the usual metric of peak crowding in SOT, the effective saturation is independent of the average minimum resolution and can be determined using traditional values of resolution and peak capacity. Our most important finding is that, under a wide range of practical conditions, 1D and 2D separations of the same mixture produce almost equal numbers of observed peaks when the traditional peak capacities of the separations are the same, provided that the effective saturation and not the usual saturation is used as the measure of crowding. This is the case when peak distributions are random and when edge effects are minor. The numerical results supporting this finding can be described by empirical functions of the effective saturation, including one for the traditional peak capacity needed to separate a given fraction of mixture constituents as observed peaks. The near equality of the number of observed peaks in 1D and 2D separations based on the effective saturation is confirmed by simulations. However, this equality is compromised in 2D separations when edge effects are large. The new finding does not contradict previous predictions by SOT of differences between 1D and 2D separations at equal saturation. Indeed, the simulations reaffirm their validity. Rather, the usual metric, i.e., the saturation, is just not as simple a metric for comparing 1D and 2D separations as is the new metric, i.e., the effective saturation. We strongly recommend use of the new metric for its great simplifying effect.

Theory of the probability of total resolution in chromatograms with systematic variation of average peak spacing and peak width.

An equation is proposed for the probability that all mixture constituents are separated, when the density (i.e., average number of eluting constituents per time) and width of single-component peaks (SCPs) vary systematically. The probability Pr that m SCPs are separated is modeled as the product of the m - 1 probabilities that adjacent pairs of SCPs are separated. Pr is then expressed as the geometric mean of the probability product raised to the power of m - 1. This geometric mean is approximated by an arithmetic mean equaling the probability that adjacent SCPs are separated, as calculated from previously developed statistical overlap theory (SOT) for variable SCP density and width. The theory is tested using previously reported and current in-house simulations of isocratic chromatograms of SCPs with random differences in standard chemical potential. In such chromatograms, more SCPs elute at short times than long times, and their widths are less at short times than long times. The average difference between 179 previously reported and currently predicted values of 100 x Pr is about 0.6, when 100 x Pr > 5. The theory requires numerical computation of one integral but can be approximated by an analytic equation for SOT probabilities close to one. For SCPs having retention times exceeding twice the void time, this equation simplifies to a previous SOT expression, with the gradient peak capacity replaced by the isocratic peak capacity. The versatility of the Pr theory is tested using three other models of chromatograms, in which the density and width of SCPs vary. The Pr predictions agree with simulation for all three models.

Dominance of the tiger: The displacement of Aedes aegypti by Aedes albopictus in parts of the Torres Strait, Australia

Most of the inhabited islands in the Torres Strait region of Australia have experienced dengue outbreaks transmitted by Aedes aegypti at various times since at least the 1890s. However, another potential dengue vector, Aedes albopictus, the Asian tiger mosquito, was detected for the first time in 2005 and it expanded across most of the Torres Strait within a few years. In 2016, a survey of container-inhabiting mosquitoes was conducted in all island communities and Ae. aegypti was undetectable on most of the islands which the species had previously occupied, and had been replaced by Ae. albopictus. It is suspected that competitive displacement was responsible for the changes in species distribution. Aedes aegypti was only detected on Boigu Island and Thursday Island. Recent dengue outbreaks in the Torres Strait have apparently been driven by both Ae. albopictus and Ae. aegypti. The findings have major implications on management of dengue outbreaks in the region.

Identification of rare-disease genes using blood transcriptome sequencing and large control cohorts.

It is estimated that 350 million individuals worldwide suffer from rare diseases, which are predominantly caused by mutation in a single gene1. The current molecular diagnostic rate is estimated at 50%, with whole-exome sequencing (WES) among the most successful approaches2-5. For patients in whom WES is uninformative, RNA sequencing (RNA-seq) has shown diagnostic utility in specific tissues and diseases6-8. This includes muscle biopsies from patients with undiagnosed rare muscle disorders6,9, and cultured fibroblasts from patients with mitochondrial disorders7. However, for many individuals, biopsies are not performed for clinical care, and tissues are difficult to access. We sought to assess the utility of RNA-seq from blood as a diagnostic tool for rare diseases of different pathophysiologies. We generated whole-blood RNA-seq from 94 individuals with undiagnosed rare diseases spanning 16 diverse disease categories. We developed a robust approach to compare data from these individuals with large sets of RNA-seq data for controls (n = 1,594 unrelated controls and n = 49 family members) and demonstrated the impacts of expression, splicing, gene and variant filtering strategies on disease gene identification. Across our cohort, we observed that RNA-seq yields a 7.5% diagnostic rate, and an additional 16.7% with improved candidate gene resolution.