Is catalytic activity of chaperones a selectable trait for the emergence of heat shock response?

Although heat shock response is ubiquitous in bacterial cells, the underlying physical chemistry behind heat shock response remains poorly understood. To study the response of cell populations to heat shock we employ a physics-based ab initio model of living cells where protein biophysics (i.e., folding and protein-protein interactions in crowded cellular environments) and important aspects of proteins homeostasis are coupled with realistic population dynamics simulations. By postulating a genotype-phenotype relationship we define a cell division rate in terms of functional concentrations of proteins and protein complexes, whose Boltzmann stabilities of folding and strengths of their functional interactions are exactly evaluated from their sequence information. We compare and contrast evolutionary dynamics for two models of chaperon action. In the active model, foldase chaperones function as nonequilibrium machines to accelerate the rate of protein folding. In the passive model, holdase chaperones form reversible complexes with proteins in their misfolded conformations to maintain their solubility. We find that only cells expressing foldase chaperones are capable of genuine heat shock response to the increase in the amount of unfolded proteins at elevated temperatures. In response to heat shock, cells' limited resources are redistributed differently for active and passive models. For the active model, foldase chaperones are overexpressed at the expense of downregulation of high abundance proteins, whereas for the passive model; cells react to heat shock by downregulating their high abundance proteins, as their low abundance proteins are upregulated.

The influence of selection for protein stability on dN/dS estimations

Understanding the relative contributions of various evolutionary processes-purifying selection, neutral drift, and adaptation-is fundamental to evolutionary biology. A common metric to distinguish these processes is the ratio of nonsynonymous to synonymous substitutions (i.e., dN/dS) interpreted from the neutral theory as a null model. However, from biophysical considerations, mutations have non-negligible effects on the biophysical properties of proteins such as folding stability. In this work, we investigated how stability affects the rate of protein evolution in phylogenetic trees by using simulations that combine explicit protein sequences with associated stability changes. We first simulated myoglobin evolution in phylogenetic trees with a biophysically realistic approach that accounts for 3D structural information and estimates of changes in stability upon mutation. We then compared evolutionary rates inferred directly from simulation to those estimated using maximum-likelihood (ML) methods. We found that the dN/dS estimated by ML methods (ωML) is highly predictive of the per gene dN/dS inferred from the simulated phylogenetic trees. This agreement is strong in the regime of high stability where protein evolution is neutral. At low folding stabilities and under mutation-selection balance, we observe deviations from neutrality (per gene dN/dS > 1 and dN/dS < 1). We showed that although per gene dN/dS is robust to these deviations, ML tests for positive selection detect statistically significant per site dN/dS > 1. Altogether, we show how protein biophysics affects the dN/dS estimations and its subsequent interpretation. These results are important for improving the current approaches for detecting positive selection.

A simulated intermediate state for folding and aggregation provides insights into ΔN6 β2-microglobulin amyloidogenic behavior

A major component of ex vivo amyloid plaques of patients with dialysis-related amyloidosis (DRA) is a cleaved variant of β₂-microglobulin (ΔN6) lacking the first six N-terminal residues. Here we perform a computational study on ΔN6, which provides clues to understand the amyloidogenicity of the full-length β₂-microglobulin. Contrary to the wild-type form, ΔN6 is able to efficiently nucleate fibrillogenesis in vitro at physiological pH. This behavior is enhanced by a mild acidification of the medium such as that occurring in the synovial fluid of DRA patients. Results reported in this work, based on molecular simulations, indicate that deletion of the N-terminal hexapeptide triggers the formation of an intermediate state for folding and aggregation with an unstructured strand A and a native-like core. Strand A plays a pivotal role in aggregation by acting as a sticky hook in dimer assembly. This study further predicts that the detachment of strand A from the core is maximized at pH 6.2 resulting into higher aggregation efficiency. The structural mapping of the dimerization interface suggests that Tyr10, His13, Phe30 and His84 are hot-spot residues in ΔN6 amyloidogenesis.

Dialysis-related amyloidosis (DRA) is a conformational disease that affects individuals undergoing long-term haemodialysis. In DRA the progressive accumulation of protein human β₂-microglobulin (Hβ₂m) in the osteoarticular system, followed by its assembly into amyloid fibrils, eventually leads to tissue erosion and destruction. Disclosing the aggregation mechanism of Hβ₂m under physiologically relevant conditions represents a major challenge due to the inability of the protein to efficiently nucleate fibrillogenesis in vitro at physiological pH. On the other hand, ΔN6, a truncated variant of Hβ₂m, which is also a major component of ex vivo amyloid deposits extracted from DRA patients, is able to efficiently form amyloid fibrils de novo in physiological conditions. This amyloidogenic behavior is dramatically enhanced in a slightly more acidic pH (6.2) compatible with the mild acidification that occurs in the synovial fluid of DRA patients. In this work, an innovative three-stage methodological approach, relying on an array of molecular simulations, spanning different levels of resolution is used to investigate the initial stage of the de novo aggregation mechanism of ΔN6 in a physiologically relevant pH range. We identify an intermediate state for folding and aggregation, whose potential to dimerize is enhanced at pH 6.2. Our results provide rationalizations for previous experimental observations and new insights into the molecular basis of DRA.

Influenza A H1N1 pandemic strain evolution--divergence and the potential for antigenic drift variants

The emergence of a novel A(H1N1) strain in 2009 was the first influenza pandemic of the genomic age, and unprecedented surveillance of the virus provides the opportunity to better understand the evolution of influenza. We examined changes in the nucleotide coding regions and the amino acid sequences of the hemagglutinin (HA), neuraminidase (NA), and nucleoprotein (NP) segments of the A(H1N1)pdm09 strain using publicly available data. We calculated the nucleotide and amino acid hamming distance from the vaccine strain A/California/07/2009 for each sequence. We also estimated Pepitope-a measure of antigenic diversity based on changes in the epitope regions-for each isolate. Finally, we compared our results to A(H3N2) strains collected over the same period. Our analysis found that the mean hamming distance for the HA protein of the A(H1N1)pdm09 strain increased from 3.6 (standard deviation [SD]: 1.3) in 2009 to 11.7 (SD: 1.0) in 2013, while the mean hamming distance in the coding region increased from 7.4 (SD: 2.2) in 2009 to 28.3 (SD: 2.1) in 2013. These trends are broadly similar to the rate of mutation in H3N2 over the same time period. However, in contrast to H3N2 strains, the rate of mutation accumulation has slowed in recent years. Our results are notable because, over the course of the study, mutation rates in H3N2 similar to that seen with A(H1N1)pdm09 led to the emergence of two antigenic drift variants. However, while there has been an H1N1 epidemic in North America this season, evidence to date indicates the vaccine is still effective, suggesting the epidemic is not due to the emergence of an antigenic drift variant. Our results suggest that more research is needed to understand how viral mutations are related to vaccine effectiveness so that future vaccine choices and development can be more predictive.

Structural and energetic determinants of co-translational folding

We performed extensive lattice Monte Carlo simulations of ribosome-bound stalled nascent chains (RNCs) to explore the relative roles of native topology and non-native interactions in co-translational folding of small proteins. We found that the formation of a substantial part of the native structure generally occurs towards the end of protein synthesis. However, multi-domain structures, which are rich in local interactions, are able to develop gradually during chain elongation, while those with proximate chain termini require full protein synthesis to fold. A detailed assessment of the conformational ensembles populated by RNCs with different lengths reveals that the directionality of protein synthesis has a fine-tuning effect on the probability to populate low-energy conformations. In particular, if the participation of non-native interactions in folding energetics is mild, the formation of native-like conformations is majorly determined by the properties of the contact map around the tethering terminus. Likewise, a pair of RNCs differing by only 1-2 residues can populate structurally well-resolved low energy conformations with significantly different probabilities. An interesting structural feature of these low-energy conformations is that, irrespective of native structure, their non-native interactions are always long-ranged and marginally stabilizing. A comparison between the conformational spectra of RNCs and chain fragments folding freely in the bulk reveals drastic changes amongst the two set-ups depending on the native structure. Furthermore, they also show that the ribosome may enhance (up to 20%) the population of low energy conformations for chains folding to native structures dominated by local interactions. In contrast, a RNC folding to a non-local topology is forced to remain largely unstructured but can attain low energy conformations in bulk.

Catalysis of protein folding by chaperones accelerates evolutionary dynamics in adapting cell populations

Although molecular chaperones are essential components of protein homeostatic machinery, their mechanism of action and impact on adaptation and evolutionary dynamics remain controversial. Here we developed a physics-based ab initio multi-scale model of a living cell for population dynamics simulations to elucidate the effect of chaperones on adaptive evolution. The 6-loci genomes of model cells encode model proteins, whose folding and interactions in cellular milieu can be evaluated exactly from their genome sequences. A genotype-phenotype relationship that is based on a simple yet non-trivially postulated protein-protein interaction (PPI) network determines the cell division rate. Model proteins can exist in native and molten globule states and participate in functional and all possible promiscuous non-functional PPIs. We find that an active chaperone mechanism, whereby chaperones directly catalyze protein folding, has a significant impact on the cellular fitness and the rate of evolutionary dynamics, while passive chaperones, which just maintain misfolded proteins in soluble complexes have a negligible effect on the fitness. We find that by partially releasing the constraint on protein stability, active chaperones promote a deeper exploration of sequence space to strengthen functional PPIs, and diminish the non-functional PPIs. A key experimentally testable prediction emerging from our analysis is that down-regulation of chaperones that catalyze protein folding significantly slows down the adaptation dynamics.

Molecular chaperones or heat-shock proteins are essential components of protein homeostatic machinery in all three domains of life, whose role is not only to prevent protein aggregation but also catalyze the protein folding process by decreasing the energetic barrier for folding. Importantly, chaperones have often been implicated as phenotypic capacitors since they buffer the deleterious effects of mutations, promote genetic diversity, and thus speed up adaptive evolution. Here we explore computationally the consequences of chaperone activity in cytoplasm via long-time evolutionary dynamics simulations. We use a 6-loci multi scale model of cell populations, where the fitness of each cell is determined from its genome, based on statistical mechanical principles of protein folding and protein-protein interactions. We find that by catalyzing protein folding chaperones buffer the deleterious effect of mutations on folding stability and thus open up a sequence space for efficient and simultaneous optimization of multiple molecular traits determining the cellular fitness. As a result, chaperones dramatically accelerate adaptation dynamics.

Contribution of selection for protein folding stability in shaping the patterns of polymorphisms in coding regions

The patterns of polymorphisms in genomes are imprints of the evolutionary forces at play in nature. In particular, polymorphisms have been extensively used to infer the fitness effects of mutations and their dynamics of fixation. However, the role and contribution of molecular biophysics to these observations remain unclear. Here, we couple robust findings from protein biophysics, enzymatic flux theory, the selection against the cytotoxic effects of protein misfolding, and explicit population dynamics simulations in the polyclonal regime. First, we recapitulate results on the dynamics of clonal interference and on the shape of the DFE, thus providing them with a molecular and mechanistic foundation. Second, we predict that if evolution is indeed under the dynamic equilibrium of mutation-selection balance, the fraction of stabilizing and destabilizing mutations is almost equal among single-nucleotide polymorphisms segregating at high allele frequencies. This prediction is proven true for polymorphisms in the human coding region. Overall, our results show how selection for protein folding stability predominantly shapes the patterns of polymorphisms in coding regions.

Highly abundant proteins favor more stable 3D structures in yeast

To understand the variation of protein sequences in nature, we need to reckon with evolutionary constraints that are biophysical, cellular, and ecological. Here, we show that under the global selection against protein misfolding, there exists a scaling among protein folding stability, protein cellular abundance, and effective population size. The specific scaling implies that the several-orders-of-magnitude range of protein abundances in the cell should leave imprints on extant protein structures, a prediction that is supported by our structural analysis of the yeast proteome.

A one-shot germinal center model under protein structural stability constraints

The germinal center reaction is the process by which low-affinity B cells evolve into potent, immunoglobulin-secreting plasma and memory B cells. Since the recycling hypothesis was created, experimental studies have both tracked movement of a small population of B cells from the light zone into the dark zone, supporting the recycling model, and parallel to the light zone-dark zone interface, indicating a one-way trajectory. We present a novel, sequence-based ab initio model of protein stability and protein interactions. Our model contains a dark zone region of clonal expansion and somatic hypermutation and a light zone site of antigenic selection. We show not only that a one-shot model is sufficient to achieve biologically-realistic rates of affinity growth, population dynamics, and silent:non-silent mutation ratios in the complementary determining region and framework region of antibodies, but also that a stochastic recycling program with or without realistic constraints on the structural stabilities of GC antibodies cannot produce biologically-observed affinity growth, population dynamics or silent:non-silent mutation profiles. The effect of recycling erases affinity gains made by potent antibodies cycling back from the light zone and causes B cells to pool in the dark zone under high replication rates.

Positively selected sites in cetacean myoglobins contribute to protein stability

Since divergence ∼50 Ma ago from their terrestrial ancestors, cetaceans underwent a series of adaptations such as a ∼10-20 fold increase in myoglobin (Mb) concentration in skeletal muscle, critical for increasing oxygen storage capacity and prolonging dive time. Whereas the O2-binding affinity of Mbs is not significantly different among mammals (with typical oxygenation constants of ∼0.8-1.2 µM(-1)), folding stabilities of cetacean Mbs are ∼2-4 kcal/mol higher than for terrestrial Mbs. Using ancestral sequence reconstruction, maximum likelihood and bayesian tests to describe the evolution of cetacean Mbs, and experimentally calibrated computation of stability effects of mutations, we observe accelerated evolution in cetaceans and identify seven positively selected sites in Mb. Overall, these sites contribute to Mb stabilization with a conditional probability of 0.8. We observe a correlation between Mb folding stability and protein abundance, suggesting that a selection pressure for stability acts proportionally to higher expression. We also identify a major divergence event leading to the common ancestor of whales, during which major stabilization occurred. Most of the positively selected sites that occur later act against other destabilizing mutations to maintain stability across the clade, except for the shallow divers, where late stability relaxation occurs, probably due to the shorter aerobic dive limits of these species. The three main positively selected sites 66, 5, and 35 undergo changes that favor hydrophobic folding, structural integrity, and intra-helical hydrogen bonds.

A universal trend among proteomes indicates an oily last common ancestor

Despite progresses in ancestral protein sequence reconstruction, much needs to be unraveled about the nature of the putative last common ancestral proteome that served as the prototype of all extant lifeforms. Here, we present data that indicate a steady decline (oil escape) in proteome hydrophobicity over species evolvedness (node number) evident in 272 diverse proteomes, which indicates a highly hydrophobic (oily) last common ancestor (LCA). This trend, obtained from simple considerations (free from sequence reconstruction methods), was corroborated by regression studies within homologous and orthologous protein clusters as well as phylogenetic estimates of the ancestral oil content. While indicating an inherent irreversibility in molecular evolution, oil escape also serves as a rare and universal reaction-coordinate for evolution (reinforcing Darwin's principle of Common Descent), and may prove important in matters such as (i) explaining the emergence of intrinsically disordered proteins, (ii) developing composition- and speciation-based "global" molecular clocks, and (iii) improving the statistical methods for ancestral sequence reconstruction.

Protein quality control acts on folding intermediates to shape the effects of mutations on organismal fitness

What are the molecular properties of proteins that fall on the radar of protein quality control (PQC)? Here we mutate the E. coli's gene encoding dihydrofolate reductase (DHFR) and replace it with bacterial orthologous genes to determine how components of PQC modulate fitness effects of these genetic changes. We find that chaperonins GroEL/ES and protease Lon compete for binding to molten globule intermediate of DHFR, resulting in a peculiar symmetry in their action: overexpression of GroEL/ES and deletion of Lon both restore growth of deleterious DHFR mutants and most of the slow-growing orthologous DHFR strains. Kinetic steady-state modeling predicts and experimentation verifies that mutations affect fitness by shifting the flux balance in cellular milieu between protein production, folding, and degradation orchestrated by PQC through the interaction with folding intermediates.

Collapse of unfolded proteins in a mixture of denaturants

Both urea and guanidinium chloride (GdmCl) are frequently used as protein denaturants. Given that proteins generally adopt extended or unfolded conformations in either aqueous urea or GdmCl, one might expect that the unfolded protein chains will remain or become further extended due to the addition of another denaturant. However, a collapse of denatured proteins is revealed using atomistic molecular dynamics simulations when a mixture of denaturants is used. Both hen egg-white lysozyme and protein L are found to undergo collapse in the denaturant mixture. The collapse of the protein conformational ensembles is accompanied by a decreased solubility and increased non-native self-interactions of hydrophobic residues in the urea/GdmCl mixture. The increase of non-native interactions rather than the native contacts indicates that the proteins experience a simple collapse transition from the fully denatured states. During the protein collapse, the relatively stronger denaturant GdmCl displays a higher tendency to be absorbed onto the protein surface due to their stronger electrostatic interactions with proteins. At the same time, urea molecules also accumulate near the protein surface, resulting in an enhanced "local crowding" for the protein near its first solvation shell. This rearrangement of denaturants near the protein surface and crowded local environment induce the protein collapse, mainly by burying their hydrophobic residues. These findings from molecular simulations are then further explained by a simple analytical model based on statistical mechanics.

Mutation induced extinction in finite populations: lethal mutagenesis and lethal isolation

Reproduction is inherently risky, in part because genomic replication can introduce new mutations that are usually deleterious toward fitness. This risk is especially severe for organisms whose genomes replicate "semi-conservatively," e.g. viruses and bacteria, where no master copy of the genome is preserved. Lethal mutagenesis refers to extinction of populations due to an unbearably high mutation rate (U), and is important both theoretically and clinically, where drugs can extinguish pathogens by increasing their mutation rate. Previous theoretical models of lethal mutagenesis assume infinite population size (N). However, in addition to high U, small N can accelerate extinction by strengthening genetic drift and relaxing selection. Here, we examine how the time until extinction depends jointly on N and U. We first analytically compute the mean time until extinction (τ) in a simplistic model where all mutations are either lethal or neutral. The solution motivates the definition of two distinct regimes: a survival phase and an extinction phase, which differ dramatically in both how τ scales with N and in the coefficient of variation in time until extinction. Next, we perform stochastic population-genetics simulations on a realistic fitness landscape that both (i) features an epistatic distribution of fitness effects that agrees with experimental data on viruses and (ii) is based on the biophysics of protein folding. More specifically, we assume that mutations inflict fitness penalties proportional to the extent that they unfold proteins. We find that decreasing N can cause phase transition-like behavior from survival to extinction, which motivates the concept of "lethal isolation." Furthermore, we find that lethal mutagenesis and lethal isolation interact synergistically, which may have clinical implications for treating infections. Broadly, we conclude that stably folded proteins are only possible in ecological settings that support sufficiently large populations.

Is there an en route folding intermediate for Cold shock proteins?

Cold shock proteins (Csps) play an important role in cold shock response of a diverse number of organisms ranging from bacteria to humans. Numerous studies of the Csp from various species showed that a two-state folding mechanism is conserved and the transition state (TS) appears to be very compact. However, the atomic details of the folding mechanism of Csp remain unclear. This study presents the folding mechanism of Csp in atomic detail using an all-atom Go model-based simulations. Our simulations predict that there may exist an en route intermediate, in which β strands 1-2-3 are well ordered and the contacts between β1 and β4 are almost developed. Such an intermediate might be too unstable to be detected in the previous fluorescence energy transfer experiments. The transition state ensemble has been determined from the P(fold) analysis and the TS appears even more compact than the intermediate state.

Thermodynamics and kinetics of the hairpin ribozyme from atomistic folding/unfolding simulations

We report a set of atomistic folding/unfolding simulations for the hairpin ribozyme using a Monte Carlo algorithm. The hairpin ribozyme folds in solution and catalyzes self-cleavage or ligation via a specific two-domain structure. The minimal active ribozyme has been studied extensively, showing stabilization of the active structure by cations and dynamic motion of the active structure. Here, we introduce a simple model of tertiary-structure formation that leads to a phase diagram for the RNA as a function of temperature and tertiary-structure strength. We then employ this model to capture many folding/unfolding events and to examine the transition-state ensemble (TSE) of the RNA during folding to its active "docked" conformation. The TSE is compact but with few tertiary interactions formed, in agreement with single-molecule dynamics experiments. To compare with experimental kinetic parameters, we introduce a novel method to benchmark Monte Carlo kinetic parameters to docking/undocking rates collected over many single molecular trajectories. We find that topology alone, as encoded in a biased potential that discriminates between secondary and tertiary interactions, is sufficient to predict the thermodynamic behavior and kinetic folding pathway of the hairpin ribozyme. This method should be useful in predicting folding transition states for many natural or man-made RNA tertiary structures.

The ensemble folding kinetics of the FBP28 WW domain revealed by an all-atom Monte Carlo simulation in a knowledge-based potential

In this work, we apply a detailed all-atom model with a transferable knowledge-based potential to study the folding kinetics of Formin-Binding protein, FBP28, which is a canonical three-stranded β-sheet WW domain. Replica exchange Monte Carlo simulations starting from random coils find native-like (Cα RMSD of 2.68 Å) lowest energy structure. We also study the folding kinetics of FBP28 WW domain by performing a large number of ab initio Monte Carlo folding simulations. Using these trajectories, we examine the order of formation of two β-hairpins, the folding mechanism of each individual β-hairpin, and transition state ensemble (TSE) of FBP28 WW domain and compare our results with experimental data and previous computational studies. To obtain detailed structural information on the folding dynamics viewed as an ensemble process, we perform a clustering analysis procedure based on graph theory. Further, a rigorous P(fold) analysis is used to obtain representative samples of the TSEs showing good quantitative agreement between experimental and simulated Φ values. Our analysis shows that the turn structure between first and second β strands is a partially stable structural motif that gets formed before entering the TSE in FBP28 WW domain and there exist two major pathways for the folding of FBP28 WW domain, which differ in the order and mechanism of hairpin formation.

Non-native interactions play an effective role in protein folding dynamics

Systematic Monte Carlo simulations of simple lattice models show that the final stage of protein folding is an ordered process where native contacts get locked (i.e., the residues come into contact and remain in contact for the duration of the folding process) in a well-defined order. The detailed study of the folding dynamics of protein-like sequences designed as to exhibit different contact energy distributions, as well as different degrees of sequence optimization (i.e., participation of non-native interactions in the folding process), reveals significant differences in the corresponding locking scenarios--the collection of native contacts and their average locking times, which are largely ascribable to the dynamics of non-native contacts. Furthermore, strong evidence for a positive role played by non-native contacts at an early folding stage was also found. Interestingly, for topologically simple target structures, a positive interplay between native and non-native contacts is observed also toward the end of the folding process, suggesting that non-native contacts may indeed affect the overall folding process. For target models exhibiting clear two-state kinetics, the relation between the nucleation mechanism of folding and the locking scenario is investigated. Our results suggest that the stabilization of the folding transition state can be achieved through the establishment of a very small network of native contacts that are the first to lock during the folding process.

Thermal adaptation of viruses and bacteria

A previously established multiscale population genetics model posits that fitness can be inferred from the physical properties of proteins under the physiological assumption that a loss of stability by any protein confers the lethal phenotype to an organism. Here, we develop this model further by positing that replication rate (fitness) of a bacterial or viral strain directly depends on the copy number of folded proteins, which determine its replication rate. Using this model, and both numerical and analytical approaches, we studied the adaptation process of bacteria and viruses at varied environmental temperatures. We found that a broad distribution of protein stabilities observed in the model and in experiment is the key determinant of thermal response for viruses and bacteria. Our results explain most of the earlier experimental observations: the striking asymmetry of thermal response curves; the absence of evolutionary tradeoff, which was expected but not found in experiments; correlation between denaturation temperature for several protein families and the optimal growth temperature of their carrier organisms; and proximity of bacterial or viral optimal growth temperatures to their evolutionary temperatures. Our theory quantitatively and with high accuracy described thermal response curves for 35 bacterial species using, for each species, only two adjustable parameters-the number of rate-determining genes and the energy barrier for metabolic reactions.

De novo design: balancing novelty and confined chemical space

IMPORTANCE OF THE FIELD

De novo drug design serves as a tool for the discovery of new ligands for macromolecular targets as well as optimization of known ligands. Recently developed tools aim to address the multi-objective nature of drug design in an unprecedented manner.

AREAS COVERED IN THIS REVIEW

This article discusses recent advances in de novo drug design programs and accessory programs used to evaluate compounds post-generation.

WHAT THE READER WILL GAIN

The reader is introduced to the challenges inherent in de novo drug design and will become familiar with current trends in de novo design. Furthermore, the reader will be better prepared to assess the value of a tool, and be equipped to design more elegant tools in the future.

TAKE HOME MESSAGE

De novo drug design can assist in the efficient discovery of new compounds with a high affinity for a given target. The inclusion of existing chemoinformatic methods with current structure-based de novo design tools provides a means of enhancing the therapeutic value of these generated compounds.

Optimality of mutation and selection in germinal centers

The population dynamics theory of B cells in a typical germinal center could play an important role in revealing how affinity maturation is achieved. However, the existing models encountered some conflicts with experiments. To resolve these conflicts, we present a coarse-grained model to calculate the B cell population development in affinity maturation, which allows a comprehensive analysis of its parameter space to look for optimal values of mutation rate, selection strength, and initial antibody-antigen binding level that maximize the affinity improvement. With these optimized parameters, the model is compatible with the experimental observations such as the ∼100-fold affinity improvements, the number of mutations, the hypermutation rate, and the “all or none” phenomenon. Moreover, we study the reasons behind the optimal parameters. The optimal mutation rate, in agreement with the hypermutation rate in vivo, results from a tradeoff between accumulating enough beneficial mutations and avoiding too many deleterious or lethal mutations. The optimal selection strength evolves as a balance between the need for affinity improvement and the requirement to pass the population bottleneck. These findings point to the conclusion that germinal centers have been optimized by evolution to generate strong affinity antibodies effectively and rapidly. In addition, we study the enhancement of affinity improvement due to B cell migration between germinal centers. These results could enhance our understanding of the functions of germinal centers.

The antibodies in our immune system could efficiently improve their abilities in recognizing new antigens. This is done with the help of proliferation, mutation and selection of B cells which carry antibodies, but we have difficulties in developing a quantitative description of this adaptation process which is consistent with the various aspects of experimental observations. Based on the knowledge from experiments, here we present a theoretical model to calculate the numbers of B cells with different antigen recognizing abilities all the time, and look for the best possible design that improves the antigen recognizing ability most efficiently. We find that the best possible design is consistent with the experimental observations, pointing to the conclusion that the immune system has been optimized in evolution. We then study the trade-offs leading to the optimization of the design. The results will not only improve our understanding of the functions in immune system, but also reveal the design principles behind the details. In addition, the study enhances our understanding of the population dynamics in evolution.