Benchmarking Tree and Ancestral Sequence Inference for B Cell Receptor Sequences

Inferring B Cell Phylogenies from Paired H and L Chain BCR Sequences with Dowser

Abs are vital to human immune responses and are composed of genetically variable H and L chains. These structures are initially expressed as BCRs. BCR diversity is shaped through somatic hypermutation and selection during immune responses. This evolutionary process produces B cell clones, cells that descend from a common ancestor but differ by mutations. Phylogenetic trees inferred from BCR sequences can reconstruct the history of mutations within a clone. Until recently, BCR sequencing technologies separated H and L chains, but advancements in single-cell sequencing now pair H and L chains from individual cells. However, it is unclear how these separate genes should be combined to infer B cell phylogenies. In this study, we investigated strategies for using paired H and L chain sequences to build phylogenetic trees. We found that incorporating L chains significantly improved tree accuracy and reproducibility across all methods tested. This improvement was greater than the difference between tree-building methods and persisted even when mixing bulk and single-cell sequencing data. However, we also found that many phylogenetic methods estimated significantly biased branch lengths when some L chains were missing, such as when mixing single-cell and bulk BCR data. This bias was eliminated using maximum likelihood methods with separate branch lengths for H and L chain gene partitions. Thus, we recommend using maximum likelihood methods with separate H and L chain partitions, especially when mixing data types. We implemented these methods in the R package Dowser: https://dowser.readthedocs.io.

B cell phylogenetics in the single cell era

The widespread availability of single-cell RNA sequencing (scRNA-seq) has led to the development of new methods for understanding immune responses. Single-cell transcriptome data can now be paired with B cell receptor (BCR) sequences. However, RNA from BCRs cannot be analyzed like most other genes because BCRs are genetically diverse within individuals. In humans, BCRs are shaped through recombination followed by mutation and selection for antigen binding. As these processes co-occur with cell division, B cells can be studied using phylogenetic trees representing the mutations within a clone. B cell trees can link experimental timepoints, tissues, or cellular subtypes. Here, we review the current state and potential of how B cell phylogenetics can be combined with single-cell data to understand immune responses.

TRIBAL: Tree Inference of B cell Clonal Lineages

B cells are a critical component of the adaptive immune system, responsible for producing antibodies that help protect the body from infections and foreign substances. Single cell RNA-sequencing (scRNA-seq) has allowed for both profiling of B cell receptor (BCR) sequences and gene expression. However, understanding the adaptive and evolutionary mechanisms of B cells in response to specific stimuli remains a significant challenge in the field of immunology. We introduce a new method, TRIBAL, which aims to infer the evolutionary history of clonally related B cells from scRNA-seq data. The key insight of TRIBAL is that inclusion of isotype data into the B cell lineage inference problem is valuable for reducing phylogenetic uncertainty that arises when only considering the receptor sequences. Consequently, the TRIBAL inferred B cell lineage trees jointly capture the somatic mutations introduced to the B cell receptor during affinity maturation and isotype transitions during class switch recombination. In addition, TRIBAL infers isotype transition probabilities that are valuable for gaining insight into the dynamics of class switching. Via in silico experiments, we demonstrate that TRIBAL infers isotype transition probabilities with the ability to distinguish between direct versus sequential switching in a B cell population. This results in more accurate B cell lineage trees and corresponding ancestral sequence and class switch reconstruction compared to competing methods. Using real-world scRNA-seq datasets, we show that TRIBAL recapitulates expected biological trends in a model affinity maturation system. Furthermore, the B cell lineage trees inferred by TRIBAL were equally plausible for the BCR sequences as those inferred by competing methods but yielded lower entropic partitions for the isotypes of the sequenced B cell. Thus, our method holds the potential to further advance our understanding of vaccine responses, disease progression, and the identification of therapeutic antibodies.

Inferring B cell phylogenies from paired heavy and light chain BCR sequences with Dowser

Antibodies are vital to human immune responses and are composed of genetically variable heavy and light chains. These structures are initially expressed as B cell receptors (BCRs). BCR diversity is shaped through somatic hypermutation and selection during immune responses. This evolutionary process produces B cell clones, cells that descend from a common ancestor but differ by mutations. Phylogenetic trees inferred from BCR sequences can reconstruct the history of mutations within a clone. Until recently, BCR sequencing technologies separated heavy and light chains, but advancements in single cell sequencing now pair heavy and light chains from individual cells. However, it is unclear how these separate genes should be combined to infer B cell phylogenies. In this study, we investigated strategies for using paired heavy and light chain sequences to build phylogenetic trees. We found incorporating light chains significantly improved tree accuracy and reproducibility across all methods tested. This improvement was greater than the difference between tree building methods and persisted even when mixing bulk and single cell sequencing data. However, we also found that many phylogenetic methods estimated significantly biased branch lengths when some light chains were missing, such as when mixing single cell and bulk BCR data. This bias was eliminated using maximum likelihood methods with separate branch lengths for heavy and light chain gene partitions. Thus, we recommend using maximum likelihood methods with separate heavy and light chain partitions, especially when mixing data types. We implemented these methods in the R package Dowser: https://dowser.readthedocs.io.

Reconstructing B cell lineage trees with minimum spanning tree and genotype abundances

B cell receptor (BCR) genes exposed to an antigen undergo somatic hypermutations and Darwinian antigen selection, generating a large BCR-antibody diversity. This process, known as B cell affinity maturation, increases antibody affinity, forming a specific B cell lineage that includes the unmutated ancestor and mutated variants. In a B cell lineage, cells with a higher antigen affinity will undergo clonal expansion, while those with a lower affinity will not proliferate and probably be eliminated. Therefore, cellular (genotype) abundance provides a valuable perspective on the ongoing evolutionary process. Phylogenetic tree inference is often used to reconstruct B cell lineage trees and represents the evolutionary dynamic of BCR affinity maturation. However, such methods should process B-cell population data derived from experimental sampling that might contain different cellular abundances. There are a few phylogenetic methods for tracing the evolutionary events occurring in B cell lineages; best-performing solutions are time-demanding and restricted to analysing a reduced number of sequences, while time-efficient methods do not consider cellular abundances. We propose ClonalTree, a low-complexity and accurate approach to construct B-cell lineage trees that incorporates genotype abundances into minimum spanning tree (MST) algorithms. Using both simulated and experimental data, we demonstrate that ClonalTree outperforms MST-based algorithms and achieves a comparable performance to a method that explores tree-generating space exhaustively. Furthermore, ClonalTree has a lower running time, being more convenient for building B-cell lineage trees from high-throughput BCR sequencing data, mainly in biomedical applications, where a lower computational time is appreciable. It is hundreds to thousands of times faster than exhaustive approaches, enabling the analysis of a large set of sequences within minutes or seconds and without loss of accuracy. The source code is freely available at github.com/julibinho/ClonalTree.

A scalable model for simulating multi-round antibody evolution and benchmarking of clonal tree reconstruction methods

Affinity maturation (AM) of B cells through somatic hypermutations (SHMs) enables the immune system to evolve to recognize diverse pathogens. The accumulation of SHMs leads to the formation of clonal lineages of antibody-secreting b cells that have evolved from a common naïve B cell. Advances in high-throughput sequencing have enabled deep scans of B cell receptor repertoires, paving the way for reconstructing clonal trees. However, it is not clear if clonal trees, which capture microevolutionary time scales, can be reconstructed using traditional phylogenetic reconstruction methods with adequate accuracy. In fact, several clonal tree reconstruction methods have been developed to fix supposed shortcomings of phylogenetic methods. Nevertheless, no consensus has been reached regarding the relative accuracy of these methods, partially because evaluation is challenging. Benchmarking the performance of existing methods and developing better methods would both benefit from realistic models of clonal lineage evolution specifically designed for emulating B cell evolution. In this paper, we propose a model for modeling B cell clonal lineage evolution and use this model to benchmark several existing clonal tree reconstruction methods. Our model, designed to be extensible, has several features: by evolving the clonal tree and sequences simultaneously, it allows modeling selective pressure due to changes in affinity binding; it enables scalable simulations of large numbers of cells; it enables several rounds of infection by an evolving pathogen; and, it models building of memory. In addition, we also suggest a set of metrics for comparing clonal trees and measuring their properties. Our results show that while maximum likelihood phylogenetic reconstruction methods can fail to capture key features of clonal tree expansion if applied naively, a simple post-processing of their results, where short branches are contracted, leads to inferences that are better than alternative methods.

Inference of B cell clonal families using heavy/light chain pairing information

Next generation sequencing of B cell receptor (BCR) repertoires has become a ubiquitous tool for understanding the antibody-mediated immune response: it is now common to have large volumes of sequence data coding for both the heavy and light chain subunits of the BCR. However, until the recent development of high throughput methods of preserving heavy/light chain pairing information, these samples contained no explicit information on which heavy chain sequence pairs with which light chain sequence. One of the first steps in analyzing such BCR repertoire samples is grouping sequences into clonally related families, where each stems from a single rearrangement event. Many methods of accomplishing this have been developed, however, none so far has taken full advantage of the newly-available pairing information. This information can dramatically improve clustering performance, especially for the light chain. The light chain has traditionally been challenging for clonal family inference because of its low diversity and consequent abundance of non-clonal families with indistinguishable naive rearrangements. Here we present a method of incorporating this pairing information into the clustering process in order to arrive at a more accurate partition of the data into clonally related families. We also demonstrate two methods of fixing imperfect pairing information, which may allow for simplified sample preparation and increased sequencing depth. Finally, we describe several other improvements to the partis software package.

Echidna: integrated simulations of single-cell immune receptor repertoires and transcriptomes

Motivation

Single-cell sequencing now enables the recovery of full-length immune receptor repertoires [B cell receptor (BCR) and T cell receptor (TCR) repertoires], in addition to gene expression information. The feature-rich datasets produced from such experiments require extensive and diverse computational analyses, each of which can significantly influence the downstream immunological interpretations, such as clonal selection and expansion. Simulations produce validated standard datasets, where the underlying generative model can be precisely defined and furthermore perturbed to investigate specific questions of interest. Currently, there is no tool that can be used to simulate single-cell datasets incorporating immune receptor repertoires and gene expression.

Results

We developed Echidna, an R package that simulates immune receptors and transcriptomes at single-cell resolution with user-tunable parameters controlling a wide range of features such as clonal expansion, germline gene usage, somatic hypermutation, transcriptional phenotypes and spatial location. Echidna can additionally simulate time-resolved B cell evolution, producing mutational networks with complex selection histories incorporating class-switching and B cell subtype information. We demonstrated the benchmarking potential of Echidna by simulating clonal lineages and comparing the known simulated networks with those inferred from only the BCR sequences as input. Finally, we simulated immune repertoire information onto existing spatial transcriptomic experiments, thereby generating novel datasets that could be used to develop and integrate methods to profile clonal selection in a spatially resolved manner. Together, Echidna provides a framework that can incorporate experimental data to simulate single-cell immune repertoires to aid software development and bioinformatic benchmarking of clonotyping, phylogenetics, transcriptomics and machine learning strategies.

Availability and implementation

The R package and code used in this manuscript can be found at github.com/alexyermanos/echidna and also in the R package Platypus (Yermanos et al., 2021). Installation instructions and the vignette for Echidna is described in the Platypus Computational Ecosystem (https://alexyermanos.github.io/Platypus/index.html). Publicly available data and corresponding sample accession numbers can be found in Supplementary Tables S2 and S3.

Supplementary information

Supplementary data are available at Bioinformatics Advances online.

Antibody repertoire sequencing analysis

High-throughput sequencing for B cell receptor (BCR) repertoire provides useful insights for the adaptive immune system. With the continuous development of the BCR-seq technology, many efforts have been made to develop methods for analyzing the ever-increasing BCR repertoire data. In this review, we comprehensively outline different BCR repertoire library preparation protocols and summarize three major steps of BCR-seq data analysis, i. e., V(D)J sequence annotation, clonal phylogenetic inference, and BCR repertoire profiling and mining. Different from other reviews in this field, we emphasize background intuition and the statistical principle of each method to help biologists better understand it. Finally, we discuss data mining problems for BCR-seq data and with a highlight on recently emerging multiple-sample analysis.

Phylogenetic analysis of migration, differentiation, and class switching in B cells

B cells undergo rapid mutation and selection for antibody binding affinity when producing antibodies capable of neutralizing pathogens. This evolutionary process can be intermixed with migration between tissues, differentiation between cellular subsets, and switching between functional isotypes. B cell receptor (BCR) sequence data has the potential to elucidate important information about these processes. However, there is currently no robust, generalizable framework for making such inferences from BCR sequence data. To address this, we develop three parsimony-based summary statistics to characterize migration, differentiation, and isotype switching along B cell phylogenetic trees. We use simulations to demonstrate the effectiveness of this approach. We then use this framework to infer patterns of cellular differentiation and isotype switching from high throughput BCR sequence datasets obtained from patients in a study of HIV infection and a study of food allergy. These methods are implemented in the R package dowser, available at https://dowser.readthedocs.io.

B cells produce high affinity antibodies through an evolutionary process of mutation and selection during adaptive immune responses. Migration between tissues, differentiation to cellular subtypes, and switching between different antibody isotypes can be important factors in shaping the role B cells play in response to infection, autoimmune disease, and allergies. B cell receptor (BCR) sequence data has the potential to elucidate important information about these processes. However, there is currently no robust, generalizable framework for making such inferences from BCR sequence data. Here, we develop three parsimony-based summary statistics to characterize migration, differentiation, and isotype switching along B cell phylogenetic trees. We confirm the effectiveness of our approach using simulations and further apply our method to data from patients with HIV and allergy. Our methods are released in the R package dowser: https://dowser.readthedocs.io.

Adaptive Immune Receptor Repertoire (AIRR) Community Guide to Repertoire Analysis

Adaptive immune receptor repertoires (AIRRs) are rich with information that can be mined for insights into the workings of the immune system. Gene usage, CDR3 properties, clonal lineage structure, and sequence diversity are all capable of revealing the dynamic immune response to perturbation by disease, vaccination, or other interventions. Here we focus on a conceptual introduction to the many aspects of repertoire analysis and orient the reader toward the uses and advantages of each. Along the way, we note some of the many software tools that have been developed for these investigations and link the ideas discussed to chapters on methods provided elsewhere in this volume.

Human B cell lineages associated with germinal centers following influenza vaccination are measurably evolving

The poor efficacy of seasonal influenza virus vaccines is often attributed to pre-existing immunity interfering with the persistence and maturation of vaccine-induced B cell responses. We previously showed that a subset of vaccine-induced B cell lineages are recruited into germinal centers (GCs) following vaccination, suggesting that affinity maturation of these lineages against vaccine antigens can occur. However, it remains to be determined whether seasonal influenza vaccination stimulates additional evolution of vaccine-specific lineages, and previous work has found no significant increase in somatic hypermutation among influenza-binding lineages sampled from the blood following seasonal vaccination in humans. Here, we investigate this issue using a phylogenetic test of measurable immunoglobulin sequence evolution. We first validate this test through simulations and survey measurable evolution across multiple conditions. We find significant heterogeneity in measurable B cell evolution across conditions, with enrichment in primary response conditions such as HIV infection and early childhood development. We then show that measurable evolution following influenza vaccination is highly compartmentalized: while lineages in the blood are rarely measurably evolving following influenza vaccination, lineages containing GC B cells are frequently measurably evolving. Many of these lineages appear to derive from memory B cells. We conclude from these findings that seasonal influenza virus vaccination can stimulate additional evolution of responding B cell lineages, and imply that the poor efficacy of seasonal influenza vaccination is not due to a complete inhibition of vaccine-specific B cell evolution.

Ymir: A 3D structural affinity model for multi-epitope vaccine simulations

Vaccine development is challenged by the hierarchy of immunodominance between target antigen epitopes and the emergence of antigenic variants by pathogen mutation. The strength and breadth of antibody responses relies on selection and mutation in the germinal center and on the structural similarity between antigens. Computational methods for assessing the breadth of germinal center responses to multivalent antigens are critical to speed up vaccine development. Yet, such methods have poorly reflected the 3D antigen structure and antibody breadth. Here, we present Ymir, a new 3D-lattice-based framework that calculates in silico antibody-antigen affinities. Key physiological properties naturally emerge from Ymir such as affinity jumps, cross-reactivity, and differential epitope accessibility. We validated Ymir by replicating known features of germinal center dynamics. We show that combining antigens with mutated but structurally related epitopes enhances vaccine breadth. Ymir opens a new avenue for understanding vaccine potency based on the structural relationship between vaccine antigens.

Using B cell receptor lineage structures to predict affinity

We are frequently faced with a large collection of antibodies, and want to select those with highest affinity for their cognate antigen. When developing a first-line therapeutic for a novel pathogen, for instance, we might look for such antibodies in patients that have recovered. There exist effective experimental methods of accomplishing this, such as cell sorting and baiting; however they are time consuming and expensive. Next generation sequencing of B cell receptor (BCR) repertoires offers an additional source of sequences that could be tapped if we had a reliable method of selecting those coding for the best antibodies. In this paper we introduce a method that uses evolutionary information from the family of related sequences that share a naive ancestor to predict the affinity of each resulting antibody for its antigen. When combined with information on the identity of the antigen, this method should provide a source of effective new antibodies. We also introduce a method for a related task: given an antibody of interest and its inferred ancestral lineage, which branches in the tree are likely to harbor key affinity-increasing mutations? We evaluate the performance of these methods on a wide variety of simulated samples, as well as two real data samples. These methods are implemented as part of continuing development of the partis BCR inference package, available at https://github.com/psathyrella/partis. Comments Please post comments or questions on this paper as new issues at https://git.io/Jvxkn.

Survival analysis of DNA mutation motifs with penalized proportional hazards

Antibodies, an essential part of our immune system, develop through an intricate process to bind a wide array of pathogens. This process involves randomly mutating DNA sequences encoding these antibodies to find variants with improved binding, though mutations are not distributed uniformly across sequence sites. Immunologists observe this nonuniformity to be consistent with "mutation motifs", which are short DNA subsequences that affect how likely a given site is to experience a mutation. Quantifying the effect of motifs on mutation rates is challenging: a large number of possible motifs makes this statistical problem high dimensional, while the unobserved history of the mutation process leads to a nontrivial missing data problem. We introduce an ℓ 1-penalized proportional hazards model to infer mutation motifs and their effects. In order to estimate model parameters, our method uses a Monte Carlo EM algorithm to marginalize over the unknown ordering of mutations. We show that our method performs better on simulated data compared to current methods and leads to more parsimonious models. The application of proportional hazards to mutation processes is, to our knowledge, novel and formalizes the current methods in a statistical framework that can be easily extended to analyze the effect of other biological features on mutation rates.

Kappa chain maturation helps drive rapid development of an infant HIV-1 broadly neutralizing antibody lineage

HIV-infected infants develop broadly neutralizing plasma responses with more rapid kinetics than adults, suggesting the ontogeny of infant responses could better inform a path to achievable vaccine targets. Here we reconstruct the developmental lineage of BF520.1, an infant-derived HIV-specific broadly neutralizing antibody (bnAb), using computational methods developed specifically for this purpose. We find that the BF520.1 inferred naive precursor binds HIV Env. We also show that heterologous cross-clade neutralizing activity evolved in the infant within six months of infection and that, ultimately, only 2% SHM is needed to achieve the full breadth of the mature antibody. Mutagenesis and structural analyses reveal that, for this infant bnAb, substitutions in the kappa chain were critical for activity, particularly in CDRL1. Overall, the developmental pathway of this infant antibody includes features distinct from adult antibodies, including several that may be amenable to better vaccine responses.