ClusTrast: a short read de novo transcript isoform assembler guided by clustered contigs

BACKGROUND

Transcriptome assembly from RNA-sequencing data in species without a reliable reference genome has to be performed de novo, but studies have shown that de novo methods often have inadequate ability to reconstruct transcript isoforms. We address this issue by constructing an assembly pipeline whose main purpose is to produce a comprehensive set of transcript isoforms.

RESULTS

We present the de novo transcript isoform assembler ClusTrast, which takes short read RNA-seq data as input, assembles a primary assembly, clusters a set of guiding contigs, aligns the short reads to the guiding contigs, assembles each clustered set of short reads individually, and merges the primary and clusterwise assemblies into the final assembly. We tested ClusTrast on real datasets from six eukaryotic species, and showed that ClusTrast reconstructed more expressed known isoforms than any of the other tested de novo assemblers, at a moderate reduction in precision. For recall, ClusTrast was on top in the lower end of expression levels (<15% percentile) for all tested datasets, and over the entire range for almost all datasets. Reference transcripts were often (35-69% for the six datasets) reconstructed to at least 95% of their length by ClusTrast, and more than half of reference transcripts (58-81%) were reconstructed with contigs that exhibited polymorphism, measuring on a subset of reliably predicted contigs. ClusTrast recall increased when using a union of assembled transcripts from more than one assembly tool as primary assembly.

CONCLUSION

We suggest that ClusTrast can be a useful tool for studying isoforms in species without a reliable reference genome, in particular when the goal is to produce a comprehensive transcriptome set with polymorphic variants.

Erythroid Differentiation Enhances RNA Mis-Splicing in SF3B1-Mutant Myelodysplastic Syndromes with Ring Sideroblasts

Myelodysplastic syndromes with ring sideroblasts (MDS-RS) commonly develop from hematopoietic stem cells (HSC) bearing mutations in the splicing factor SF3B1 (SF3B1mt). Direct studies into MDS-RS pathobiology have been limited by a lack of model systems that fully recapitulate erythroid biology and RS development and the inability to isolate viable human RS. Here, we combined successful direct RS isolation from patient samples, high-throughput multiomics analysis of cells encompassing the SF3B1mt stem-erythroid continuum, and functional assays to investigate the impact of SF3B1mt on erythropoiesis and RS accumulation. The isolated RS differentiated, egressed into the blood, escaped traditional nonsense-mediated decay (NMD) mechanisms, and leveraged stress-survival pathways that hinder wild-type hematopoiesis through pathogenic GDF15 overexpression. Importantly, RS constituted a contaminant of magnetically enriched CD34+ cells, skewing bulk transcriptomic data. Mis-splicing in SF3B1mt cells was intensified by erythroid differentiation through accelerated RNA splicing and decreased NMD activity, and SF3B1mt led to truncations in several MDS-implicated genes. Finally, RNA mis-splicing induced an uncoupling of RNA and protein expression, leading to critical abnormalities in proapoptotic p53 pathway genes. Overall, this characterization of erythropoiesis in SF3B1mt RS provides a resource for studying MDS-RS and uncovers insights into the unexpectedly active biology of the "dead-end" RS.

SIGNIFICANCE

Ring sideroblast isolation combined with state-of-the-art multiomics identifies survival mechanisms underlying SF3B1-mutant erythropoiesis and establishes an active role for erythroid differentiation and ring sideroblasts themselves in SF3B1-mutant myelodysplastic syndrome pathogenesis.

Single-cell multi-omics identifies chronic inflammation as a driver of TP53-mutant leukemic evolution

Understanding the genetic and nongenetic determinants of tumor protein 53 (TP53)-mutation-driven clonal evolution and subsequent transformation is a crucial step toward the design of rational therapeutic strategies. Here we carry out allelic resolution single-cell multi-omic analysis of hematopoietic stem/progenitor cells (HSPCs) from patients with a myeloproliferative neoplasm who transform to TP53-mutant secondary acute myeloid leukemia (sAML). All patients showed dominant TP53 'multihit' HSPC clones at transformation, with a leukemia stem cell transcriptional signature strongly predictive of adverse outcomes in independent cohorts, across both TP53-mutant and wild-type (WT) AML. Through analysis of serial samples, antecedent TP53-heterozygous clones and in vivo perturbations, we demonstrate a hitherto unrecognized effect of chronic inflammation, which suppressed TP53 WT HSPCs while enhancing the fitness advantage of TP53-mutant cells and promoted genetic evolution. Our findings will facilitate the development of risk-stratification, early detection and treatment strategies for TP53-mutant leukemia, and are of broad relevance to other cancer types.

Genotyping and population characteristics of the China Kadoorie Biobank

The China Kadoorie Biobank (CKB) is a population-based prospective cohort of >512,000 adults recruited from 2004 to 2008 from 10 geographically diverse regions across China. Detailed data from questionnaires and physical measurements were collected at baseline, with additional measurements at three resurveys involving ∼5% of surviving participants. Analyses of genome-wide genotyping, for >100,000 participants using custom-designed Axiom arrays, reveal extensive relatedness, recent consanguinity, and signatures reflecting large-scale population movements from recent Chinese history. Systematic genome-wide association studies of incident disease, captured through electronic linkage to death and disease registries and to the national health insurance system, replicate established disease loci and identify 14 novel disease associations. Together with studies of candidate drug targets and disease risk factors and contributions to international genetics consortia, these demonstrate the breadth, depth, and quality of the CKB data. Ongoing high-throughput omics assays of collected biosamples and planned whole-genome sequencing will further enhance the scientific value of this biobank.

Cone-setting in spruce is regulated by conserved elements of the age-dependent flowering pathway

Reproductive phase change is well characterized in angiosperm model species, but less studied in gymnosperms. We utilize the early cone-setting acrocona mutant to study reproductive phase change in the conifer Picea abies (Norway spruce), a gymnosperm. The acrocona mutant frequently initiates cone-like structures, called transition shoots, in positions where wild-type P. abies always produces vegetative shoots. We collect acrocona and wild-type samples, and RNA-sequence their messenger RNA (mRNA) and microRNA (miRNA) fractions. We establish gene expression patterns and then use allele-specific transcript assembly to identify mutations in acrocona. We genotype a segregating population of inbred acrocona trees. A member of the SQUAMOSA BINDING PROTEIN-LIKE (SPL) gene family, PaSPL1, is active in reproductive meristems, whereas two putative negative regulators of PaSPL1, miRNA156 and the conifer specific miRNA529, are upregulated in vegetative and transition shoot meristems. We identify a mutation in a putative miRNA156/529 binding site of the acrocona PaSPL1 allele and show that the mutation renders the acrocona allele tolerant to these miRNAs. We show co-segregation between the early cone-setting phenotype and trees homozygous for the acrocona mutation. In conclusion, we demonstrate evolutionary conservation of the age-dependent flowering pathway and involvement of this pathway in regulating reproductive phase change in the conifer P. abies.

The MDM2 Inhibitor Navtemadlin Arrests Mouse Melanoma Growth In Vivo and Potentiates Radiotherapy

The tumor suppressor protein p53 is mutated in close to 50% of human tumors and is dysregulated in many others, for instance by silencing or loss of p14^ARF. Under steady-state conditions, the two E3 ligases MDM2/MDM4 interact with and inhibit the transcriptional activity of p53. Inhibition of p53-MDM2/4 interaction to reactivate p53 in tumors with wild-type (WT) p53 has therefore been considered a therapeutic strategy. Moreover, studies indicate that p53 reactivation may synergize with radiation and increase tumor immunogenicity. In vivo studies of most MDM2 inhibitors have utilized immunodeficient xenograft mouse models, preventing detailed studies of action of these molecules on the immune response. The mouse melanoma cell line B16-F10 carries functional, WT p53 but does not express the MDM2 regulator p19^ARF. In this study, we tested a p53-MDM2 protein-protein interaction inhibitor, the small molecule Navtemadlin, which is currently being tested in phase II clinical trials. Using mass spectrometry-based proteomics and imaging flow cytometry, we identified specific protein expression patterns following Navtemadlin treatment of B16-F10 melanoma cells compared with their p53 CRISPR-inactivated control cells. In vitro, Navtemadlin induced a significant, p53-dependent, growth arrest but little apoptosis in B16-F10 cells. When combined with radiotherapy, Navtemadlin showed synergistic effects and increased apoptosis. In vivo, Navtemadlin treatment significantly reduced the growth of B16-F10 melanoma cells implanted in C57Bl/6 mice. Our data highlight the utility of a syngeneic B16-F10 p53^+/+ mouse melanoma model for assessing existing and novel p53-MDM2/MDM4 inhibitors and in identifying new combination therapies that can efficiently eliminate tumors in vivo.

Significance

The MDM2 inhibitor Navtemadlin arrests mouse tumor growth and potentiates radiotherapy. Our results support a threshold model for apoptosis induction that requires a high, prolonged p53 signaling for cancer cells to become apoptotic.

Identifying the Common Genetic Basis of Antidepressant Response

Background

Antidepressants are a first-line treatment for depression. However, only a third of individuals experience remission after the first treatment. Common genetic variation, in part, likely regulates antidepressant response, yet the success of previous genome-wide association studies has been limited by sample size. This study performs the largest genetic analysis of prospectively assessed antidepressant response in major depressive disorder to gain insight into the underlying biology and enable out-of-sample prediction.

Methods

Genome-wide analysis of remission (n remit = 1852, n nonremit = 3299) and percentage improvement (n = 5218) was performed. Single nucleotide polymorphism-based heritability was estimated using genome-wide complex trait analysis. Genetic covariance with eight mental health phenotypes was estimated using polygenic scores/AVENGEME. Out-of-sample prediction of antidepressant response polygenic scores was assessed. Gene-level association analysis was performed using MAGMA and transcriptome-wide association study. Tissue, pathway, and drug binding enrichment were estimated using MAGMA.

Results

Neither genome-wide association study identified genome-wide significant associations. Single nucleotide polymorphism-based heritability was significantly different from zero for remission (h 2 = 0.132, SE = 0.056) but not for percentage improvement (h 2 = -0.018, SE = 0.032). Better antidepressant response was negatively associated with genetic risk for schizophrenia and positively associated with genetic propensity for educational attainment. Leave-one-out validation of antidepressant response polygenic scores demonstrated significant evidence of out-of-sample prediction, though results varied in external cohorts. Gene-based analyses identified ETV4 and DHX8 as significantly associated with antidepressant response.

Conclusions

This study demonstrates that antidepressant response is influenced by common genetic variation, has a genetic overlap schizophrenia and educational attainment, and provides a useful resource for future research. Larger sample sizes are required to attain the potential of genetics for understanding and predicting antidepressant response.

Efficient ancestry and mutation simulation with msprime 1.0

Stochastic simulation is a key tool in population genetics, since the models involved are often analytically intractable and simulation is usually the only way of obtaining ground-truth data to evaluate inferences. Because of this, a large number of specialized simulation programs have been developed, each filling a particular niche, but with largely overlapping functionality and a substantial duplication of effort. Here, we introduce msprime version 1.0, which efficiently implements ancestry and mutation simulations based on the succinct tree sequence data structure and the tskit library. We summarize msprime’s many features, and show that its performance is excellent, often many times faster and more memory efficient than specialized alternatives. These high-performance features have been thoroughly tested and validated, and built using a collaborative, open source development model, which reduces duplication of effort and promotes software quality via community engagement.

Reciprocal interactions between tumour cell populations enhance growth and reduce radiation sensitivity in prostate cancer

Intratumoural heterogeneity (ITH) contributes to local recurrence following radiotherapy in prostate cancer. Recent studies also show that ecological interactions between heterogeneous tumour cell populations can lead to resistance in chemotherapy. Here, we evaluated whether interactions between heterogenous populations could impact growth and response to radiotherapy in prostate cancer. Using mixed 3D cultures of parental and radioresistant populations from two prostate cancer cell lines and a predator-prey mathematical model to investigate various types of ecological interactions, we show that reciprocal interactions between heterogeneous populations enhance overall growth and reduce radiation sensitivity. The type of interaction influences the time of regrowth after radiation, and, at the population level, alters the survival and cell cycle of each population without eliminating either one. These interactions can arise from oxygen constraints and from cellular cross-talk that alter the tumour microenvironment. These findings suggest that ecological-type interactions are important in radiation response and could be targeted to reduce local recurrence.

Functional Parameters Derived from Magnetic Resonance Imaging Reflect Vascular Morphology in Preclinical Tumors and in Human Liver Metastases

Purpose: Tumor vessels influence the growth and response of tumors to therapy. Imaging vascular changes in vivo using dynamic contrast-enhanced MRI (DCE-MRI) has shown potential to guide clinical decision making for treatment. However, quantitative MR imaging biomarkers of vascular function have not been widely adopted, partly because their relationship to structural changes in vessels remains unclear. We aimed to elucidate the relationships between vessel function and morphology in vivo Experimental Design: Untreated preclinical tumors with different levels of vascularization were imaged sequentially using DCE-MRI and CT. Relationships between functional parameters from MR (iAUC, K trans, and BATfrac) and structural parameters from CT (vessel volume, radius, and tortuosity) were assessed using linear models. Tumors treated with anti-VEGFR2 antibody were then imaged to determine whether antiangiogenic therapy altered these relationships. Finally, functional-structural relationships were measured in 10 patients with liver metastases from colorectal cancer.Results: Functional parameters iAUC and K trans primarily reflected vessel volume in untreated preclinical tumors. The relationships varied spatially and with tumor vascularity, and were altered by antiangiogenic treatment. In human liver metastases, all three structural parameters were linearly correlated with iAUC and K trans For iAUC, structural parameters also modified each other's effect.Conclusions: Our findings suggest that MR imaging biomarkers of vascular function are linked to structural changes in tumor vessels and that antiangiogenic therapy can affect this link. Our work also demonstrates the feasibility of three-dimensional functional-structural validation of MR biomarkers in vivo to improve their biological interpretation and clinical utility. Clin Cancer Res; 24(19); 4694-704. ©2018 AACR.

Genome-wide association analyses identify 44 risk variants and refine the genetic architecture of major depression

Major depressive disorder (MDD) is a common illness accompanied by considerable morbidity, mortality, costs, and heightened risk of suicide. We conducted a genome-wide association meta-analysis based in 135,458 cases and 344,901 controls and identified 44 independent and significant loci. The genetic findings were associated with clinical features of major depression and implicated brain regions exhibiting anatomical differences in cases. Targets of antidepressant medications and genes involved in gene splicing were enriched for smaller association signal. We found important relationships of genetic risk for major depression with educational attainment, body mass, and schizophrenia: lower educational attainment and higher body mass were putatively causal, whereas major depression and schizophrenia reflected a partly shared biological etiology. All humans carry lesser or greater numbers of genetic risk factors for major depression. These findings help refine the basis of major depression and imply that a continuous measure of risk underlies the clinical phenotype.

Integrative Analysis of Three RNA Sequencing Methods Identifies Mutually Exclusive Exons of MADS-Box Isoforms During Early Bud Development in Picea abies

Recent efforts to sequence the genomes and transcriptomes of several gymnosperm species have revealed an increased complexity in certain gene families in gymnosperms as compared to angiosperms. One example of this is the gymnosperm sister clade to angiosperm TM3-like MADS-box genes, which at least in the conifer lineage has expanded in number of genes. We have previously identified a member of this sub-clade, the conifer gene DEFICIENS AGAMOUS LIKE 19 (DAL19), as being specifically upregulated in cone-setting shoots. Here, we show through Sanger sequencing of mRNA-derived cDNA and mapping to assembled conifer genomic sequences that DAL19 produces six mature mRNA splice variants in Picea abies. These splice variants use alternate first and last exons, while their four central exons constitute a core region present in all six transcripts. Thus, they are likely to be transcript isoforms. Quantitative Real-Time PCR revealed that two mutually exclusive first DAL19 exons are differentially expressed across meristems that will form either male or female cones, or vegetative shoots. Furthermore, mRNA in situ hybridization revealed that two mutually exclusive last DAL19 exons were expressed in a cell-specific pattern within bud meristems. Based on these findings in DAL19, we developed a sensitive approach to transcript isoform assembly from short-read sequencing of mRNA. We applied this method to 42 putative MADS-box core regions in P. abies, from which we assembled 1084 putative transcripts. We manually curated these transcripts to arrive at 933 assembled transcript isoforms of 38 putative MADS-box genes. 152 of these isoforms, which we assign to 28 putative MADS-box genes, were differentially expressed across eight female, male, and vegetative buds. We further provide evidence of the expression of 16 out of the 38 putative MADS-box genes by mapping PacBio Iso-Seq circular consensus reads derived from pooled sample sequencing to assembled transcripts. In summary, our analyses reveal the use of mutually exclusive exons of MADS-box gene isoforms during early bud development in P. abies, and we find that the large number of identified MADS-box transcripts in P. abies results not only from expansion of the gene family through gene duplication events but also from the generation of numerous splice variants.

An automated method measures variability in P-glycoprotein and ABCG2 densities across brain regions and brain matter

Changes in P-glycoprotein and ABCG2 densities may play a role in amyloid-beta accumulation in Alzheimer's disease. However, previous studies report conflicting results from different brain regions, without correcting for changes in vessel density. We developed an automated method to measure transporter density exclusively within the vascular space, thereby correcting for vessel density. We then examined variability in transporter density across brain regions, matter, and disease using two cohorts of post-mortem brains from Alzheimer's disease patients and age-matched controls. Changes in transporter density were also investigated in capillaries near plaques and on the mRNA level. P-glycoprotein density varied with brain region and matter, whereas ABCG2 density varied with brain matter. In temporal cortex, P-glycoprotein density was 53% lower in Alzheimer's disease samples than in controls, and was reduced by 35% in capillaries near plaque deposits within Alzheimer's disease samples. ABCG2 density was unaffected in Alzheimer's disease. No differences were detected at the transcript level. Our study indicates that region-specific changes in transporter densities can occur globally and locally near amyloid-beta deposits in Alzheimer's disease, providing an explanation for conflicting results in the literature. When differences in region and matter are accounted for, changes in density can be reproducibly measured using our automated method.

Amyloid arthropathy associated with multiple myeloma: polyarthritis without synovial infiltration of CD20+ or CD38+ cells

OBJECTIVES

To describe histological, immunohistochemical and ultrastructural features of synovial biopsies of amyloid arthropathy associated with multiple myeloma (MM).

METHODS

Synovial biopsies from affected joints of two patients with MM and amyloid arthropathy were examined with light and electron microscopy, and immunohistochemically for expression of CD3, CD8, CD20, CD38, CD68, Ki-67 and vWF. Results were compared to values from osteoarthritis (OA, n = 26), rheumatoid arthritis (RA, n = 24) and normal (n = 15) synovial membranes.

RESULTS

There was no or only mild lining hyperplasia. Vascular density was not elevated, and there were few Ki-67+ proliferating cells in the stroma. The Krenn synovitis score classified one specimen as "low-grade" and one as "high-grade" synovitis. CD68+ and CD3+ cells were the predominant mononuclear inflammatory cells, whereas CD20+ and CD38+ cells were absent from both synovial membrane and synovial fluid sediment. Electron microscopy demonstrated amyloid phagocytosis by synovial macrophages. In hierarchical clustering the two amyloid arthropathy specimens were more closely related to OA than to RA or normal synovium.

CONCLUSIONS

This first detailed immunohistological analysis of MM-associated amyloid arthropathy suggests that it is a chronic synovitis that evolves despite the loss of humoral immunity seen in advanced MM. Instead, amyloid phagocytosis by synovial macrophages likely triggers and perpetuates local disease.

Sequence diversity of Pan troglodytes subspecies and the impact of WFDC6 selective constraints in reproductive immunity

Recent efforts have attempted to describe the population structure of common chimpanzee, focusing on four subspecies: Pan troglodytes verus, P. t. ellioti, P. t. troglodytes, and P. t. schweinfurthii. However, few studies have pursued the effects of natural selection in shaping their response to pathogens and reproduction. Whey acidic protein (WAP) four-disulfide core domain (WFDC) genes and neighboring semenogelin (SEMG) genes encode proteins with combined roles in immunity and fertility. They display a strikingly high rate of amino acid replacement (d_N/d_S), indicative of adaptive pressures during primate evolution. In human populations, three signals of selection at the WFDC locus were described, possibly influencing the proteolytic profile and antimicrobial activities of the male reproductive tract. To evaluate the patterns of genomic variation and selection at the WFDC locus in chimpanzees, we sequenced 17 WFDC genes and 47 autosomal pseudogenes in 68 chimpanzees (15 P. t. troglodytes, 22 P. t. verus, and 31 P. t. ellioti). We found a clear differentiation of P. t. verus and estimated the divergence of P. t. troglodytes and P. t. ellioti subspecies in 0.173 Myr; further, at the WFDC locus we identified a signature of strong selective constraints common to the three subspecies in WFDC6—a recent paralog of the epididymal protease inhibitor EPPIN. Overall, chimpanzees and humans do not display similar footprints of selection across the WFDC locus, possibly due to different selective pressures between the two species related to immune response and reproductive biology.

Reproduction and immunity-driven natural selection in the human WFDC locus

The whey acidic protein (WAP) four-disulfide core domain (WFDC) locus located on human chromosome 20q13 spans 19 genes with WAP and/or Kunitz domains. These genes participate in antimicrobial, immune, and tissue homoeostasis activities. Neighboring SEMG genes encode seminal proteins Semenogelin 1 and 2 (SEMG1 and SEMG2). WFDC and SEMG genes have a strikingly high rate of amino acid replacement (dN/dS), indicative of responses to adaptive pressures during vertebrate evolution. To better understand the selection pressures acting on WFDC genes in human populations, we resequenced 18 genes and 54 noncoding segments in 71 European (CEU), African (YRI), and Asian (CHB + JPT) individuals. Overall, we identified 484 single-nucleotide polymorphisms (SNPs), including 65 coding variants (of which 49 are nonsynonymous differences). Using classic neutrality tests, we confirmed the signature of short-term balancing selection on WFDC8 in Europeans and a signature of positive selection spanning genes PI3, SEMG1, SEMG2, and SLPI. Associated with the latter signal, we identified an unusually homogeneous-derived 100-kb haplotype with a frequency of 88% in Asian populations. A putative candidate variant targeted by selection is Thr56Ser in SEMG1, which may alter the proteolytic profile of SEMG1 and antimicrobial activities of semen. All the well-characterized genes residing in the WDFC locus encode proteins that appear to have a role in immunity and/or fertility, two processes that are often associated with adaptive evolution. This study provides further evidence that the WFDC and SEMG loci have been under strong adaptive pressure within the short timescale of modern humans.

Easy detection of GFP-positive mutants following forward mutations at specific gene locus in cultured human cells

We have generated a new mutation assay system using HT1080 human fibrosarcoma cells, which consists of a combination of tetracycline-operator dependent GFP gene (TetO-EGFP) and tetracycline repressor (TetR) genes, where the expression of GFP gene is under strict control of TetR protein, and the TetR gene is located within the endogenous HPRT gene. In this system, any inactivating mutation at the TetR gene or large deletions including the gene itself results in high expression of GFP gene (>200-fold increase) in the cells, which can be readily scored not only by a flow cytometer but also under a fluorescent microscope. With this new cell line, we show that the spontaneous mutation rate at the TetR locus was 2.8-3.4×10(-6)/cell division, slightly lower than the rate at the endogenous HPRT gene of HT1080 cells, and has a dose response to X rays as a mutagen. We also isolated variant clones with elevated spontaneous mutation rate (i.e., genetically unstable cells) following X irradiation. Spontaneous GFP-positive mutants were predominantly base-change mutations at the TetR gene while those obtained after X irradiation often contained large deletions which spanned up to 6Mb. The results indicate that the bacterial TetR/TetO regulatory units work extremely well as a mutation detection system in human cells, and any part of the human genome may be tested for mutation sensitivity following targeted insertion of the TetR gene in a stably expressing gene.

Balancing selection maintains a form of ERAP2 that undergoes nonsense-mediated decay and affects antigen presentation

A remarkable characteristic of the human major histocompatibility complex (MHC) is its extreme genetic diversity, which is maintained by balancing selection. In fact, the MHC complex remains one of the best-known examples of natural selection in humans, with well-established genetic signatures and biological mechanisms for the action of selection. Here, we present genetic and functional evidence that another gene with a fundamental role in MHC class I presentation, endoplasmic reticulum aminopeptidase 2 (ERAP2), has also evolved under balancing selection and contains a variant that affects antigen presentation. Specifically, genetic analyses of six human populations revealed strong and consistent signatures of balancing selection affecting ERAP2. This selection maintains two highly differentiated haplotypes (Haplotype A and Haplotype B), with frequencies 0.44 and 0.56, respectively. We found that ERAP2 expressed from Haplotype B undergoes differential splicing and encodes a truncated protein, leading to nonsense-mediated decay of the mRNA. To investigate the consequences of ERAP2 deficiency on MHC presentation, we correlated surface MHC class I expression with ERAP2 genotypes in primary lymphocytes. Haplotype B homozygotes had lower levels of MHC class I expressed on the surface of B cells, suggesting that naturally occurring ERAP2 deficiency affects MHC presentation and immune response. Interestingly, an ERAP2 paralog, endoplasmic reticulum aminopeptidase 1 (ERAP1), also shows genetic signatures of balancing selection. Together, our findings link the genetic signatures of selection with an effect on splicing and a cellular phenotype. Although the precise selective pressure that maintains polymorphism is unknown, the demonstrated differences between the ERAP2 splice forms provide important insights into the potential mechanism for the action of selection.

It has long been known that the extremely high levels of genetic diversity present in the major histocompatibility locus (MHC) are due to balancing selection, a type of natural selection that maintains advantageous genetic diversity in populations. The MHC encodes for molecules required for a type of antigen presentation that mediates detection of infected and cancerous cells by the immune system; the genetic diversity of the MHC thus ensures an adequate response to the wide variety of pathogens that humans encounter. Here, we show that other genes involved in the same antigen-presentation pathway are also subject to balancing selection in humans. Specifically, we show that balancing selection acts to maintain two forms of the endoplasmic reticulum aminopeptidase 2 gene (ERAP2), which encodes a protein also involved in antigen presentation. Although the two ERAP2 forms are present in a similar frequency (close to 0.5), they are associated with differences with respect to the levels of MHC molecules on the cell surface of immune cells. In summary, our findings show that natural selection maintains variants of ERAP2 that affect immune surveillance; they also establish ERAP2 as one of the few examples of balancing selection in humans where the selected variant, its functional consequences, and its influence in interpersonal diversity are known.