Complexities of recapitulating polygenic effects in natural populations: replication of genetic effects on wing shape in artificially selected and wild-caught populations of Drosophila melanogaster.

Identifying the genetic architecture of complex traits is important to many geneticists, including those interested in human disease, plant and animal breeding, and evolutionary genetics. Advances in sequencing technology and statistical methods for genome-wide association studies have allowed for the identification of more variants with smaller effect sizes, however, many of these identified polymorphisms fail to be replicated in subsequent studies. In addition to sampling variation, this failure to replicate reflects the complexities introduced by factors including environmental variation, genetic background, and differences in allele frequencies among populations. Using Drosophila melanogaster wing shape, we ask if we can replicate allelic effects of polymorphisms first identified in a genome-wide association studies in three genes: dachsous, extra-macrochaete, and neuralized, using artificial selection in the lab, and bulk segregant mapping in natural populations. We demonstrate that multivariate wing shape changes associated with these genes are aligned with major axes of phenotypic and genetic variation in natural populations. Following seven generations of artificial selection along the dachsous shape change vector, we observe genetic differentiation of variants in dachsous and genomic regions containing other genes in the hippo signaling pathway. This suggests a shared direction of effects within a developmental network. We also performed artificial selection with the extra-macrochaete shape change vector, which is not a part of the hippo signaling network, but showed a largely shared direction of effects. The response to selection along the emc vector was similar to that of dachsous, suggesting that the available genetic diversity of a population, summarized by the genetic (co)variance matrix (G), influenced alleles captured by selection. Despite the success with artificial selection, bulk segregant analysis using natural populations did not detect these same variants, likely due to the contribution of environmental variation and low minor allele frequencies, coupled with small effect sizes of the contributing variants.

An experimental model of anti-PD-1 resistance exhibits activation of TGFß and Notch pathways and is sensitive to local mRNA immunotherapy.

Immune checkpoint blockade elicits durable anti-cancer responses in the clinic, however a large proportion of patients do not benefit from treatment. Several mechanisms of innate and acquired resistance to checkpoint blockade have been defined and include mutations of MHC I and IFNγ signaling pathways. However, such mutations occur in a low frequency of patients and additional mechanisms have yet to be elucidated. In an effort to better understand acquired resistance to checkpoint blockade, we generated a mouse tumor model exhibiting in vivo resistance to anti-PD-1 antibody treatment. MC38 tumors acquired resistance to PD-1 blockade following serial in vivo passaging. Lack of sensitivity to PD-1 blockade was not attributed to dysregulation of PD-L1 or ß2M expression, as both were expressed at similar levels in parental and resistant cells. Similarly, IFNγ signaling and antigen processing and presentation pathways were functional in both parental and resistant cell lines. Unbiased gene expression analysis was used to further characterize potential resistance mechanisms. RNA-sequencing revealed substantial differences in global gene expression, with tumors resistant to anti-PD-1 displaying a marked reduction in expression of immune-related genes relative to parental MC38 tumors. Indeed, resistant tumors exhibited reduced immune infiltration across multiple cell types, including T and NK cells. Pathway analysis revealed activation of TGFß and Notch signaling in anti-PD-1 resistant tumors, and activation of these pathways was associated with poorer survival in human cancer patients. While pharmacological inhibition of TGFß and Notch in combination with PD-1 blockade decelerated tumor growth, a local mRNA-based immunotherapy potently induced regression of resistant tumors, resulting in complete tumor remission, and resensitized tumors to treatment with anti-PD-1. Overall, this study describes a novel anti-PD-1 resistant mouse tumor model and underscores the role of two well-defined signaling pathways in response to immune checkpoint blockade. Furthermore, our data highlights the potential of intratumoral mRNA therapy in overcoming acquired resistance to PD-1 blockade.

The lethal sex gap: COVID-19.

While Coronavirus disease 2019 (COVID-19), caused by severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2), is disrupting lives across the globe for everyone, it has a more devastating impact on the health of older adults, especially that of older men. This pandemic has highlighted the crucial importance of considering an individual's age and biological sex in the clinic in addition to other confounding diseases (Kuchel, G.A, J Am Geriatr Soc, 67, 203, 2019, Tannenbaum, C., Nature, 575 451-458, 2009) As an interdisciplinary team of scientists in immunology, hematology, genomics, bioinformatics, and geriatrics, we have been studying how age and sex shape the human immune system. Herein we reflect on how our recent findings on the alterations of the immune system in aging might contribute to our current understanding of COVID-19 infection rate and disease risk.

V-SVA: an R Shiny application for detecting and annotating hidden sources of variation in single-cell RNA-seq data.

SUMMARY: Single-cell RNA-sequencing (scRNA-seq) technology enables studying gene expression programs from individual cells. However, these data are subject to diverse sources of variation, including 'unwanted' variation that needs to be removed in downstream analyses (e.g. batch effects) and 'wanted' or biological sources of variation (e.g. variation associated with a cell type) that needs to be precisely described. Surrogate variable analysis (SVA)-based algorithms, are commonly used for batch correction and more recently for studying 'wanted' variation in scRNA-seq data. However, interpreting whether these variables are biologically meaningful or stemming from technical reasons remains a challenge. To facilitate the interpretation of surrogate variables detected by algorithms including IA-SVA, SVA or ZINB-WaVE, we developed an R Shiny application [Visual Surrogate Variable Analysis (V-SVA)] that provides a web-browser interface for the identification and annotation of hidden sources of variation in scRNA-seq data. This interactive framework includes tools for discovery of genes associated with detected sources of variation, gene annotation using publicly available databases and gene sets, and data visualization using dimension reduction methods. AVAILABILITY AND IMPLEMENTATION: The V-SVA Shiny application is publicly hosted at https://vsva.jax.org/ and the source code is freely available at https://github.com/nlawlor/V-SVA. CONTACT: leed13@miamioh.edu or duygu.ucar@jax.org. SUPPLEMENTARY INFORMATION: Supplementary data are available at Bioinformatics online.

Sexual-dimorphism in human immune system aging.

Differences in immune function and responses contribute to health- and life-span disparities between sexes. However, the role of sex in immune system aging is not well understood. Here, we characterize peripheral blood mononuclear cells from 172 healthy adults 22-93 years of age using ATAC-seq, RNA-seq, and flow cytometry. These data reveal a shared epigenomic signature of aging including declining naïve T cell and increasing monocyte and cytotoxic cell functions. These changes are greater in magnitude in men and accompanied by a male-specific decline in B-cell specific loci. Age-related epigenomic changes first spike around late-thirties with similar timing and magnitude between sexes, whereas the second spike is earlier and stronger in men. Unexpectedly, genomic differences between sexes increase after age 65, with men having higher innate and pro-inflammatory activity and lower adaptive activity. Impact of age and sex on immune phenotypes can be visualized at https://immune-aging.jax.org to provide insights into future studies.

A Multivariate Genome-Wide Association Study of Wing Shape in Drosophila melanogaster.

Due to the complexity of genotype-phenotype relationships, simultaneous analyses of genomic associations with multiple traits will be more powerful and informative than a series of univariate analyses. However, in most cases, studies of genotype-phenotype relationships have been analyzed only one trait at a time. Here, we report the results of a fully integrated multivariate genome-wide association analysis of the shape of the Drosophila melanogaster wing in the Drosophila Genetic Reference Panel. Genotypic effects on wing shape were highly correlated between two different laboratories. We found 2396 significant SNPs using a 5% false discovery rate cutoff in the multivariate analyses, but just four significant SNPs in univariate analyses of scores on the first 20 principal component axes. One quarter of these initially significant SNPs retain their effects in regularized models that take into account population structure and linkage disequilibrium. A key advantage of multivariate analysis is that the direction of the estimated phenotypic effect is much more informative than a univariate one. We exploit this fact to show that the effects of knockdowns of genes implicated in the initial screen were on average more similar than expected under a null model. A subset of SNP effects were replicable in an unrelated panel of inbred lines. Association studies that take a phenomic approach, considering many traits simultaneously, are an important complement to the power of genomics.

Multiomic Profiling Identifies cis-Regulatory Networks Underlying Human Pancreatic ß Cell Identity and Function.

EndoC-ßH1 is emerging as a critical human ß cell model to study the genetic and environmental etiologies of ß cell (dys)function and diabetes. Comprehensive knowledge of its molecular landscape is lacking, yet required, for effective use of this model. Here, we report chromosomal (spectral karyotyping), genetic (genotyping), epigenomic (ChIP-seq and ATAC-seq), chromatin interaction (Hi-C and Pol2 ChIA-PET), and transcriptomic (RNA-seq and miRNA-seq) maps of EndoC-ßH1. Analyses of these maps define known (e.g., PDX1 and ISL1) and putative (e.g., PCSK1 and mir-375) ß cell-specific transcriptional cis-regulatory networks and identify allelic effects on cis-regulatory element use. Importantly, comparison with maps generated in primary human islets and/or ß cells indicates preservation of chromatin looping but also highlights chromosomal aberrations and fetal genomic signatures in EndoC-ßH1. Together, these maps, and a web application we created for their exploration, provide important tools for the design of experiments to probe and manipulate the genetic programs governing ß cell identity and (dys)function in diabetes.

BiFET: sequencing Bias-free transcription factor Footprint Enrichment Test.

Transcription factor (TF) footprinting uncovers putative protein-DNA binding via combined analyses of chromatin accessibility patterns and their underlying TF sequence motifs. TF footprints are frequently used to identify TFs that regulate activities of cell/condition-specific genomic regions (target loci) in comparison to control regions (background loci) using standard enrichment tests. However, there is a strong association between the chromatin accessibility level and the GC content of a locus and the number and types of TF footprints that can be detected at this site. Traditional enrichment tests (e.g. hypergeometric) do not account for this bias and inflate false positive associations. Therefore, we developed a novel post-processing method, Bias-free Footprint Enrichment Test (BiFET), that corrects for the biases arising from the differences in chromatin accessibility levels and GC contents between target and background loci in footprint enrichment analyses. We applied BiFET on TF footprint calls obtained from EndoC-ßH1 ATAC-seq samples using three different algorithms (CENTIPEDE, HINT-BC and PIQ) and showed BiFET's ability to increase power and reduce false positive rate when compared to hypergeometric test. Furthermore, we used BiFET to study TF footprints from human PBMC and pancreatic islet ATAC-seq samples to show its utility to identify putative TFs associated with cell-type-specific loci.

Type 2 Diabetes-Associated Genetic Variants Regulate Chromatin Accessibility in Human Islets.

Type 2 diabetes (T2D) is a complex disorder in which both genetic and environmental risk factors contribute to islet dysfunction and failure. Genome-wide association studies (GWAS) have linked single nucleotide polymorphisms (SNPs), most of which are noncoding, in >200 loci to islet dysfunction and T2D. Identification of the putative causal variants and their target genes and whether they lead to gain or loss of function remains challenging. Here, we profiled chromatin accessibility in pancreatic islet samples from 19 genotyped individuals and identified 2,949 SNPs associated with in vivo cis-regulatory element use (i.e., chromatin accessibility quantitative trait loci [caQTL]). Among the caQTLs tested (n = 13) using luciferase reporter assays in MIN6 ß-cells, more than half exhibited effects on enhancer activity that were consistent with in vivo chromatin accessibility changes. Importantly, islet caQTL analysis nominated putative causal SNPs in 13 T2D-associated GWAS loci, linking 7 and 6 T2D risk alleles, respectively, to gain or loss of in vivo chromatin accessibility. By investigating the effect of genetic variants on chromatin accessibility in islets, this study is an important step forward in translating T2D-associated GWAS SNP into functional molecular consequences.

Chromatin interaction networks revealed unique connectivity patterns of broad H3K4me3 domains and super enhancers in 3D chromatin.

Broad domain promoters and super enhancers are regulatory elements that govern cell-specific functions and harbor disease-associated sequence variants. These elements are characterized by distinct epigenomic profiles, such as expanded deposition of histone marks H3K27ac for super enhancers and H3K4me3 for broad domains, however little is known about how they interact with each other and the rest of the genome in three-dimensional chromatin space. Using network theory methods, we studied chromatin interactions between broad domains and super enhancers in three ENCODE cell lines (K562, MCF7, GM12878) obtained via ChIA-PET, Hi-C, and Hi-CHIP assays. In these networks, broad domains and super enhancers interact more frequently with each other compared to their typical counterparts. Network measures and graphlets revealed distinct connectivity patterns associated with these regulatory elements that are robust across cell types and alternative assays. Machine learning models showed that these connectivity patterns could effectively discriminate broad domains from typical promoters and super enhancers from typical enhancers. Finally, targets of broad domains in these networks were enriched in disease-causing SNPs of cognate cell types. Taken together these results suggest a robust and unique organization of the chromatin around broad domains and super enhancers: loci critical for pathologies and cell-specific functions.

The chromatin accessibility signature of human immune aging stems from CD8+ T cells.

Aging is linked to deficiencies in immune responses and increased systemic inflammation. To unravel the regulatory programs behind these changes, we applied systems immunology approaches and profiled chromatin accessibility and the transcriptome in PBMCs and purified monocytes, B cells, and T cells. Analysis of samples from 77 young and elderly donors revealed a novel and robust aging signature in PBMCs, with simultaneous systematic chromatin closing at promoters and enhancers associated with T cell signaling and a potentially stochastic chromatin opening mostly found at quiescent and repressed sites. Combined analyses of chromatin accessibility and the transcriptome uncovered immune molecules activated/inactivated with aging and identified the silencing of the IL7R gene and the IL-7 signaling pathway genes as potential biomarkers. This signature is borne by memory CD8+ T cells, which exhibited an aging-related loss in binding of NF-κB and STAT factors. Thus, our study provides a unique and comprehensive approach to identifying candidate biomarkers and provides mechanistic insights into aging-associated immunodeficiency.

Leukaemia cell of origin identified by chromatin landscape of bulk tumour cells.

The precise identity of a tumour's cell of origin can influence disease prognosis and outcome. Methods to reliably define tumour cell of origin from primary, bulk tumour cell samples has been a challenge. Here we use a well-defined model of MLL-rearranged acute myeloid leukaemia (AML) to demonstrate that transforming haematopoietic stem cells (HSCs) and multipotent progenitors results in more aggressive AML than transforming committed progenitor cells. Transcriptome profiling reveals a gene expression signature broadly distinguishing stem cell-derived versus progenitor cell-derived AML, including genes involved in immune escape, extravasation and small GTPase signal transduction. However, whole-genome profiling of open chromatin reveals precise and robust biomarkers reflecting each cell of origin tested, from bulk AML tumour cell sampling. We find that bulk AML tumour cells exhibit distinct open chromatin loci that reflect the transformed cell of origin and suggest that open chromatin patterns may be leveraged as prognostic signatures in human AML.

QuIN: A Web Server for Querying and Visualizing Chromatin Interaction Networks.

UNLABELLED: Recent studies of the human genome have indicated that regulatory elements (e.g. promoters and enhancers) at distal genomic locations can interact with each other via chromatin folding and affect gene expression levels. Genomic technologies for mapping interactions between DNA regions, e.g., ChIA-PET and HiC, can generate genome-wide maps of interactions between regulatory elements. These interaction datasets are important resources to infer distal gene targets of non-coding regulatory elements and to facilitate prioritization of critical loci for important cellular functions. With the increasing diversity and complexity of genomic information and public ontologies, making sense of these datasets demands integrative and easy-to-use software tools. Moreover, network representation of chromatin interaction maps enables effective data visualization, integration, and mining. Currently, there is no software that can take full advantage of network theory approaches for the analysis of chromatin interaction datasets. To fill this gap, we developed a web-based application, QuIN, which enables: 1) building and visualizing chromatin interaction networks, 2) annotating networks with user-provided private and publicly available functional genomics and interaction datasets, 3) querying network components based on gene name or chromosome location, and 4) utilizing network based measures to identify and prioritize critical regulatory targets and their direct and indirect interactions. AVAILABILITY: QuIN's web server is available at http://quin.jax.org QuIN is developed in Java and JavaScript, utilizing an Apache Tomcat web server and MySQL database and the source code is available under the GPLV3 license available on GitHub: https://github.com/UcarLab/QuIN/.

The tandem duplicator phenotype as a distinct genomic configuration in cancer.

Next-generation sequencing studies have revealed genome-wide structural variation patterns in cancer, such as chromothripsis and chromoplexy, that do not engage a single discernable driver mutation, and whose clinical relevance is unclear. We devised a robust genomic metric able to identify cancers with a chromotype called tandem duplicator phenotype (TDP) characterized by frequent and distributed tandem duplications (TDs). Enriched only in triple-negative breast cancer (TNBC) and in ovarian, endometrial, and liver cancers, TDP tumors conjointly exhibit tumor protein p53 (TP53) mutations, disruption of breast cancer 1 (BRCA1), and increased expression of DNA replication genes pointing at rereplication in a defective checkpoint environment as a plausible causal mechanism. The resultant TDs in TDP augment global oncogene expression and disrupt tumor suppressor genes. Importantly, the TDP strongly correlates with cisplatin sensitivity in both TNBC cell lines and primary patient-derived xenografts. We conclude that the TDP is a common cancer chromotype that coordinately alters oncogene/tumor suppressor expression with potential as a marker for chemotherapeutic response.

Complex constraints on allometry revealed by artificial selection on the wing of Drosophila melanogaster.

Precise exponential scaling with size is a fundamental aspect of phenotypic variation. These allometric power laws are often invariant across taxa and have long been hypothesized to reflect developmental constraints. Here we test this hypothesis by investigating the evolutionary potential of an allometric scaling relationship in drosophilid wing shape that is nearly invariant across 111 species separated by at least 50 million years of evolution. In only 26 generations of artificial selection in a population of Drosophila melanogaster, we were able to drive the allometric slope to the outer range of those found among the 111 sampled species. This response was rapidly lost when selection was suspended. Only a small proportion of this reversal could be explained by breakup of linkage disequilibrium, and direct selection on wing shape is also unlikely to explain the reversal, because the more divergent wing shapes produced by selection on the allometric intercept did not revert. We hypothesize that the reversal was instead caused by internal selection arising from pleiotropic links to unknown traits. Our results also suggest that the observed selection response in the allometric slope was due to a component expressed late in larval development and that variation in earlier development did not respond to selection. Together, these results are consistent with a role for pleiotropic constraints in explaining the remarkable evolutionary stability of allometric scaling.

Linkage Disequilibrium and Inversion-Typing of the Drosophila melanogaster Genome Reference Panel.

We calculated the linkage disequilibrium between all pairs of variants in the Drosophila Genome Reference Panel with minor allele count ≥5. We used r(2) ≥ 0.5 as the cutoff for a highly correlated SNP. We make available the list of all highly correlated SNPs for use in association studies. Seventy-six percent of variant SNPs are highly correlated with at least one other SNP, and the mean number of highly correlated SNPs per variant over the whole genome is 83.9. Disequilibrium between distant SNPs is also common when minor allele frequency (MAF) is low: 37% of SNPs with MAF < 0.1 are highly correlated with SNPs more than 100 kb distant. Although SNPs within regions with polymorphic inversions are highly correlated with somewhat larger numbers of SNPs, and these correlated SNPs are on average farther away, the probability that a SNP in such regions is highly correlated with at least one other SNP is very similar to SNPs outside inversions. Previous karyotyping of the DGRP lines has been inconsistent, and we used LD and genotype to investigate these discrepancies. When previous studies agreed on inversion karyotype, our analysis was almost perfectly concordant with those assignments. In discordant cases, and for inversion heterozygotes, our results suggest errors in two previous analyses or discordance between genotype and karyotype. Heterozygosities of chromosome arms are, in many cases, surprisingly highly correlated, suggesting strong epsistatic selection during the inbreeding and maintenance of the DGRP lines.

The Measurement of Local Variation in Shape.

Geometric morphometrics comprises tools for measuring and analyzing shape as captured by an entire set of landmark configurations. Many interesting questions in evolutionary, genetic, and developmental research, however, are only meaningful at a local level, where a focus on "parts" or "traits" takes priority over properties of wholes. To study variational properties of such traits, current approaches partition configurations into subsets of landmarks which are then studied separately. This approach is unable to fully capture both variational and spatial characteristics of these subsets because interpretability of shape differences is context-dependent. Landmarks omitted from a partition usually contain information about that partition's shape. We present an interpolation-based approach that can be used to model shape differences at a local, infinitesimal level as a function of information available globally. This approach belongs in a large family of methods that see shape differences as continuous "fields" spanning an entire structure, for which landmarks serve as reference parameters rather than as data. We show, via analyses of simulated and real data, how interpolation models provide a more accurate representation of regional shapes than partitioned data. A key difference of this interpolation approach from current morphometric practice is that one must assume an explicit interpolation model, which in turn implies a particular kind of behavior of the regions between landmarks. This choice presents novel methodological challenges, but also an opportunity to incorporate and test biomechanical models that have sought to explain tissue-level processes underlying the generation of morphological shape.

Constraint and opportunity: the genetic basis and evolution of modularity in the cichlid mandible.

Modular variation, whereby the relative degree of connectivity varies within a system, is thought to evolve through a process of selection that favors the integration of certain traits and the decoupling of others. In this way, modularity may facilitate the pace of evolution and determine evolvability. Alternatively, conserved patterns of modularity may act to constrain the rate and direction of evolution by preventing certain functions from evolving. A comprehensive understanding of the potential interplay between these phenomena will require knowledge of the inheritance and the genetic basis of modularity. Here we explore these ideas in the cichlid mandible by investigating patterns of modularity at the clade and species levels and through the introduction of a new approach, the individual level. Specifically, we assessed patterns of covariation in Lake Malawi cichlid species that employ alternate "biting" and "suction-feeding" modes of feeding and in a hybrid cross between these two ecotypes. Across the suction-feeding clade, patterns of modularity were largely conserved and reflected a functionally based pattern. In contrast, the biting species displayed a pattern of modularity that more closely matched developmental modules. The pattern of modularity present in our F2 population was very similar to the pattern exhibited by the biter, suggesting a role for dominant inheritance. We demonstrate that our individual-level metric of modularity (IMM) is a valid quantitative trait that has a nonlinear relationship with shape. IMMs for each model were used as quantitative characters to map quantitative trait loci (QTL) that underlie modularity. Our QTL analysis offers new insights into the genetic basis of modularity in these fishes that may eventually lead to the discovery of the genetic processes that delineate particular modules. In all, our findings suggest that modularity is both a constraining and an evolvable force in cichlid evolution, as distinct patterns occur between species and variation exists among individuals.

A statistical framework for testing modularity in multidimensional data.

Modular variation of multivariate traits results from modular distribution of effects of genetic and epigenetic interactions among those traits. However, statistical methods rarely detect truly modular patterns, possibly because the processes that generate intramodular associations may overlap spatially. Methodologically, this overlap may cause multiple patterns of modularity to be equally consistent with observed covariances. To deal with this indeterminacy, the present study outlines a framework for testing a priori hypotheses of modularity in which putative modules are mathematically represented as multidimensional subspaces embedded in the data. Model expectations are computed by subdividing the data into arrays of variables, and intermodular interactions are represented by overlapping arrays. Covariance structures are thus modeled as the outcome of complex and nonorthogonal intermodular interactions. This approach is demonstrated by analyzing mandibular modularity in nine rodent species. A total of 620 models are fit to each species, and the most strongly supported are heuristically modified to improve their fit. Five modules common to all species are identified, which approximately map to the developmental modules of the mandible. Within species, these modules are embedded within larger "super-modules," suggesting that these conserved modules act as building blocks from which covariation patterns are built.