Many-molecule encapsulation by an icosahedral shell.

We computationally study how an icosahedral shell assembles around hundreds of molecules. Such a process occurs during the formation of the carboxysome, a bacterial microcompartment that assembles around many copies of the enzymes ribulose 1,5-bisphosphate carboxylase/ oxygenase and carbonic anhydrase to facilitate carbon fixation in cyanobacteria. Our simulations identify two classes of assembly pathways leading to encapsulation of many-molecule cargoes. In one, shell assembly proceeds concomitantly with cargo condensation. In the other, the cargo first forms a dense globule; then, shell proteins assemble around and bud from the condensed cargo complex. Although the model is simplified, the simulations predict intermediates and closure mechanisms not accessible in experiments, and show how assembly can be tuned between these two pathways by modulating protein interactions. In addition to elucidating assembly pathways and critical control parameters for microcompartment assembly, our results may guide the reengineering of viruses as nanoreactors that self-assemble around their reactants.

Introduction

Encapsulation is a hallmark of biology. A cell must co-localize high concentrations of enzymes and reactants to perform the reactions that sustain life, and it must safely store genetic material to ensure long-term viability. While lipid-based organelles primarily fulfill these functions in eukaryotes, self-assembling protein shells take the lead in simpler organisms. For example, viruses surround their genomes with a protein capsid, while bacteria use large icosahedral shells known as bacterial microcompartments (BMCs) to sequester the enzymes and reactions responsible for particular metabolic pathways (Kerfeld et al., 2010; Axen et al., 2014; Shively et al., 1998; Bobik et al., 1999; Erbilgin et al., 2014; Petit et al., 2013; Price and Badger, 1991; Shively et al., 1973; Shively et al., 1973; Kerfeld and Erbilgin, 2015). Within diverse bacteria, BMC functions have been linked to bacterial growth, carbon fixation, symbiosis, or pathogenesis (Kerfeld and Erbilgin, 2015). Other protein-based compartments are found in bacteria and archea (e.g. encapsulins (Sutter et al., 2008) and gas vesicles (Pfeifer, 2012; Sutter et al., 2008)) and even eukaryotes (e.g. vault particles (Kickhoefer et al., 1998)), while some viruses may assemble around lipidic globules (Lindenbach and Rice, 2013; Faustino et al., 2014). Thus, understanding the factors that control microcompartment assembly and encapsulation is a central question in modern cell biology. From the perspectives of synthetic biology and nanoscience, there is great interest in reengineering BMCs or viruses as nanoreactors that spontaneously encapsulate enzymes and reagents in vitro (e.g. Luque et al., 2014; Douglas and Young, 1998; Rurup et al., 2014; Patterson et al., 2014; Patterson et al., 2012; Zhu et al., 2014; Rhee et al., 2011; Rurup et al., 2014; Wörsdörfer et al., 2012; Comas-Garcia et al., 2014), or as customizable organelles that assemble around a programmable set of core enzymes in vivo, introducing capabilities such as carbon fixation or biofuel production into bacteria or other organisms (e.g. Kerfeld and Erbilgin, 2015; Bonacci et al., 2012; Parsons et al., 2010; Choudhary et al., 2012; Lassila et al., 2014). However, the principles controlling such co-assembly processes have yet to be established, and it is not clear how to design systems to maximize encapsulation.

In this article we take a step toward this goal, by developing theoretical and computational models that describe the dynamical encapsulation of hundreds of cargo molecules by self-assembling icosahedral shells. Although our models are general, we are motivated by recent experiments on a type of BMC known as the carboxysome (Kerfeld et al., 2010; Schmid et al., 2006; Iancu et al., 2007; Tanaka et al., 2008). Carboxysomes are large (40–400 nm), roughly icosahedral shells that encapsulate a dense complex of the enzyme ribulose-1,5-bisphosphate carboxylase/oxygenase (RuBisCO) and other proteins to facilitate the Calvin-Bensen-Bassham cycle in autotrophic bacteria (Price and Badger, 1991; Shively et al., 1973; Shively et al., 1973; Iancu et al., 2007; 2010; Kerfeld et al., 2010; Tanaka et al., 2008). Recently, striking microscopy experiments visualized carboxysome shells assembling on and budding from procarboxysomes (the condensed complex of RuBisCO and other proteins found in the interior of carboxysomes) (Cameron et al., 2013; Chen et al., 2013). Genomic analysis suggests that many BMCs with diverse functions assemble via similar pathways (Cameron et al., 2013; Kerfeld and Erbilgin, 2015). However, the mechanisms of budding and pinch-off to close the shell remain incompletely understood because of the small size and transient nature of assembly intermediates. Moreover, experiments suggest that carboxysomes (another form of carboxysome) assemble by a different mechanism, in which shell assembly encapsulates an initially diffuse pool of RuBisCO (Iancu et al., 2010; Cai et al., 2015). The factors determining which of these assembly pathways occurs are unknown.

BMC assembly is driven by a complex interplay of interactions among the proteins forming the external shell and the interior cargo. It is difficult, with experiments alone, to parse these interactions for those mechanisms and factors that critically influence assembly pathways, especially due to the lack of an in vitro assembly system. Models which can correlate individual factors to their effect on assembly are therefore an important complement to experiments.

Previous experimental and theoretical studies of encapsulation by icosahedral shells, e.g. the assembly of viral capsids around their nucleic acid genomes (e.g. Hu and Shklovskii, 2007; Kivenson and Hagan, 2010; Elrad and Hagan, 2010; Perlmutter et al., 2013; 2014; Mahalik and Muthukumar, 2012; Zhang et al., 2013; Zhang and Linse, 2013; Hagan, 2008; Devkota et al., 2009; Dixit et al., 2006; Borodavka et al., 2012; Dykeman et al., 2013; 2014; Zlotnick et al., 2013; Johnson et al., 2004; Patel et al., 2015; Cadena-Nava et al., 2012; Comas-Garcia et al., 2012; 2014; Garmann et al., 2014a; 2014b; Malyutin and Dragnea, 2013), have demonstrated that the structure of the cargo can strongly influence assembly pathways and products. However, BMCs assemble around a cargo which is topologically different from a nucleic acid — a fluid complex comprising many, noncovalently linked molecules. We demonstrate here that changing the cargo topology leads to new assembly pathways and different critical control parameters.

We present phase diagrams and analysis of dynamical simulation trajectories showing how the thermodynamics, assembly pathways, and emergent structures depend on the interactions among shell proteins and cargo molecules. Within distinct parameter ranges, we observe two classes of assembly pathways, which resemble those suggested for respectively or carboxysomes. We find that tunability of cargo loading is a key functional difference between the two classes of pathways. Shells assembled around a diffuse cargo can be varied from empty (containing almost no cargo) to completely full, whereas assembly around a condensed, procarboxysome-like complex invariably produces full shells. While we find that the encapsulated cargo becomes ordered due to confinement, complete crystalline order in the globule before encapsulation inhibits budding. We discuss these results in the context of recent observations on carboxysome assembly, and their implications for engineering BMCs, viruses or drug delivery vehicles that assemble around a fluid cargo (e.g. Refs. [Kerfeld and Erbilgin, 2015; Parsons et al., 2010; Choudhary et al., 2012; Lassila et al., 2014; Luque et al., 2014; Douglas and Young, 1998; Rurup et al., 2014; Patterson et al., 2014; Patterson et al., 2012; Zhu et al., 2014; Rhee et al., 2011; Rurup et al., 2014; Wörsdörfer et al., 2012]).

Results

Our model system is motivated by icosahedral viral capsids and BMCs (Tanaka et al., 2008; Kerfeld et al., 2010). Since icosahedral symmetry can accommodate at most 60 identical subunits, formation of large icosahedral structures requires subunits to assemble into different local environments. The subunits can be grouped into pentamers and hexamers, with 12 pentamers at the icosahedron vertices and the remaining subunits in hexamers. Viruses typically assemble from small oligomers of the capsid protein, which we refer to as the basic assembly unit (Hagan, 2014). Recent AFM experiments demonstrated that hexamers are the basic assembly unit during the assembly of BMC shell facets (Sutter et al., 2016), and the carboxysome major shell proteins crystallize as pentamers and hexamers (Tanaka et al., 2008). Motivated by these observations, our model considers two basic assembly units, one a pentamer and the other a hexamer, with interactions designed so that the lowest energy structure corresponds to a truncated icosahedron with 12 pentamers and 20 hexamers (Figure 1). While BMCs generally have more hexamers, our model is intended to explore the general principles of assembly around a fluid cargo rather than model a specific system. Further details of the model and a thermodynamic analysis are given in section 3 and the appendices.

Figure 1.

Description of the model.

(A) Each shell subunit contains ‘Attractors’ (green circles) on the perimeter, a ‘Top’ (tan circle, ‘T’) in the center above the plane, and a ‘Bottom’ (purple circle, ‘B’ below the plane). (B) Interactions between complementary Attractors drive subunit dimerization, with the Top-Top repulsions (tan arrow) tuned to favor the subunit-subunit angle in a complete shell. Complementary pairs of attractors are indicated by green arrows in (A) for the pentamer-hexamer interface and in (B) for the hexamer-hexamer interface. (C) Bottom psuedoatoms bind cargo molecules (terra cotta circles, ‘C’), while excluder atoms (blue and brown pseudoatoms in (D)) placed in the plane of the pentagon experience excluded volume interactions with the cargo. (D) The positions of excluder atoms in the lowest energy shell geometry, a truncated icosahedron with 12 pentamers (blue) and 20 hexamers (brown).

DOI: http://dx.doi.org/10.7554/eLife.14078.003

To understand how assembly around multiple cargo molecules depends on the relative strengths of interactions between components, we performed dynamical simulations as a function of the parameters controlling shell subunit-subunit (), shell subunit-cargo (), and cargo-cargo () interaction strengths. All energy values are given in units of the thermal energy, . We focus on parameters for which shell subunit-subunit interactions are too weak to drive assembly in the absence of cargo (). Except where mentioned otherwise, the cargo diameter is set equal to the circumradius of a shell subunit.

For the simulated density of cargo particles, the phase behavior (in the absence of shells) corresponds to a vapor at , liquid-vapor phase coexistence for  (the phase coexistence boundary is slightly below ), and a solid phase at . We find that tuning through phase coexistence dramatically alters the typical assembly process. Strong cargo interactions () drive formation of a globule followed by assembly and budding of a shell, such as observed for carboxysomes (Figure 2A, Simulation Video 1), while under weak interactions () shell assembly usually proceeds in concert with cargo encapsulation (Figure 2B, Simulation Video 2), as suggested for assembly of carboxysomes. We now elaborate on these classes of assembly pathways, and how the resulting assembly products depend on parameter values.

Figure 2.

Snapshots illustrating typical assembly trajectories.

(A) Multi-step assembly involving an amorphous globule of cargo and shell subunits. (B) Single-step assembly, in which shell assembly drives local cargo condensation. and (C) when shell-cargo interactions are too weak to condense the cargo. The values of the cargo-cargo (), shell subunit-cargo (), and subunit-subunit () interaction strengths are listed above each panel (all energies are in units of the thermal energy ), and the time (in units of timesteps) is noted below each image. The color scheme here and throughout the manuscript is: Red=Cargo, Blue=Pentagon Excluder, Brown=Hexagon Excluder. Attractor and Bottom pseudoatoms are omitted to aid visibility. Videos of assembly trajectories are included below.

DOI: http://dx.doi.org/10.7554/eLife.14078.004

(A) Globule-mediated assembly under stronger shell-cargo interactions, leading to the simultaneous formation of multiple globules. Parameters are indicated in the figure, and the initial conditions are random (no pre-equilibration). (B) Assembly around a pre-equilibrated cargo globule (see text and Figure 3—figure supplement 1 for an explanation of this alternative initial condition).

DOI: http://dx.doi.org/10.7554/eLife.14078.005

Motivated by the octomer structure of RuBisCO, we constructed model cargo particles, each containing eight pseudoatoms positioned at the corners of a cube, with center-to-center spacing of . Each pseudoatom interacts with other particles in the system as described for cargo particles in appendix 1. Parameters are (A) , , , and (B) , , , (the magnitude of the shell-cargo and cargo-cargo interaction strengths is reduced because multiple cargo pseudoatoms participate in each interaction). To keep the cargo volume fraction equal to those in other simulations in this manuscript, we set the number of octomer cargo particles to 132 (each has a volume of ). Other parameters are as described in the Model section.

DOI: http://dx.doi.org/10.7554/eLife.14078.006

Video 1.

Animation of a typical simulation showing assembly around a cargo globule.

Parameters are , , and .

DOI: http://dx.doi.org/10.7554/eLife.14078.007

Video 2.

Animation of a typical simulation showing simultaneous assembly and cargo condensation.

Parameters are , , and .

DOI: http://dx.doi.org/10.7554/eLife.14078.008

Figure 3—figure supplement 1.

The distribution of assembly outcomes in Figure 3A is shown as a function of for indicated values of and .

Ten simulations were performed at each set of parameter values. Representative snapshots corresponding to each outcome are shown in Figure 3B. Simulations were performed for timesteps.

DOI: http://dx.doi.org/10.7554/eLife.14078.012

Assembly and budding from a cargo globule

We begin by discussing assembly behavior when the cargo-cargo interactions are strong enough to drive equilibrium phase coexistence (). Near the phase boundary () a system of pure cargo particles is metastable on the timescales we simulate. However, for , adding shell subunits drives nucleation of a cargo globule with shell subunits adsorbed on the surface. The subsequent fate of the globule depends on parameter values; typical simulation end-states are shown as a function of parameter values in Figure 3. For moderate interaction strengths () the globule grows to a large size, typically containing at least twice the cargo molecules that can be packaged within a complete shell. Adsorbed shell subunits then reversibly associate to form ordered clusters. Once a cluster acquires enough inter-subunit interactions to be a stable nucleus, it grows by coagulation of additional subunits or other adsorbed clusters. For the parameter set corresponding to Figure 2A, nucleation is fast in comparison to cluster growth, and thus two nuclei grow simultaneously. The last three images show the system immediately preceding and following detachment of the lower shell. Missing only one of its 32 subunits, the shell is connected to the remainder of the droplet only by a narrow neck of cargo. Insertion of the final subunit breaks the neck and completes shell detachment. The complete shell contains 120–130 cargo particles, which is slighty above random close packing ( particles) but below fcc density ( particles, see appendix 1.2).

Figure 3.

Results of assembly around a cargo globule.

(A) The most frequently observed assembly outcome is overlaid on a color map of the theoretical free energy density difference (Equation (3)) between assembled shells and the unassembled globule. Results are plotted against the shell-cargo adsorption strength and the shell-shell interaction strength for indicated values of the cargo-cargo interaction strength . (B) Representative snapshots of the predominant assembly outcomes shown in (A).

DOI: http://dx.doi.org/10.7554/eLife.14078.009

DOI: http://dx.doi.org/10.7554/eLife.14078.010

The sizes of each cargo globule and shell assemblage, and associations between shell assemblages and cargo globules, were determined by clustering. The outcome was then categorized according to the criteria listed in this table.

DOI: http://dx.doi.org/10.7554/eLife.14078.011

Ten simulations were performed at each set of parameter values. Representative snapshots corresponding to each outcome are shown in Figure 3B. Simulations were performed for timesteps.

DOI: http://dx.doi.org/10.7554/eLife.14078.012

The most frequently observed assembly outcome is overlaid on a color map of the theoretical free energy density difference (Equation (3)) between assembled shells and the unassembled globule. Results are plotted against the shell-cargo adsorption strength and the shell-shell interaction strength for indicated values of the cargo-cargo interaction strength . Outcomes are defined as in Figure 3B of the main text. The outcome of each simulation for this figure is listed in Figure 3—figure supplement 2—source data 1.

DOI: http://dx.doi.org/10.7554/eLife.14078.013

DOI: http://dx.doi.org/10.7554/eLife.14078.014

(A) The mean number of cargo molecules encapsulated by shells assembled in dynamics simulations for . The results are averaged over all complete shells (for any ) assembled at each value of , the error bars indicate 95% confidence intervals. (B) The equilibrium number of cargo particles packaged in shells as a function of the shell-cargo and cargo-cargo interaction strengths. The equilibrium cargo loading was calculated by performing simulations initialized with a pre-assembled shell, for which the excluders on one subunit were made permeable to cargo particles. We then performed two simulations at each parameter set, each of length timesteps.

DOI: http://dx.doi.org/10.7554/eLife.14078.015

List of all simulation outcomes for Figures 3A,5A.

Criteria used to categorize assembly outcomes.

List of all simulation outcomes for Figure 3—figure supplement 1—2.

Increasing the shell-shell interaction strength drives faster shell assembly and closure, thus limiting the size of the globule before budding. For the largest interaction strength we simulated () the globule typically does not exceed the size of a single shell, and multiple globules nucleate within the simulation box (Figure 2—figure supplement 1). This observation could place an upper bound on shell-shell interaction strengths, since multiple nucleation events were rare in the carboxysome assembly experiments (Cameron et al., 2013) (however, we discuss potential complicating factors within the cellular environment below). To quantify the relationship between assembly mechanism and parameter values, we calculate an assembly order parameter, defined as the maximum number of unassembled subunits adsorbed onto a globule during an assembly trajectory. The order parameter is shown as a function of the interaction strengths in Figure 4. For and we observe large values of the order parameter (e.g. , the red and yellow regions in Figure 4), which indicate formation of a large amorphous globule consisent with the procarboxysome precursor to carboxysome shell assembly (Cameron et al., 2013).

Figure 2—figure supplement 1.

Snapshots from additional trajectories, including a trajectory with a pre-equilibrated cargo globule.

(A) Globule-mediated assembly under stronger shell-cargo interactions, leading to the simultaneous formation of multiple globules. Parameters are indicated in the figure, and the initial conditions are random (no pre-equilibration). (B) Assembly around a pre-equilibrated cargo globule (see text and Figure 3—figure supplement 1 for an explanation of this alternative initial condition).

DOI: http://dx.doi.org/10.7554/eLife.14078.005

Figure 4.

Dependence of assembly pathway on shell-cargo and shell-shell interaction strength.

The assembly order parameter, defined as the maximum number of unassembled shell subunits adsorbed on a globule at any point during a trajectory, is shown as a function of and for indicated values of the cargo-cargo interaction . Large numbers of adsorbed unassembled subunits () indicate the two step assembly mechanism (Figure 2A), whereas smaller values correspond to simultaneous assembly and cargo condensation (Figure 2B).

DOI: http://dx.doi.org/10.7554/eLife.14078.016

DOI: http://dx.doi.org/10.7554/eLife.14078.017

Assembly order parameter values for and .

Other assembly products

Outside of the optimal parameter ranges, we observe several classes of alternative outcomes. Overly weak shell-shell interactions fail to drive assembly. For and the cargo vapor phase is metastable, and the system remains ‘Unnucleated’ (with no cargo globule) on simulated timescales (we discuss alternative initial conditions below). Stronger cargo-cargo or shell-cargo interactions result in unassembled ‘Globules’, where a cargo globule forms but the shell subunits on its surface fail to nucleate. As increases, we observe assembly on the globule, leading either to complete shells or two classes of incomplete assembly. In the first incomplete case, ‘Attached’, one or more shells almost reaches completion, but fails to detach from the droplet within simulated timescales. ‘Attached’ configurations occur for low , when the subunit-cargo interaction does not provide a strong enough driving force for the last subunit(s) to penetrate the cargo and close the shell. Overly strong interactions drive the other class of incomplete assembly: ‘Over-nucleated/Malformed’, in which an excess of partially assembled shells deplete the system of free subunits before any shells are completed. In this regime it is also common to observe malformed structures, in which defects become trapped within growing shells.

As the cargo-cargo interaction increases (), multiple effects narrow the parameter range that leads to complete assembly and detachment. Firstly, cargo globules nucleate rapidly at multiple locations within the simulation box, increasing the likelihood of the ‘Over-nucleated’ outcome. Secondly, the threshold value of required for cargo penetration increases, resulting in ‘Attached’ shells over a wider parameter range. We also observe a configuration we refer to as ‘Stalled’, in which shell assembly fails to penetrate the globule surface (and thus does not even proceed to the attached stage). The latter is especially prevalent for , when the cargo crystallizes even in the absence of shell encapsulation. For both ‘Attached’ and ‘Stalled’ configurations, regardless of the initial number of nucleation events, we typically observe coarsening into a large globule.

Simultaneous shell assembly and cargo condensation

For the cargo forms an equilibrium vapor phase in the absence of shell subunits. However, above threshold values of and , the diffuse cargo molecules drive nucleation of shell assembly. The subsequent assembly pathway depends sensitively on the shell-cargo interaction strength. For low (Figure 2C), assembly captures only a few cargo molecules, leading to complete, but nearly empty shells. For larger (Figure 2B, and Simulation Video 2), the shell-cargo interactions drive local condensation of cargo molecules. Shell assembly and cargo complexation then proceed in concert, resembling the mechanism proposed for assembly of -carboxysomes (Iancu et al., 2010). Thus, tuning the shell-cargo interaction dramatically affects cargo loading, with a sharp transition from empty to filled shells around . This transition closely tracks the equilibrium filling fraction (Figure 5C), measured by simulating a complete shell made permeable to cargo molecules. This effect is comparable to the condensation of water vapor below its dew point inside of hydrophilic cavities. In contrast, assembly around a globule only generates full shells.

Figure 5.

Results of assembly around a cargo with weak interactions ().

(A) The most frequently observed assembly outcome as a function of and . The distribution of outcomes for is shown in Figure 3—figure supplement 2, and a data file containing the outcome for each trial at each parameter set is included (Figure 3—source data 1). (B) Representative snapshots for the outcomes shown in (A). The complete shell outcomes are shown with the excluders rendered opaque (left) and transparent (right) to enable visualizing the encapsulated cargo. (C) The number of cargo molecules encapsulated by shells assembled in dynamics simulations (red symbols) is compared to the results of equilibrium simulations (black line). The dynamics results are averaged over all complete shells (for any ) assembled at each value of , the error bars indicate 95% confidence intervals. Most simulations were performed for timesteps; simulations with , , and exhibited partially assembled shells at timesteps, and were continued up to timesteps.

DOI: http://dx.doi.org/10.7554/eLife.14078.018

(A) The color map shows the fraction of subunits in assembled shells () obtained from numerically solving Equation (1–2), using the parameter values described in appendix 2, as a function of and for assembly around gas phase cargo (left) and assembly and budding from a pre-equilibrated globule (right). The white circles overlaid on the plots quantify the fraction of dynamical simulations that led to at least one well-formed capsid (defined as shells containing 12 pentamers and 20 hexamers, each interacting with respectively 5 or 6 neighbors). The size of each white circle is proportional to the yield obtained from dynamical simulations at that parameter set, with the largest circle corresponding to 100%. (B) The color map shows the simulation result for the fraction of subunits in of any type of assemblage (defined as any assemblage comprising 10 or more subunits) as a function of and for (left) and (right). We see that the theoretical prediction of the onset of assembly roughly corresponds to the boundary between assembly and no assembly in the simulations, except that the simulation boundary is seen at slightly higher parameter values in all cases, and for the simulation boundary slopes upward with more sharply than the theoretical prediction. Both of these differences can be attributed to prohibitive nucleation barriers which arise for parameter values near the threshold equilibrium values. As discussed in the main text, decreasing reduces the ability for a partially assembled shells to condense cargo molecules, leading to longer nucleation timescales and hence a wider range of between the equilibrium threshold for assembly and the threshold for observing nucleation within our simulation timescale.

DOI: http://dx.doi.org/10.7554/eLife.14078.019

(A) The fraction of trajectories resulting in assembly of a complete shell is shown as a function of for indicated cargo diameters () and . Each data point corresponds to 10 independent simulations. (B) The maximum number of cargo particles encapsulated into a complete shell for each diameter in the simulations shown in (A). (C) Cutaway view of assembled shells corresponding to each data point in (B). Further information. For most results shown in this article, we set the size of the cargo to be commensurate with the size of the shell subunits. This is qualitatively consistent with BMCs; e.g., for carboxysomes the diameter of a RuBisCO holoenzyme is about twice the circum-diameter of a hexamer or pentamer. To investigate how sensitive cargo encapsulation is to the ratio of cargo and shell subunits sizes, we performed additional simulations with cargo diameters in the range , where is the cargo diameter in Equation (A5). In these simulations we maintained a constant cargo volume fraction and box size, so the number of simulated cargo particles varies inversely with the cargo volume. As shown in (A), assembly can accommodate such variations in the cargo diameter, but the yields and robustness to variations in diminish as varies from 1. This may suggest that commensurate shell subunits and cargo sizes are optimal for encapsulation; however, further exploration is required to determine whether varying other parameters such as the cargo volume fraction or the length scale of the subunit-cargo interaction would change this result.

DOI: http://dx.doi.org/10.7554/eLife.14078.020

Assembly of full shells (by either pathway, Figure 2A or Figure 2B) is typically about two orders of magnitude faster than assembly of empty shells (Figure 2C). This disparity demonstrates the key role that the cargo plays in promoting shell association, during all stages of assembly. Cargo molecules initially promote shell nucleation by stabilizing interactions among small, sub-nucleated clusters. Then, the presence of a condensed globule provides a large cross-section for adsorption of additional subunits, significantly enhancing the flux of subunits to the partial capsid, thus increasing its growth rate. The condensed cargo particularly facilitates insertion of the last few subunits, which are significantly hindered by steric interactions, as noted previously for simulations of empty virus capsids (Nguyen et al., 2007).

Figure 5A shows how the products of assembly around cargo with weak interactions depends on parameters. While moderate parameter values lead to complete assembly, overly weak and (lower left region of Figure 5A) prevent shell nucleation, leading to the ‘Unnucleated’ outcome. In the limit of large but weak the shell-cargo interaction stabilizes small disordered globules ( cargo particles, lower right region of Figure 5A), while under strong subunit and weak cargo interactions (, ) shells nucleate but cannot condense the cargo, leading to the complete but slow assembly just discussed. As for assembly around a globule, overly strong interactions lead to overnucleation and malformed shells. However, the predominant mode of malformation is now shell collapse. Because the cargo is below its dew point, the locally condensed globule leads to a negative pressure on the shell subunits, which can flatten the shell and thus prevent closure of a symmetric shell.

Thermodynamic model

The simple free energy model (Equations (1–2)) reproduces the threshold parameter values required for shell assembly with no adjustable parameters (color map in Figure 3). Since it is an equilibrium model and only considers the free energy difference between complete and unassembled configurations, it cannot distinguish between parameter values that lead to complete assembly or kinetic traps at the long but finite simulation times. However, the thermodynamic calculation does suggest that the simulations resulting in ‘Attached’ shells would eventually reach completion on a longer timescale. We do not show in Figure 5A because the globule is always less favorable than assembled shells for , but the yield of well-formed shells in our simulations roughly follows the prediction of the equilibrium theory (Figure 5—figure supplement 1).

Figure 5—figure supplement 1.

Assembly yields calculated by simulation and theory.

(A) The color map shows the fraction of subunits in assembled shells () obtained from numerically solving Equation (1–2), using the parameter values described in appendix 2, as a function of and for assembly around gas phase cargo (left) and assembly and budding from a pre-equilibrated globule (right). The white circles overlaid on the plots quantify the fraction of dynamical simulations that led to at least one well-formed capsid (defined as shells containing 12 pentamers and 20 hexamers, each interacting with respectively 5 or 6 neighbors). The size of each white circle is proportional to the yield obtained from dynamical simulations at that parameter set, with the largest circle corresponding to 100%. (B) The color map shows the simulation result for the fraction of subunits in of any type of assemblage (defined as any assemblage comprising 10 or more subunits) as a function of and for (left) and (right). We see that the theoretical prediction of the onset of assembly roughly corresponds to the boundary between assembly and no assembly in the simulations, except that the simulation boundary is seen at slightly higher parameter values in all cases, and for the simulation boundary slopes upward with more sharply than the theoretical prediction. Both of these differences can be attributed to prohibitive nucleation barriers which arise for parameter values near the threshold equilibrium values. As discussed in the main text, decreasing reduces the ability for a partially assembled shells to condense cargo molecules, leading to longer nucleation timescales and hence a wider range of between the equilibrium threshold for assembly and the threshold for observing nucleation within our simulation timescale.

DOI: http://dx.doi.org/10.7554/eLife.14078.019

Effects of varying other parameters or initial conditions

To investigate whether the results described above depend on assumptions within our model, we performed several sets of additional simulations. Firstly, we performed simulations in which the ratio between cargo diameter in shell subunit size was varied. As shown in Figure 5—figure supplement 2, assembly is most robust for our default cargo diameter (for which the model was parameterized), but productive assembly occurs for cargo diameters varied over a factor of four. Secondly, we performed assembly simulations with anisotropic cargo molecules with a shape motivated by the octomer structure of the RuBisCO holoenzyme (Figure 2—figure supplement 2).

Figure 5—figure supplement 2.

The effect of varying cargo diameter on assembly.

(A) The fraction of trajectories resulting in assembly of a complete shell is shown as a function of for indicated cargo diameters () and . Each data point corresponds to 10 independent simulations. (B) The maximum number of cargo particles encapsulated into a complete shell for each diameter in the simulations shown in (A). (C) Cutaway view of assembled shells corresponding to each data point in (B). Further information. For most results shown in this article, we set the size of the cargo to be commensurate with the size of the shell subunits. This is qualitatively consistent with BMCs; e.g., for carboxysomes the diameter of a RuBisCO holoenzyme is about twice the circum-diameter of a hexamer or pentamer. To investigate how sensitive cargo encapsulation is to the ratio of cargo and shell subunits sizes, we performed additional simulations with cargo diameters in the range , where is the cargo diameter in Equation (A5). In these simulations we maintained a constant cargo volume fraction and box size, so the number of simulated cargo particles varies inversely with the cargo volume. As shown in (A), assembly can accommodate such variations in the cargo diameter, but the yields and robustness to variations in diminish as varies from 1. This may suggest that commensurate shell subunits and cargo sizes are optimal for encapsulation; however, further exploration is required to determine whether varying other parameters such as the cargo volume fraction or the length scale of the subunit-cargo interaction would change this result.

DOI: http://dx.doi.org/10.7554/eLife.14078.020

Figure 2—figure supplement 2.

Snapshots from assembly trajectories around anisotropic cargo particles, for (A) strong cargo-cargo interactions leading to two-step, globule-mediated assembly, and (B) weak cargo-cargo interactions leading to simultaneous assembly and cargo condensation.

Motivated by the octomer structure of RuBisCO, we constructed model cargo particles, each containing eight pseudoatoms positioned at the corners of a cube, with center-to-center spacing of . Each pseudoatom interacts with other particles in the system as described for cargo particles in appendix 1. Parameters are (A) , , , and (B) , , , (the magnitude of the shell-cargo and cargo-cargo interaction strengths is reduced because multiple cargo pseudoatoms participate in each interaction). To keep the cargo volume fraction equal to those in other simulations in this manuscript, we set the number of octomer cargo particles to 132 (each has a volume of ). Other parameters are as described in the Model section.

DOI: http://dx.doi.org/10.7554/eLife.14078.006

Thirdly, we performed a set of simulations in which we pre-equilibrated the cargo globule before introducing shell subunits into the system (Figure 3—figure supplement 2, Simulation Video 3). Investigating this alternative initial condition was motivated by the fact that RuBisCO is present in the cell before induction of the carboxysome gene in the experiments of Ref. (Cameron et al., 2013), and by the observation that multiple carboxysomes bud sequentially in time from a single procarboxysome. For the results are very similar to those obtained without pre-equilibrating the cargo. However, for , successful assembly and detachment is limited to more narrow ranges of shell-shell and shell-cargo interaction strengths than in Figure 3, due to an increased prevalence of ‘Attached’ and ‘Stalled’ configurations. The latter are particularly common for , when the cargo forms a hexagonally close packed crystal which strongly resists deformation by shell protein assembly.

Figure 3—figure supplement 2.

Results of assembly around a pre-equilibrated cargo globule.

The most frequently observed assembly outcome is overlaid on a color map of the theoretical free energy density difference (Equation (3)) between assembled shells and the unassembled globule. Results are plotted against the shell-cargo adsorption strength and the shell-shell interaction strength for indicated values of the cargo-cargo interaction strength . Outcomes are defined as in Figure 3B of the main text. The outcome of each simulation for this figure is listed in Figure 3—figure supplement 2—source data 1.

DOI: http://dx.doi.org/10.7554/eLife.14078.013

DOI: http://dx.doi.org/10.7554/eLife.14078.014

Video 3.

Animation of a simulation with a pre-equilibrated cargo globule.

Parameters are , , and .

DOI: http://dx.doi.org/10.7554/eLife.14078.021

Taken together, the results from both assembly protocols (Figure 3 and Figure 3—figure supplement 2) suggest that moderate effective cargo-cargo interactions are most consistent with the observations of shell assembly and budding in Refs. (Cameron et al., 2013; Chen et al., 2013). Such interactions are strong enough to drive cargo globule formation, but malleable enough to allow shell assembly to deform and eventually sever intra-globule interactions.

Organization of encapsulated cargo

Studies of assembled carboxysomes report varying degrees of order for the encapsulated cargo, ranging from none to paracrystalline order (Iancu et al., 2007; 2010; Kaneko et al., 2006; Schmid et al., 2006). We therefore studied the relationship between cargo order and interaction parameters using equilibrium simulations (see Figure 6 and Figure 6—figure supplement 1). Below , we do not observe true fcc order of the encapsulated cargo. However, for all parameters leading to significant filling, even those well below the cargo liquid-vapor transition, the cargo becomes organized in concentric layers (Figure 6). We observe similar cargo organizations within shells which have budded from cargo globules in dynamical simulations. These results demonstrate that ordering of the cargo does not require crystallinity of the initial globule. Moreover, the magnitude of ordering increases with cargo loading, but, for fixed loading, is essentially independent of the cargo-shell interaction strength . We observe ordering within filled shells due to confinement, even if even if is set to 0 (Figure 6—figure supplement 1), as previously noted by Iancu et al. (Iancu et al., 2007).

Figure 6.

Order of the encapsulated cargo.

The spherically averaged density of cargo molecules inside a shell is shown as a function of radius for (A) and (B) for indicated values of the cargo-shell adhesion strength , measured in equilibrium simulations. The density of the encapsulated cargo ranges from below random close packing to near hexagonal close packing density as and are increased (see Figure 3—figure supplement 3). A snapshot of cargo inside the shell is shown in Figure 5—figure supplement 2. The raw data for this figure is provided in Figure 6—source data 1.

DOI: http://dx.doi.org/10.7554/eLife.14078.022

DOI: http://dx.doi.org/10.7554/eLife.14078.023

The spherically averaged density distribution is shown as a function of distance from the shell center, for simulations in which a preset number of cargo molecules are trapped within a complete shell, with the cargo-shell attraction turned off (). The value of corresponding to each curve is given in the legend, and the value of the subunit-shell energy is shown above each plot. These simulations were each run for timesteps. The simulations shown in Figure 6 were also run with a complete shell; however, one excluder was rendered permeable to cargo molecules allowing the number of encapsulated cargo molecules to equilibrate. Those simulations were also each run for timesteps. Raw data for this figure is provided in Figure 6—figure supplement 1—source data 1.

DOI: http://dx.doi.org/10.7554/eLife.14078.024

DOI: http://dx.doi.org/10.7554/eLife.14078.025

Figure 6—figure supplement 1.

Ordering of the encapsulated cargo is primarily driven by confinement, not adhesion to the inner surface of the shell. 

The spherically averaged density distribution is shown as a function of distance from the shell center, for simulations in which a preset number of cargo molecules are trapped within a complete shell, with the cargo-shell attraction turned off (). The value of corresponding to each curve is given in the legend, and the value of the subunit-shell energy is shown above each plot. These simulations were each run for timesteps. The simulations shown in Figure 6 were also run with a complete shell; however, one excluder was rendered permeable to cargo molecules allowing the number of encapsulated cargo molecules to equilibrate. Those simulations were also each run for timesteps. Raw data for this figure is provided in Figure 6—figure supplement 1—source data 1.

DOI: http://dx.doi.org/10.7554/eLife.14078.024

DOI: http://dx.doi.org/10.7554/eLife.14078.025

Figure 3—figure supplement 3.

The number of cargo particles packaged as a function of parameters.

(A) The mean number of cargo molecules encapsulated by shells assembled in dynamics simulations for . The results are averaged over all complete shells (for any ) assembled at each value of , the error bars indicate 95% confidence intervals. (B) The equilibrium number of cargo particles packaged in shells as a function of the shell-cargo and cargo-cargo interaction strengths. The equilibrium cargo loading was calculated by performing simulations initialized with a pre-assembled shell, for which the excluders on one subunit were made permeable to cargo particles. We then performed two simulations at each parameter set, each of length timesteps.

DOI: http://dx.doi.org/10.7554/eLife.14078.015

Raw data for Figure 6.

Raw data for Figure 6—figure supplement 1.

Discussion

We have described an equilibrium theory and a dynamical computational model for the assembly of shells around a fluid cargo. Our simulations show that assembly can proceed by two classes of pathways: (i) a multi-step process in which the cargo forms a dense globule, followed by adsorption, assembly, and budding of shell proteins, or (ii) single-step assembly, with simultaneous aggregation of cargo molecules and shell assembly. This result demonstrates that the minimal interactions included in our model are sufficient to drive both classes of assembly pathways, suggesting that they are a generic feature of assembly around a fluid cargo. Moreover, while we cannot rule out the existence of active mechanisms in biological examples such as carboxysomes, our model demonstrates that the same interactions which drive assembly of shells can also drive budding from and closure around an amorphous globule of cargo.

Our results suggest bounds on the relative strengths of interactions that drive BMC assembly in cells. The decisive control parameter determining the assembly pathway is the cohesive energy between cargo molecules, which could arise through direct cargo-cargo interactions or be mediated by auxiliary proteins (Cameron et al., 2013). Relatively weak cargo interactions lead to single-step assembly pathways, while stronger interactions favor formation of the cargo-shell globule. However, the strength of cargo-shell and shell-shell interactions also play a role. Strong shell-shell interactions cause assembly to proceed rapidly during globule formation, limiting the size of the globule. Moreover, if a large globule is already present (e.g. due to time-dependent protein concentrations within a cell), strong interactions tend to drive malformed assemblies. We find that an important functional difference between the two classes of assembly pathways is control over the amount of packaged cargo. While the multi-step assembly pathways always generate a shell filled with cargo molecules, shells assembling around a diffuse cargo can be tuned from nearly empty to completely full by controlling the strength of cargo-shell interactions.

These results have implications for reengineering BMCs to encapsulate new core enzymes. Recent works demonstrated that protein cargos can be targeted to BMCs via encapsulation peptides that mediate cargo-shell interactions. However, packaged amounts were much lower than for native core enzymes (Parsons et al., 2010; Choudhary et al., 2012; Lassila et al., 2014). Our simulations show that both cargo-shell and cargo-cargo interactions (direct or mediated) must be controlled to assemble full shells.

We also find that a general equilibrium theory describes the ranges of parameter values for which assembly occurs. However, the dynamical simulations demonstrate that, at finite timescales, there is a rich variety of assembly morphologies. Formation of ordered, full shells requires a delicate balance of cargo-cargo, cargo-shell, and shell-shell interactions, all of which must be on the order . This constraint is consistent with previous studies on viruses and other assembly systems, which found that formation of ordered states requires multiple, cooperative weak interactions between subunits (Hagan, 2014; Whitelam and Jack, 2015). Outside of optimal parameter regimes, the simulations predict alternative outcomes, ranging from no assembly to various alternative trapped intermediates, with the morphology depending on which interaction is strongest. We find that assembly is least robust to parameter variations when the cargo crystallizes before shell assembly. The assembling shell is unable to deform or penetrate the cargo complex, leading to defect-riddled, non-budded complexes. Within the limits of our simplified model, this observation suggests that procarboxysome complexes are at least partially fluid prior to successful shell assembly. Moreover, we find that observations of ordered cargo within assembled shells may be explained by packing constraints.

An important limitation of the present study is that the model interactions are specific to the shell geometry shown in Figure 1 (containing 20 hexamers) because alternating edges on hexagonal subunits have attractive interactions only with pentagonal subunits. In reality BMCs contain many more hexamers (formed from multiple protein sequences) and thus must include a greater range of hexamer-hexamer interactions. Extension of the model to allow for this possibility would allow consideration of two important questions: (1) The mechanism controlling insertion of the 12 pentagons required for a closed shell topology. (2) The relationship between assembly pathway and BMC size polydispersity. In particular, experiments suggest that -carboxysomes are more polydisperse than -carboxysomes (Price and Badger, 1991; Shively et al., 1973; Shively et al., 1973; Iancu et al., 2007; 2010; Kerfeld et al., 2010; Tanaka et al., 2008). We speculate that in the case of assembly around vapor-phase cargo, the size of the assembling shell will be primarily dictated by the preferred shell protein curvature and thus relatively uniform. However, during assembly around a condensed globule, the shell protein interactions could be strained to accommodate a globule which is larger or smaller than the preferred curvature, causing the shell size to depend on a complex balance of intermolecular interaction strengths and variables such as the local RuBisCO concentration.

Our model is minimal, intended to elucidate general principles of assembly around a fluid cargo, and thus may apply to diverse systems including prokaryotic microcompartments, viruses, and engineered delivery vehicles. The predicted trends for how assembly mechanisms and morphologies vary with control parameters can be experimentally tested by microscopy experiments. Such testing will be most straightforward in vitro (e.g. Luque et al., 2014; Douglas and Young, 1998; Rurup et al., 2014; Patterson et al., 2014; Patterson et al., 2012; Zhu et al., 2014; Rhee et al., 2011; Rurup et al., 2014; Wörsdörfer et al., 2012), where subunit-subunit interactions can be tuned by varying solution conditions and the stoichiometries of shell and cargo species can be readily varied. While there is currently no BMC assembly system starting from purified components, our findings can be tested in vivo by mutations which alter known protein binding interfaces, or by altering expression levels of RuBisCO or carboxysome proteins.

We anticipate that our model can serve as a qualitative guide for understanding how such multicomponent complexes assemble in natural systems, or to reengineer them for new applications. More broadly, our results demonstrate that the properties of encapsulated cargo, such as its topology, geometry and interaction strengths, strongly influence assembly pathways and morphologies.

Materials and methods

Computational model

Shell subunits

We have adapted a model for virus assembly (Perlmutter et al., 2013; 2014; Perlmutter and Hagan, 2015a; Wales, 2005; Fejer et al., 2009; Johnston et al., 2010; Ruiz-Herrero and Hagan, 2015) to describe assembly of an icosahedral shell around a fluid cargo. Each subunit contains ‘Attractors’ on its perimeter that mediate subunit-subunit attractions (as in Ruiz-Herrero and Hagan, 2015). Attractor interactions are specific – complementary pairs of Attractors (see Figure 1A,B and appendix 1) have short-range interactions (modeled by a Morse potential), whereas non-complementary pairs have no interactions. A repulsive interaction between pairs of ‘Top’ (type ‘T’) pseudoatoms favors the correct subunit-subunit angle. The ‘Bottom’ (type ‘B’) pseudoatoms mediate short-ranged subunit-cargo attractions (e.g. due to interactions with shell ‘encapsulation peptides’ (Kinney et al., 2012; Cameron et al., 2013; Fan et al., 2010)), represented by a Morse potential. We also add a layer of ‘Excluders’ in the plane of the ‘Top’ pseudoatoms, which represent subunit-cargo excluded volume interactions. The strengths of subunit-subunit and subunit-cargo attractions are parameterized by potential well depths and respectively (appendix 1).

Cargo

As a minimal representation of globular proteins, the cargo is modeled as spherical particles which interact via an attractive Lennard-Jones (LJ) potential, with well-depth . The attractions implicitly model hydrophobic and screened electrostatics interactions between cargo molecules, as well as effective cargo-cargo interactions mediated by auxiliary proteins (e.g. the carboxysome protein CcmM (Cameron et al., 2013)).

Simulations

We simulated assembly dynamics using the Langevin dynamics algorithm in HOOMD (a software package that uses GPUs to perform highly efficient dynamics simulations [Anderson et al., 2008]) and periodic boundary conditions to represent a bulk system. The subunits are modeled as rigid bodies (Nguyen et al., 2011). The simulations were performed using a set of fundamental units (URL. http://codeblue.umich.edu/hoomd-blue/doc/page_units.html), with defined as the circumradius of the pentagonal subunit (the cargo diameter is also set to 1 ). Unless specified otherwise, each simulation contained enough subunits to form four complete shells (48 pentamers and 80 hexamers) and 611 cargo particles (a shell typically encapsulates 120–130 cargo particles) in a cubic box with side length . The simulation time step was in dimensionless time units, and dynamics was performed for timesteps unless mentioned otherwise.

We performed two sets of simulations, using different initial conditions. In the first, simulations were initialized by introducing cargo particles and shell subunits simultaneously with random positions and orientations (except avoiding high-energy overlaps). The second set of initial conditions was motivated by the possibility that the cargo globule could form before shell subunits reach sufficient concentrations within the cell to undergo assembly. To model this situation, we pre-equilibrated the cargo by performing a long simulation with only cargo particles present. Shell subunits were then introduced with random positions and orientations (excluding high-energy overlaps). For , the assembly simulations thus began with a cargo globule already present. For the two protocols are equivalent, since no globule forms during cargo equilibration.

Sample sizes

To cover the largest range of parameter space possible given the computational expense associated with each simulation, we performed 5 independent simulations at most parameter sets. To assess statistical error and to estimate the distribution of different assembly outcomes, we performed 10 independent trials for one value of at each value of and . We also performed additional simulations at parameter sets for which 5 trials did not result in a majority outcome, or when necessary to obtain better statistics on the number of encapsidated cargo particles. Based on these results, performing additional simulations at other parameter values would not qualitatively change our results. (It would increase the statistical accuracy of estimated boundaries between different outcomes; however, these boundaries correspond to crossovers rather than sharp transitions.)

Thermodynamics of assembly around a fluid cargo

To complement the finite-time simulations, we have developed a general thermodynamic description of assembly around a fluid cargo. We consider shells composed of species , with subunits of species in a complete shell, which encapsulates cargo molecules (the index 0 refers to cargo molecules henceforth). Assembly occurs from a dilute solution of cargo molecules with density , shell subunits with density for each species, and the density of assembled, full shells as . These are in equilibrium with a globule containing cargo molecules and subunits for each species . We assume that, due to the asymmetric nature of the shell-cargo interaction, the shell subunits reside at the exterior of the globule (as we observe in our simulations). The globule containing unassembled shell subunits thus resembles a spherical microemulsion droplet (Safran, 1994). Minimizing the total free energy (see appendix 2) gives:

where is the interaction free energy of the assembled shell and are the chemical potentials of free cargo molecules and shell subunits, given by , with a standard state volume and the globule composition given by

with as the globule free energy.

(1) – (2) are the general equilibrium description for a system of assembling shells with a disordered-phase intermediate; application to a specific system requires specifying the forms of and . In appendix 2 we specify these equations for our computational model, allowing us to compare the equilibrium calculation with simulation results, using no free parameters.

To compare the relative stabilities of the globule and assembled shells, we also calculate the free energy difference

where the first term on the right-hand side is the minimized free energy for a system containing shells and free subunits but no globule, while the second term corresponds to the minimized free energy for a system containing subunits and the globule, but no assembly.

In the interests of transparency, eLife includes the editorial decision letter and accompanying author responses. A lightly edited version of the letter sent to the authors after peer review is shown, indicating the most substantive concerns; minor comments are not usually included.

Thank you for submitting your work entitled "Many-molecule encapsulation by an icosahedral shell" for consideration by eLife. Your article has been reviewed by three peer reviewers, and the evaluation has been overseen by a Reviewing Editor and Naama Barkai as the Senior Editor.

The reviewers have discussed the reviews with one another and the Reviewing Editor has drafted this decision to help you prepare a revised submission.

The following individuals involved in review of your submission have agreed to reveal their identity: Charles Knobler and Avinoam Ben-Shau (peer reviewers).

The manuscript describes novel simulations of encapsulation of hundreds of cargo molecules in icosahedral shells. This is relevant for the formation of the carboxysome microcompartment in bacteria. The work is very impressive from a theoretical standpoint and the tendency is to accept the paper for publication in eLife. However, following the main comment of Reviewer 1, re-enforced by Reviewer 3, there was debate about the suitability of the work for the biological audience of eLife. The overall feeling is that it is, but we suggest that the relevance to biology will be further emphasized. In producing the final draft of the manuscript the authors should also address the following issues:

1) Resolve the confusion about the assembly outcome (see comment of Reviewer 1).

2) Refer to Reviewer 1's point regarding the monodispersity differences between α- and β-carboxysomes.

3) Examine the effect of changing the diameter of the cargo on the results (Reviewers 1 and 2).

4) Explain why the cargo is allowed to equilibrate prior the addition of the shell, in contrast with the biological scenario (Reviewer 1).

5) Consider changing the names of the energy terms (Reviewer 1).

6) Explain why a Morse potential is used for the shell-shell and shell-cargo interactions while a Lennard-Jones potential is used for cargo-cargo interactions (Reviewer 1).

7) Consider adding the effect of non-isotropic interactions among the cargo molecules (Reviewer 1).

8) Cite the Sutter et al. article.

9) Address the point about making the model geometry more permissive (Reviewer 1).

10) Address the issue about Figure 4 (Reviewer 1).

11) Consider adding the analogy to the wetting behavior of binary, immiscible, liquid droplets but not at the expense of making the paper less clear to biologists (Reviewer 3).

Reviewer #1:

Bacterial microcompartments (BMCs) are capsid-like structures that function as organelles in bacteria. They are composed of both shell proteins and specific cargo proteins that self-assemble into the ~100 nm structures found in bacterial cells. Recent experimental evidence suggests that different BMCs may assemble through different mechanisms. Here, Perlmutter provide a compelling theoretical evidence of such mechanisms. The authors demonstrate that this minimal set of interactions is sufficient to drive efficient encapsulation, and they explore the phase diagrams with respect to these parameters. Two distinct paths to well-formed full shells are identified which are consistent with current hypotheses for α- and β-carboxysome assembly, respectively. The balance of energetics they identify could provide a useful conceptual framework for future work on in vitro reconstitution of existing BMCs and the engineering of synthetic ones. The manuscript is well written, and the work appears to be competently executed and analyzed.

My main concern is about the suitability of this paper for eLife and its readership. With its emphasis on statistical mechanical methods and interpretations it seems like this would be a better fit in something like the Biophysical Journal or J. Phys. Chem. B. Moreover, the work is likely to benefit from the exposure in one of these specific journals, as opposed to the more biological readership of eLife.

I have another general concern about the categorization of assembly outcomes, which I found confusing. For instance in Figure 3B 'unnucleated' is shown as a mass of cargo and shell pieces whereas in Figure 5B 'unnucleated' is a collection of separated pieces. It is also not clear to me what exactly 'vapor assembly' means because the section "Vapor Assembly" describes the pathway through which α-carboxysomes are thought to assemble. Are those complete shells filled through this way counted as 'complete' or 'vapor assembly' and what determines the cutoff, the amount of cargo? There also seems to be inconsistent use of symbols where X indicates 'Stalled' in Figure 3 but 'Malformed' in Figure 5. The Results needs a clear and organized description of each of the states as well as consistent usage of images and symbols. A table might be a good idea.

Finally, I wonder if the authors have considered the monodispersity differences between α- and β-carboxysomes. If I recall correctly, β-carboxysomes tend to be more heterogeneous between related species and even within the same cell, suggesting some fundamental difference. Does the current model provide any insight into this?

Reviewer #2:

Mike Hagan and his collaborators have carried out what to my mind are the best computer studies of viral assembly by using minimal models that capture the essential physics of the assembly process. Here they carry out a computational; study that addresses a related problem, the assembly of carboxysomes. The model and computational methods they employ are similar to those they have utilized in previous studies but with an additional degree of complexity. In all of their previous work, I believe, they have used models for the protein shell that is appropriate to a T =1 capsid, which can be made up of rigid capsomers of a single type, either triangles or pentamers. Here they employ a model that is appropriate to a T = 3 shell that is made up of rigid pentamers and hexamers. While the interactions between subunits and between subunits and cargo are modeled in a similar fashion in all of these simulations, the introduction of two different subunits considerably complicates the computations. To my untutored understanding of the nuts and bolts of computer modeling, this seems to me to be a tour de force.

The work breaks new ground as well in considering the effects of the phase behavior of the cargo on the assembly path. While the results are very interesting, there are no surprises. The way in which the assembly pathway depends on the strengths of the interactions is what one would have been expected after being informed by the previous computational studies and in vitro assembly studies of viral assembly. They are by no means trivial, however and especially useful because there are no in vitro studies of carboxysome assembly. I think that this is an excellent paper.

Reviewer #3:

In this paper Perlmutter et al. present an interesting computational study of a self-assembly process whereby many small (spherical) cargo particles are encapsulated by icosahedral shells composed of 12 pentagons and 20 hexagons, analogous to the T=3 protein shells of viral capsids. Using Brownian dynamics simulations the authors analyze the assembly ("phase") behavior of a mixture enabling the formation of up to four complete carboxysome-like particles (CLP). The simulations are accompanied by an equilibrium thermodynamic theoretical analysis. Three interaction parameters dictate the assembly behavior of the system ϵ(LJ), ϵ(sub), ϵ(Ads) representing the (magnitudes) of the interaction energies between neighboring cargo units, capsid subunits, and subunit-cargo pairs, respectively. As expected, large ϵ(LJ) leads to solidification of the cargo particles, preventing CLP formation. Similarly large ϵsub (not simulated) would lead to the formation of empty capsids. The interesting regime is when all energies are comparable, not too small and not too large (compared to kT), in which case two scenarios are possible: (i) co-assembly of the cargo and shell and (ii) aggregation of the cargo particles, followed by "carving" their aggregates to form a complete CLPs (provided subunit interactions are strong enough) – as illustrated in the relevant cases in Figures 2,3, and convincingly summarized in the Discussion.

Reading the Introduction and the Discussion I conclude that the simulations and their analysis are of considerable biological relevance and interest. However, not being a biologist myself, I cannot judge to what extent. On the other hand, as far as the simulations are concerned, my impression is that they are highly non-trivial and of high quality, of the same high standards as the simulations of viral assembly that Hagan et al. have pioneered. The complementing thermodynamic analysis adds a significant theoretical aspect to this work.

Reading the manuscript it also reminded me the (indeed, incomplete) analogy to the wetting behavior of binary, immiscible, liquid droplets. For such droplets of A and B molecules, depending on the interaction strengths ϵ(AA), ϵ(BB), ϵ(AB) complete segregation takes place when (ϵ(AA)+ ϵ(BB))/2>ϵ(AB), and ϵ(AA)> ϵ(AB), ϵ(BB)> ϵ(AB), and envelopment of A by B takes place when (ϵ(AA)+ ϵ(BB))/2>ϵ(AB) > ϵ(BB), etc. The analogy is incomplete because of the special symmetry of the subunit-subunit and the anisotropic subunit-cargo interactions, but some similarity exists and may worth being pointed out. There are few other points which could help the reader (and may already be present in the manuscript but I missed).For instance: we are told that there are 611 cargo particles of a given diameter and we know the inner volume of the capsids. It would be of interest to know what is the maximal number of close-packed spheres can be accommodated within the capsid's volume, and what is the actual density observed. One of the figures shows layering of the cargo particles – at what density relative to the maximal density in the capsid (or in fcc solid) this happens? Other numbers of interest – again, perhaps mentioned – are the critical density and interaction energy of the cargo spheres. Maybe a list of all the "limiting values" of this kind should be mentioned at the beginning of the Results section.

I wonder whether all 78 references are really needed. My guess is that the list can be more economical.

But – to summarize. I find this paper interesting, original and of high scientific quality and recommend its publication in eLife.

1) Resolve the confusion about the assembly outcome (see comment of Reviewer 1).

2) Refer to Reviewer 1's point regarding the monodispersity differences between α- and β-carboxysomes.

3) Examine the effect of changing the diameter of the cargo on the results (Reviewers 1 and 2).

4) Explain why the cargo is allowed to equilibrate prior the addition of the shell, in contrast with the biological scenario (Reviewer 1).

5) Consider changing the names of the energy terms (Reviewer 1).

6) Explain why a Morse potential is used for the shell-shell and shell-cargo interactions while a Lennard-Jones potential is used for cargo-cargo interactions (Reviewer 1).

7) Consider adding the effect of non-isotropic interactions among the cargo molecules (Reviewer 1).

8) Cite the Sutter et al. article.

9) Address the point about making the model geometry more permissive (Reviewer 1).

10) Address the issue about

11) Consider adding the analogy to the wetting behavior of binary, immiscible, liquid droplets but not at the expense of making the paper less clear to biologists (Reviewer 3).

12) Address the typos and minor points raised by reviewers 2 and 3.

We have revised our attached manuscript as suggested by the reviewers. We appreciate both the positive evaluations of our manuscript and the many thoughtful suggestions for its improvement by Prof. Knobler, Prof. Ben-Shaul and the anonymous reviewer. We believe that addressing these comments has significantly improved our presentation. It was clear from the reviewer comments that, given the complex variety of structures that can assemble in the system, we needed to do a better job of clearly presenting our assembly outcomes. We have reformatted our figures to address this. Secondly, there were a number of suggestions for interesting additional simulations (indeed many of these were already on our list of things to do). We have performed these simulations to the extent possible within the two months resubmission time table, leading to several new figures.

Reviewer #1:

Bacterial microcompartments (BMCs) are capsid-like structures that function as organelles in bacteria. They are composed of both shell proteins and specific cargo proteins that self-assemble into the ~100 nm structures found in bacterial cells. Recent experimental evidence suggests that different BMCs may assemble through different mechanisms. Here, Perlmutter provide a compelling theoretical evidence of such mechanisms. The authors demonstrate that this minimal set of interactions is sufficient to drive efficient encapsulation, and they explore the phase diagrams with respect to these parameters. Two distinct paths to well-formed full shells are identified which are consistent with current hypotheses for α- and β-carboxysome assembly, respectively. The balance of energetics they identify could provide a useful conceptual framework for future work on in vitro reconstitution of existing BMCs and the engineering of synthetic ones. The manuscript is well written, and the work appears to be competently executed and analyzed.

My main concern is about the suitability of this paper for eLife and its readership. With its emphasis on statistical mechanical methods and interpretations it seems like this would be a better fit in something like the Biophysical Journal or J. Phys. Chem. B. Moreover, the work is likely to benefit from the exposure in one of these specific journals, as opposed to the more biological readership of eLife.

While we did initially consider sending our work to a more physical sciences oriented journal, such as JACS or PRL, we felt that we would get more readership from the biology community in eLife. It has been our experience that the biology community is interested in and very capable of understanding statistical mechanics arguments when they are properly presented. Moreover, we are particularly interested in having a conversation with the biology community about the questions that our paper addresses, and eLife is an excellent platform for this.

We have revised the text in several ways to clarify the biological relevance of our study. In the introduction we now point out that, in a system as complex as BMCs, quantitative models are an essential complement to experiments. This is especially true given the lack of a complete in vitro assembly system. In the results we now explain the relationship between the biological scenario of BMC assembly and the simulation protocols we have used. We have also used an additional protocol (see below) which may connect better to the cellular situation. In the Discussion we point out aspects of our results which set bounds on the relative strengths of interactions that can drive BMC assembly in cells. We also point out limitations of the current model which should be addressed in the future to more completely describe the biological situation, such as the question of size polydispersity.

I have another general concern about the categorization of assembly outcomes, which I found confusing. For instance in

It is apparent from this comment and those of the other reviewers that the rich variety of configurations which assemble in this system was not adequately categorized and illustrated by our set of snapshots. We have now implemented a new, hopefully clearer, set of category names, we use a unique symbol for each different category across all figures, and we have added a table giving a short description for each category. We have also changed the category name that was formerly ‘vapor assembly’ to avoid confusion with the other assembly pathway.

Finally, I wonder if the authors have considered the monodispersity differences between α- and β-carboxysomes. If I recall correctly, β-carboxysomes tend to be more heterogeneous between related species and even within the same cell, suggesting some fundamental difference. Does the current model provide any insight into this?

We do think that the different assembly pathways between α- and β- carboxysomes play a role in the size-dispersity of carboxysomes. In the case of assembly around vapor-phase cargo, the size of the assembling carboxysome will be primarily dictated by the preferred carboxysome shell protein curvature. However, during assembly around a pre-formed globule, the shell protein interactions could be strained to accommodate a globule which is larger or smaller than the preferred curvature. The latter affect depends on the protein mechanics and orientational specificity of their interactions. For the present study, in order to focus on the simplest possible case, we designed model interactions which are size-specific and thus suppress polymorphism. Polymorphism is still possible in our model, but likely requires larger strain to the interactions than would be the case for more promiscuous interactions. We are currently investigating the best way to extend the model to allow for greater polymorphism. This requires significant exploration and model space however, and is thus beyond the scope of our current manuscript. We have added a paragraph in the Discussion mentioning the potential link between assembly pathway and size-polydispersity, and the caveat that the specificity of our current model design limits polydispersity.

Reviewer #3:

In this paper Perlmutter et al. present an interesting computational study of a self-assembly process whereby many small (spherical) cargo particles are encapsulated by icosahedral shells composed of 12 pentagons and 20 hexagons, analogous to the T=3 protein shells of viral capsids. Using Brownian dynamics simulations the authors analyze the assembly ("phase") behavior of a mixture enabling the formation of up to four complete carboxysome-like particles (CLP).

[…]

Reading the manuscript it also reminded me the (indeed, incomplete) analogy to the wetting behavior of binary, immiscible, liquid droplets. For such droplets of A and B molecules, depending on the interaction strengths ϵ(AA), ϵ(BB), ϵ(AB) complete segregation takes place when (ϵ(AA)+ ϵ(BB))/2>ϵ(AB), and ϵ(AA)> ϵ(AB), ϵ(BB)> ϵ(AB), and envelopment of A by B takes place when (ϵ(AA)+ ϵ(BB))/2>ϵ(AB) > ϵ(BB), etc. The analogy is incomplete because of the special symmetry of the subunit-subunit and the anisotropic subunit-cargo interactions, but some similarity exists and may worth being pointed out.

The reviewer raises an interesting point here about the role of interaction asymmetry versus like-like/likenonlike interaction strength disparity, which is likely to be of interest to physical scientists familiar with the statistical mechanics of binary mixtures. However, we do fear that this Discussion might be distracting to some biologists who are less familiar with these ideas. As a compromise, we have decided to point out the similarity between a globule with unassembled shell subunits on its surface and surfactants stabilizing a spherical microemulsion droplet, where the asymmetry of the surfactant molecules plays a key role in driving them to the interface. This is in the subheading “Thermodynamics of assembly around a fluid cargo” of the Model Section.

There are few other points which could help the reader (and may already be present in the manuscript but I missed).For instance: we are told that there are 611 cargo particles of a given diameter and we know the inner volume of the capsids. It would be of interest to know what is the maximal number of close-packed spheres can be accommodated within the capsid's volume, and what is the actual density observed. One of the figures shows layering of the cargo particles – at what density relative to the maximal density in the capsid (or in fcc solid) this happens? Other numbers of interest – again, perhaps mentioned – are the critical density and interaction energy of the cargo spheres. Maybe a list of all the "limiting values" of this kind should be mentioned at the beginning of the Results section.

We list the location of the binodal as well as the phase behavior corresponding to each of the simulated values of ϵCC at the beginning of Results. We have now added a section to the appendix in which we estimate the number of packaged cargo particles corresponding to fcc or rcp packing. We have added a plot showing the number of packaged cargo particles as a function of parameter values (Figure 3—figure supplement 3), and Figure 5C also shows the number of packaged cargo particles for vapor phase assembly. As suggested by the reviewer, we now refer to these plots in Figure 4, so that the reader is able to place the observation of layering in the context of packing density. (We observe layering even below rcp density, but, as would be expected, the degree of ordering increases as the packaged cargo approaches fcc density.)

I wonder whether all 78 references are really needed. My guess is that the list can be more economical.

There is an extensive background literature both in the field of BMCs and in relating modeling efforts. We feel it is important to point to these works, both for the benefit of the reader who is unfamiliar with BMCs and to acknowledge previous work.

Table 1.

Description of the assembly outcomes presented in Figures 3,5.

DOI: http://dx.doi.org/10.7554/eLife.14078.026

Symbol	Name	Description
▪	Complete shell (full)	Complete shell, full of cargo molecules
◆	Complete shell (empty)r	Complete shell, almost empty of cargo molecules
⚫	Attached	Nearly complete shells attached to a globule by a neck of cargo
✳	Over-nucleated/Malformed	Multiple globules, with incomplete or malformed shells on their surfaces
×	Stalled	Large globule with multiple incomplete or malformed shells on its surface
□	Globule	Cargo globule with unassembled shell subunits on its surface
⊙	Unnucleated	Diffuse subunits and cargo molecules