Open Access Article
This Open Access Article is licensed under a Creative Commons Attribution-Non Commercial 3.0 Unported Licence

Not all therapeutic antibody isotypes are equal: the case of IgM versus IgG in Pertuzumab and Trastuzumab

Firdaus Samsudin a, Joshua Yi Yeo a, Samuel Ken-En Gan *acd and Peter J. Bond *ab
aBioinformatics Institute (A*STAR), 30 Biopolis Street, #07-01 Matrix, Singapore 138671, Singapore. E-mail: peterjb@bii.a-star.edu.sg; samuelg@bii.a-star.edu.sg
bDepartment of Biological Sciences, National University of Singapore, 14 Science Drive 4, Singapore 117543, Singapore
cp53 Laboratory (A*STAR), 8A Biomedical Grove, #06-04/05 Neuros/Immunos, Singapore 138648
dExperimental Drug Development Center (A*STAR), 10 Biopolis Road Chromos #05-01, Singapore 138670

Received 19th September 2019 , Accepted 12th February 2020

First published on 12th February 2020


Abstract

The therapeutic potential of immunoglobulin M (IgM) is of considerable interest in immunotherapy due to its complement-activating and cell-agglutinating abilities. Pertuzumab and Trastuzumab are monoclonal antibodies used to treat human epidermal growth factor receptor 2 (HER2)-positive breast cancer but exhibit significantly different binding affinities as IgM when compared to its IgG isotype. Using integrative multiscale modelling and simulations of complete antibody assemblies, we show that Pertuzumab IgM is able to utilize all of its V-regions to bind multiple HER2 receptors simultaneously, while similar binding in Trastuzumab IgM is prohibited by steric clashes caused by the large globular domain of HER2. This is subsequently validated by confirming that Pertuzumab IgM inhibits proliferation in HER2 over-expressing live cells more effectively than its IgG counterpart and Trastuzumab IgM. Our study highlights the importance of understanding the molecular details of antibody–antigen interactions for the design and isotype selection of therapeutic antibodies.


Introduction

Immunoglobulin M (IgM) is the primary response antibody to combat foreign pathogens in adaptive immunity.1,2 As the first line of antibody defense, IgM tends to have lower antigen binding affinities. To compensate for this, secreted IgM forms multimeric structures (pentamers or hexamers), increasing the number of antigen binding sites for a higher overall avidity. This multimeric characteristic also confers other advantages on IgM. For example, the activation of the classical complement pathway requires the binding of multiple constant fragment (Fc) regions within close proximity, making multimeric IgM a very potent activator of the complement system.3 The large size and multivalency of IgM molecules also enable the formation of bridges between distant epitopes, such as those on different viral particles, leading to superior aggregation properties when neutralizing viral infections.4–6 Whilst most of the currently approved clinical monoclonal antibodies are of the IgG isotype, the high avidity of IgM and its effective complement activation and agglutination make IgM an attractive candidate for future immunotherapy.7

Multimeric IgM exists as either five (pentamer) or six (hexamer) subunits covalently linked to each other via disulfide bridges.8,9 Each IgM subunit is made of four polypeptide chains, namely two heavy chains containing five immunoglobulin (Ig) domains (Cμ1, Cμ2, Cμ3, Cμ4 and VH), and two light chains comprised of two Ig domains (CL and VL) (Fig. S1). A short polypeptide called the joining (J)-chain may also be involved in IgM multimer formation, and the absence of the J-chain has been suggested to favor hexamer formation.10,11 Due to the large size of the IgM pentamer and hexamer, high-resolution structural data for the entire complexes are absent. No crystal structure is currently available for the full-length monomeric IgM, let alone its pentameric or hexameric forms. Early studies based on negative-stain electron microscopy (EM) and small-angle X-ray scattering (SAXS) experiments suggested pentameric IgM to be a symmetric, star-shaped molecule with the antigen-binding fragment (Fab) regions pointing outwards.12–14 Subsequently, cryo-atomic force microscopy (cryo-AFM) data showed the IgM pentamer to be non-planar, forming a mushroom-like shape with part of the Fc domains protruding out of the plane formed by the rest of the antibody.15 A model of IgM Fc was built based on SAXS analysis, integrating structures of each of the Cμ2, Cμ3, and Cμ4 domains solved using X-ray crystallography and NMR spectroscopy.16 Low-resolution cryo-electron tomography (cryo-ET) revealed that both the Fab and Fc domains of IgM are flexible and adopt multiple conformations.17 More recently, EM images indicated that in the presence of the J-chain, the IgM pentamer exhibits an asymmetric pentagonal shape with a large grove, acting as a carrier for apoptosis inhibitors in macrophages.18

Nevertheless, the structural details of interactions between multimeric IgM and antigens remain elusive, largely due to the experimental limitations associated with studying such large complexes. The molecular basis for how IgM achieves its strong avidity is thus unclear. It is currently unknown whether all of the Fab domains in a multimeric IgM are able to bind antigens simultaneously, or whether the binding of an antigen on one Fab arm can affect the binding on another. The degree of avidity of IgM can also vary, especially when binding to different epitopes on the same antigen. For instance, the breast cancer therapeutic antibodies Pertuzumab and Trastuzumab19 showed remarkably different binding avidities to HER2 in their IgM form.20 Controlling for the different epitopes in HER2, there are significant discrepancies between the differences of Pertuzumab IgG1 and Trastuzumab IgG1 binding versus Pertuzumab IgM and Trastuzumab IgM. Compared to its monomeric IgG isotype counterpart, the equilibrium dissociation constant of Pertuzumab IgM to HER2 is around an order of magnitude higher. In contrast, Trastuzumab IgM has a similar equilibrium dissociation constant to HER2 for both IgG and IgM forms, suggesting a much weaker IgM avidity effect in the latter. A molecular-level understanding of how such distinct binding avidities arise for the same antigen is of importance for future design of therapeutic antibodies and epitope selection. Given that we have previously showed that antibody–antigen interactions can be drastically affected by small changes in the antibody light chain,21 antibody hinge,22 V-region pairing,23 and VH families,24 it may be necessary to study the whole IgM molecule using a holistic approach.25

Thus, we now report the first integrative models of full-length Pertuzumab and Trastuzumab IgM multimers, based initially on available X-ray and NMR structures for each Ig domain. The models were validated against previously published EM and cryo-AFM data, while their structural stability and dynamics were assessed using a combination of extensive atomistic and coarse-grained (CG) molecular dynamics (MD) simulations. By modelling the binding of each antibody to the HER2 extracellular domain (ECD), we found that Pertuzumab IgM is able to form a geometrically optimal complex with its HER2 epitopes that enables multivalent interactions. On the other hand, Trastuzumab IgM binding to its HER2 epitope results in steric hindrances that inhibited simultaneous binding of multiple sites. Based on the Trastuzumab and Pertuzumab models, we show that IgM binding avidity depends primarily on the spatial location of the epitope on the antigen to enable efficient interaction with multiple sites at the same time. We subsequently validated this in SKBR3 HER2 over-expressing breast cancer cell-line, where Pertuzumab IgM inhibited SKBR3 proliferation more effectively than Pertuzumab IgG1, whilst the reverse was observed for Trastuzumab. Collectively, our approach represents a novel strategy for providing guidance to epitope selection and antibody isotype selection in the future development of biologics.

Results

Constructing integrative models of Pertuzumab and Trastuzumab IgM

We first built a homology model of a single unit (protomer) for each of Pertuzumab and Trastuzumab IgM (Fig. 1A). The protomer comprises four chains: two heavy and two light. The two heavy chains are covalently linked to each other via a disulfide bond between residues C331 in the Cμ2 domain. The heavy and light chains are linked via a disulfide bridge between residues C131 in the heavy chain and C214 in the C-terminus of the light chain. The structures of individual Ig domains for both antibodies have been resolved by X-ray crystallography and NMR spectroscopy (details in the Methods section). IgM and IgE have a similar domain architecture within the Fc region with three Ig domains instead of two in the other antibody isotypes; we therefore used the structure of IgE Fc26 to model the arrangement of Cμ2, Cμ3 and Cμ4 domains. A previous attempt at modelling the IgM pentamer also utilized the structure of IgE as a template due to its high sequence homology.15 As a result, our model displayed a sharp bend in the Fc part found in IgE. The sharp bend in IgE Fc is in part influenced by an additional disulfide bridge linking the two heavy chains in the Cε2 domains. To test the stability of the bent structure in our IgM model, which lacks this extra disulfide linkage, we ran two independent 1 μs atomistic MD simulations of the Fc domain and compared the dynamics with that of the IgE Fc (Fig. S2). Comparing the asymmetric bending angles of the two heavy chains, we found that both IgM and IgE displayed a similar angle distribution, whereby one of the chains was bent at 60–90° and the other was predominantly bent at ∼120°. These angles remained far below 180° (perfectly flat Fc domain) for both IgM and IgE throughout the entire simulations, suggesting that it is likely that IgM adopts a similarly bent structure to IgE. Additionally, this bend in the Fc domain of the IgM protomer is necessary to form the non-planar conformation with a central protrusion observed in cryo-AFM images of human IgM oligomer complexes.15
image file: c9sc04722k-f1.tif
Fig. 1 The integrative model of Pertuzumab IgM. (A) Pertuzumab IgM protomer with inter-chain disulfide bridges highlighted. (B) Top view of the Pertuzumab IgM hexamer. Disulfide bonds between the Cμ3 domain of adjacent IgM subunits are coloured in red and an example is highlighted in the inset. (C) Side view of the Pertuzumab IgM hexamer to show the protruding Cμ4 domains. Disulfide linkers between the C-terminal domains are illustrated in red and an example is highlighted in the inset. The Fab domains are not shown for clarity. IgM is shown in cartoon representation; the heavy and light chains are colored green and orange, respectively.

To test the structural integrity of the full-length IgM protomer homology models, we performed 1 μs atomistic MD simulations, each system comprising ∼300[thin space (1/6-em)]000 atoms (Table 1). To assess the structural drift of each domain, we measured the root mean square displacement (RMSD) of the backbone atoms with respect to the initial structure. We found that for all Ig domains, the backbone RMSD reached a plateau after around 200 ns of 0.1–0.3 nm (Fig. S3), suggesting that each domain was structurally stable. The radius of gyration (Rg) of each Ig domain also showed little fluctuation during the 1 μs simulations. Analysis of the secondary structure further corroborated this, showing a high degree of structural preservation for all Ig domains throughout (Fig. S4). However, the backbone RMSD of the whole antibody was much higher, plateauing at ∼1 nm in two simulations, but reaching as much as ∼2.5 nm in the third. This indicates that the whole IgM protomer was structurally flexible such that each individual Ig domain could move relative to one another due to the unstructured linkers between them.

Table 1 List of simulations
System contents Resolution Number of atoms/particles Simulation time (μs)
Monomer
Pertuzumab Atomistic 288[thin space (1/6-em)]285 3 × 1
Trastuzumab Atomistic 259[thin space (1/6-em)]405 3 × 1
Pertuzumab Coarse-grained 34[thin space (1/6-em)]153 3 × 1
Trastuzumab Coarse-grained 30[thin space (1/6-em)]982 3 × 1
[thin space (1/6-em)]
Hexamer
Pertuzumab Coarse-grained 245[thin space (1/6-em)]090 3 × 10
Trastuzumab Coarse-grained 266[thin space (1/6-em)]950 3 × 10
Pertuzumab–HER2 Coarse-grained 559[thin space (1/6-em)]507 3 × 10
[thin space (1/6-em)]
Pentamer without J-chain
Pertuzumab Coarse-grained 231[thin space (1/6-em)]746 3 × 10
Trastuzumab Coarse-grained 228[thin space (1/6-em)]580 3 × 10
[thin space (1/6-em)]
Pentamer with a J-chain
Pertuzumab Coarse-grained 285[thin space (1/6-em)]720 3 × 10
Trastuzumab Coarse-grained 292[thin space (1/6-em)]745 3 × 10


After validating the structural stability of the protomer domains, we built integrative models for IgM multimeric states. IgM exists either as pentamers or hexamers, and the former may be with or without the J-chain (a small ∼130 amino acid residue protein involved in IgM assembly). IgM subunits within the oligomers are connected to one another via two disulfide linkages: one between residues C408 in the Cμ3 domains, and another between residues C569 in the C-terminal unstructured loops (Fig. S1). In pentamers with a J-chain, the J-chain replaces the position of the sixth IgM subunit and connects to neighboring IgMs via disulfide bonds to residues C569. In this study, we built representative models for all three oligomeric states of IgM (details in the Methods section) (Fig. 1B and S6). In agreement with cryo-AFM images,15 our resultant IgM oligomer models adopted a non-planar arrangement (Fig. 1C). Due to the acute bend in the Fc region of the protomer, the Cμ4 domains and the C-terminal loop protrude out of the plane defined by the rest of the Fc domains. The J-chain in our pentameric model is also elevated with respect to the rest of the antibody as it is covalently linked to the C-terminal loop. Our model is also in agreement with a previous SAXS-derived model of the IgM Fc hexamer, where the Cμ4 domain formed a hexameric ring that formed the structural core of the whole antibody.16 Whilst the IgM hexamer and pentamer without a J-chain are symmetric, the pentamer with a J-chain form an asymmetric pentagon as previously shown by single-particle negative-stain EM.18 Due to spatial constraints, the Fab domains of one IgM subunit were necessarily positioned within close proximity to adjacent subunits, which could potentially have important implications for simultaneous antigen binding and, therefore, binding avidity.

IgM Fabs are flexible and move independently

IgM pentamers and hexamers are very large protein complexes (around 8000 residues in the former and 10[thin space (1/6-em)]000 residues in the latter), meaning that simulations at atomic resolution would be prohibitively expensive. As such, we next studied the dynamics of IgM pentamers and hexamers at coarse-grained (CG) resolution. An elastic network was applied to maintain the secondary structures within individual Ig domains of each IgM subunit. This elastic network was omitted between Ig domains as our atomistic simulations showed that individual domains may move freely relative to each other. We performed 1 μs CG simulations with a single IgM subunit and found that it showed similar flexibility to the atomistic model, with comparable distributions of relative domain–domain motions (Fig. S5), showing our pseudo-atomic model to be a good approximation of the atomistic counterpart. The backbone RMSD and radius of gyration of each Ig domain throughout these simulations were also similar to the previously described atomistic simulations (Fig. S3). To build CG models for IgM multimers, we also omitted the elastic network between adjacent IgM subunits; they are thus only linked covalently via disulfide bridges. Using these CG models, we then performed three independent repeats of 10 μs simulations per system.

We first looked at the structural dynamics of the IgM multimers. Comparing the per-residue root mean square fluctuation (RMSF) values as a measure of average flexibility along the chains, we found the Fc region to be rigidified due to the disulfide bridges linking each subunit (Fig. 2 and S7). The Fab domains on the other hand were significantly more mobile with respect to their adjacent domains. This is consistent with negative-stain EM images of the monoclonal mouse IgM pentamer, in which a strong density was observed for the Fc chains, whilst the peripheral region corresponding to the Fab domains could not be observed clearly due to their presumed flexible motion.18 Principal component analysis (PCA) was used to assess collective motions for each complex. This revealed that the direction of motion of a given Fab domain is independent of its neighbors (Fig. 3), suggesting allosteric regulation between IgM subunits to be unlikely. Due to their close proximity, adjacent Fab domains could also associate with one another via non-specific but long-timescale interactions, which reduced their mobility. These associations involved different Fab domains in different simulations, rationalizing the varying RMSF patterns of specific Fab domains observed across the three independent repeats (Fig. S7).


image file: c9sc04722k-f2.tif
Fig. 2 Flexibility of the Fab domain in hexameric IgM multimers. Root mean square fluctuation (RMSF) values were mapped onto the protein surface (left) and protein number density (right) from a 10 μs simulation for the (A) Pertuzumab IgM hexamer and (B) Trastuzumab IgM hexamer. The dotted red line shows the approximate boundary between Fc and Fab domains.

image file: c9sc04722k-f3.tif
Fig. 3 The direction of collective motion of the Fab domains during simulations based on principal component analysis (PCA). (A) The first principal motion of all backbone atoms as determined by PCA for the: Pertuzumab IgM hexamer (left); Pertuzumab IgM pentamer without a J-chain (middle); and Pertuzumab IgM pentamer with a J-chain (right), each for one of the simulation replicas. The arrows illustrate how the direction of motion of the Fab domains is independent of adjacent neighbors. (B) The same analysis performed for simulations of Trastuzumab IgM multimers.

Pertuzumab IgM and HER2 ECDs form a geometrically optimal complex

Avidity is an important factor for IgM given its ten or twelve antigen binding sites in its pentameric or hexameric forms, respectively. This contrasts with the presence of only two or four binding sites in other antibody isotypes. Although targeting the same HER2 antigen, Pertuzumab and Trastuzumab bind to different epitopes or sites on the HER2 ECD, and the effect of the different spatial arrangement of binding sites on the avidity of IgM is yet to be investigated. We therefore extended our models of the IgM hexamer and pentamer to integrate the bound HER2 ECD, based on previously determined crystal structures of monomeric Fab–HER2 complexes.27,28

Pertuzumab interacts with HER2 near the center of domain II found in the globular region of the ECD (Fig. S8). This blocks the dimerisation hairpin from fitting into its binding site on adjacent receptors. To build a model of HER2 bound Pertuzumab IgM, we performed structural alignment of each of the Fab domains in the initial integrative model of the Pertuzumab IgM hexamer with the crystal structure (Fig. 4A and B). This model revealed the exposed region, primarily on domain IV, to be small and free of steric clashes. This makes it possible to fit twelve HER2 ECDs to the IgM hexamer, one for each Fab, maximizing the multiple binding sites on IgM, expectedly increasing avidity.


image file: c9sc04722k-f4.tif
Fig. 4 The binding of HER2 ECDs to the Pertuzumab IgM hexamer. (A) Twelve HER2 ECDs (pink and purple) were aligned onto the Pertuzumab IgM hexamer initial model, one for each Fab domain, based on the previously determined binding site.28 (B) Enlarged images of one of the IgM subunits with HER2 ECDs bound. The rest of the IgM is shown in transparent representation and the other HER2 ECDs are not shown for clarity. IgM is shown in ribbon representation and coloured green and orange for heavy and light chains, respectively, and HER2 ECDs are shown in surface representation. (C) Conformations sampled by one of the HER2 ECDs during one of the 10 μs simulations. The bound Fab domain is shown in green (heavy chain) and orange (light chain), and equally spaced conformations of the ECD are shown overlaid, colored according to timeframe. Only the backbone atoms are shown for clarity. (D) Backbone RMSDs of each of the twelve Fab–HER2 complexes during the simulation showed stable binding. (E) Per-residue RMSF of the Pertuzumab IgM hexamer during one of the 10 μs simulations. HER2 molecules are shown in grey. (F) Comparison of backbone RMSDs of the Fc domain of the IgM hexamer from apo versus HER2-bound simulations. Data are shown for three independent repeats of each system.

To test the binding stability, we performed three independent 10 μs simulations of the Pertuzumab IgM hexamer complexed to multiple HER2 ECDs. In all three simulations, all twelve HER2 ECDs remained bound to the antibody; the RMSD of the Fab–HER2 complexes reached a plateau at around 0.5 nm, suggesting stable binding. While the globular region of HER2 ECDs remained relatively rigid due to the binding to IgM, the solvent exposed domain IV was more flexible (Fig. 4C and D). To understand the effect of HER2 binding on the conformational flexibility of Pertuzumab IgM, we measured its per-residue RMSF. Similar to the apo simulations, the Fab regions exhibited more motion compared to the Fc part (Fig. 4E). Interestingly, the RMSD of the Fc domain from all three HER2-bound simulations was comparable to that of the apo simulations (Fig. 4F), suggesting that HER2 binding did not affect the conformation of IgM Fc.

Steric clashes hinder multiple binding of HER2 to Trastuzumab IgM

In contrast to Pertuzumab, Trastuzumab binds HER2 on the C-terminal portion of domain IV, close to its transmembrane domain.27 This binding poses a steric barrier for optimal interaction between the transmembrane regions of neighboring HER2 receptors, which is required for signaling. As performed above, we placed the HER2 ECD onto the initial Trastuzumab IgM hexamer model using structural alignment based on crystallographic data (Fig. 5A and B). We found that the binding of one HER2 ECD to a Fab domain resulted in significant steric clashes with the adjacent Fab domain on the same IgM subunit. As Trastuzumab interacts with the small juxtamembrane region of HER2, the large globular region composed of domains I, II, and III is exposed. Therefore, for HER2 binding to occur, adjacent Fab domains must have sufficient space between them to accommodate this globular region of HER2 ECDs.
image file: c9sc04722k-f5.tif
Fig. 5 The binding of HER2 ECDs to the Trastuzumab IgM hexamer. (A) Twelve HER2 ECDs (pink and purple) aligned onto the initial model of the Trastuzumab IgM hexamer based on the previously elucidated monomeric binding site.27 (B) Enlarged images of one IgM subunit with HER2 bound. Dotted red line highlights the steric clashes between the globular domain of HER2 and the adjacent Fab domain of IgM. Only one HER2 ECD is shown in each image for clarity. IgM is shown in ribbon representation and coloured green and orange for heavy and light chains, respectively, and HER2 ECDs are shown in surface representation. (C) Side (left) and top (right) views of the crystal structure of the Trastuzumab Fab–HER2 ECD complex. Light and heavy chains are orange and green, respectively, while the binding site residues are highlighted in red. The space required to accommodate HER2 ECDs is shown on the right. (D) Distance distribution between the center of mass of the binding site residues on a Fab domain and the bottom-most (left) or the top-most (right) point on the adjacent Fab domain. Different colors depict the ten/twelve different Fab domains in each IgM. The data were extracted from three independent apo simulations of the Trastuzumab IgM hexamer, pentamer, and pentamer with a J-chain. Dotted lines indicate the distance cut-off required for HER2 binding and red arrows show the parts of the distance distributions that meet this requirement.

As the IgM Fab is flexible, it is possible that adjacent Fab domains may move apart from one another to allow HER2 binding. To approximate the likelihood of HER2 binding from our simulations, we first estimated the distance a Fab domain has to move with respect to its neighbor to accommodate HER2. This was performed by measuring the width of the exposed globular domain of HER2 from the crystal structure. We found that a clearance of around 1.2 nm and 1.5 nm on either side of HER2 binding site would allow for binding without a steric clash (Fig. 5C). We then measured distances along the same axis between adjacent Fab domains from our apo simulations to determine how frequently this configuration was achieved. Our apo simulations were 10 μs each with snapshots recorded every 1 ns, giving a total of 10[thin space (1/6-em)]000 frames. As there are twelve or ten Fabs in the IgM hexamer or pentamer, respectively, there were 300[thin space (1/6-em)]000–360[thin space (1/6-em)]000 Fab configurations per triplicate trajectory for each system. The resultant distribution of distances is shown in Fig. 5D. For the Trastuzumab IgM hexamer, the distance requirement to allow HER2 binding was only found ∼500 times, giving rise to a probability of around 1.4 × 10−3. This probability was slightly higher for simulations of IgM pentamers with and without a J chain, with a probability of around 1.7 × 10−2 and 3.7 × 10−2, respectively. Overall, these results further corroborate the low probability of simultaneous binding of HER2 ECDs to Trastuzumab IgM, helping to rationalize the weaker avidity effect compared to Pertuzumab.

Pertuzumab IgM inhibits cell proliferation more effectively than IgG1

Finally, we sought to validate our models by testing the inhibitory effectiveness of the Pertuzumab and Trastuzumab IgM and IgG1 in live cells. Thus, HER2 over-expressing breast cancer cell line SKBR3 cells were grown in the presence of Pertuzumab or Trastuzumab in their IgM or IgG1 forms, or a water control, for 7 days, before counting the viable cells. Even though Pertuzumab and Trastuzumab IgM and IgG1 inhibited cell growth compared to the control, Pertuzumab IgM inhibited cell proliferation significantly more effectively than Pertuzumab IgG1, whereas the reverse was found for Trastuzumab IgM (Fig. 6), thus supporting the hypothesized inhibitory mechanism based on analysis of our models. To rationalize the inhibitory effect of Trastuzumab IgM compared to control, we took a snapshot from the CG simulations whereby the Fab domains had moved apart from one another to allow a single binding of the HER2 ECD (Fig. 5D). Structural alignment was performed as above (Fig. S9) and we found that the binding of a single HER2 ECD to Trastuzumab IgM blocks the HER2 dimerisation interface. The large size of IgM may obstruct another HER2 receptor from approaching the IgM-bound HER2 and therefore dimerisation is expected to be further hindered. This suggests that even though Trastuzumab IgM does not have the ability to bind to multiple HER2 receptors simultaneously like Pertuzumab, the binding to a single HER2 ECD, when sterically allowed, may block receptor dimerisation due to the large size of IgM coupled to its epitope binding mode.
image file: c9sc04722k-f6.tif
Fig. 6 SKBR3 cell count assay. HER2-positive SKBR3 cells were incubated with one of four antibodies (Pertuzumab IgM, Trastuzumab IgM, Pertuzumab IgG1, and Trastuzumab IgG1) and PBS (control). Cell count (cells per mL) of the culture was performed 7 days after incubation with the antibody. Individual standard deviations were used to plot error bars, with means represented by the black circle with their respective values labelled. * represents statistical significance (p < 0.05). Individual standard deviations were used to plot the error bars.

Discussion

We set out to report full-length models of multimeric IgM antibodies that matched existing experimental data, followed by multiscale simulations to study the detailed molecular mechanism of IgM antigen binding avidity. Our models agreed with previously published structural data and showed good overall structural integrity during simulations. We found the Fab domains in multimeric IgM to be flexible whilst the Fc domains were rigid due to the interconnecting disulfide bonds between different subunits. Each Fab arm moved independently of its neighbors suggesting that long-distance allostery between different IgM subunits is unlikely. Models of HER2 binding to IgM hexamers revealed that Pertuzumab could bind to twelve copies of HER2 ECDs simultaneously, explaining the impressive avidity effect compared to its monomeric isotypes.20 Trastuzumab IgM hexamers, however, could not bind to multiple copies of HER2 due to the steric clashes caused by the large globular domain of the latter. The hypothesized inhibitory mechanism for each antibody was then validated by testing their effectiveness in live cells.

There are three possible caveats in our study that we wish to highlight. Firstly, the use of an elastic network model in our CG simulations limits possible large conformational changes in the antibodies. While the Fab arms displayed some degree of flexibility in our simulations, we did not observe any dramatic bending in the Fc part. This behaviour contrasts with cryo-ET structures of IgM, which indicated that the core Fc region could distort resulting in a bell-shaped conformation,17 but is in good agreement with negative-stain EM data.18

Secondly, we have not explored the possible effects of glycosylation on the dynamics of the Fab domains, largely due to the lack of available parameters for the CG simulations of polysaccharides. The heavy chains of IgM have a total of five potential N-linked glycosylation sites,29 of which one is on the Cμ1 domain (residue N163 on our model) that is potentially able to influence the flexibility of the Fab domains. To approximate this effect, we added monosialylated complex glycans – which are the most predominant glycan species – to this site,29 in all of the heavy chains of the atomistic model of Pertuzumab IgM bound to HER2 ECDs (Fig. S10). We found that the glycan moieties lie in between the Fab domains of neighboring IgM subunits, suggesting that they may indeed influence the Fab mobility. Earlier atomistic simulation studies showed that N-glycosylation did not trigger large changes in the overall protein structure but rather reduced protein dynamics.30 This suggests that in the presence of glycosylation, the Fab domains may be less mobile than indicated by our CG simulations. Nevertheless, our model also shows that due to the large distance between the glycosylation sites and the antigen binding sites, it is unlikely for the glycans to affect the ability of Pertuzumab IgM to bind to multiple HER2 ECDs, and therefore should not alter the overall conclusions of the paper.

The third possible caveat is related to the fact that we modelled the J-chain in our IgM pentamers as a random loop, due to the lack of structural data. Previous secondary structural predictions of the J-chain produced conflicting results31,32 and the assigned disulfide bonds were inconsistent with subsequent predictions,33 which makes modelling this domain difficult. Our simulations showed that the J-chain displayed similar rigidity to the Fc domains, likely due to the intradomain disulfide bridges. Adding secondary structural elements to the J-chain will potentially further enhance this rigidity. This rigidity resulted in the Fab domains moving independently of the Fc and J-chain, and as such, we would predict that the secondary structure of the J-chain has only a limited effect on the flexibility of the Fabs and their ability to bind to antigens.

Allosteric modulation has been an area of interest in molecular immunology.21,24,34,35 Previous experimental and computational studies of IgG and IgA have suggested that the Fc and Fab regions of an antibody might communicate with one another;23,25 for example, point mutations in the Fc region may reduce binding affinity, whilst antigen binding to the Fab region changes the conformation of the Fc domain to promote binding to its receptor.22,24,36,37 Allostery between the different subunits in a multimeric IgM, however, remains elusive. Our multiscale simulation study indicates that the individual Fab domains are free to move independently of one another and do not affect the motion of their neighbors. Furthermore, we did not observe any difference in the flexibility of the Fc region upon binding with HER2 compared to the apo systems within the timescale of our simulations, suggesting that allosteric communication between one Fab domain and another via the Fc region is unlikely to occur. As previously discussed, the use of an elastic network model also limits any large conformational changes that may be needed for allostery. Nonetheless, allosteric modulation did not need to be invoked in order to rationalize the experimentally validated HER2 inhibitory activity of each antibody.

Conclusion

Our study highlights the importance of studying the antibody–antigen binding interface in the context of a biologically realistic assembly, specifically in multivalent isotypes like IgM with Fab domains in close proximity to one another. To date, therapeutic antibodies are mostly of the bivalent isotypes, in particular IgG; however, studies on using IgM for clinical purposes are underway especially for cancer treatment.38–41 This is unsurprising given that IgM is immunologically superior to other isotypes, vis-à-vis antigen agglutination and complement activation. The dramatic increase in overall antigen binding strength due to avidity effects is another compelling reason to use IgM. Previous experimental measurement of binding affinities involving IgM and other monomeric isotypes, however, showed that the overall increase in binding avidity in the former is not always significant compared to the latter.20 This suggests that the addition of antigen binding sites in multimeric IgM would not guarantee an improved overall binding. Indeed, our modelling and simulations show that Trastuzumab IgM is not capable of utilizing all of its Fab domains for HER2 binding due to steric clashes, which explains its poorer avidity compared to Pertuzumab IgM. Key to this steric barrier is the Pertuzumab–HER2 binding interface, located in the juxtamembrane region of the latter.27 We postulate that in the presence of the cell membrane and other membrane proteins in vivo, Trastuzumab IgM may be severely hindered from accessing HER2 juxtamembrane regions due to its large size, rendering it a less effective choice compared to the smaller IgG, as shown by the SKBR3 assay, where the Trastuzumab IgM was significantly less effective in inhibiting SKBR3 proliferation than Trastuzumab IgG1. On the other hand, Pertuzumab IgM is able to bind HER2 in all of its antigen binding sites without any steric clashes due to the binding interface being positioned at the large globular domain of HER2.28 This capability is especially beneficial to elicit multiple simultaneous interactions in cancer cells overexpressing HER2 receptors. We do note that this multiple binding will also be dependent upon the density of HER2 receptors on the membrane surface of cancer cells, as well as other additional constraints such as cell shape and membrane undulations. Modelling the binding of IgM antibodies to full-length HER2 receptors within a realistic membrane environment is beyond the scope of the current study but is a potential avenue to explore in the future, given the significant inhibitory effects displayed by Pertuzumab IgM in the SKBR3 assay. Our results shed light on the crucial role of understanding antibody–antigen interactions at the molecular level for more effective antibody isotype selection, as well as epitope selection for biologics development.

Experimental and computational details

Modelling IgM hexamers and pentamers

The homology models of Pertuzumab and Trastuzumab IgM protomers were constructed using Modeller version 9.21.42 The crystal structures of the Fab domains of Pertuzumab and Trastuzumab in complexes with human HER2 ECDs (PDB: 1S78 (ref. 28) and 1N8Z,27 respectively) were used as the template for the Fab region, whilst the Fc region was modelled based on the atomic structures of its three domains: Cμ2 (PDB: 4JVU), Cμ3 (PDB: 4BA8) and Cμ4 (PDB: 4JVW) previously determined using X-ray crystallography and NMR spectroscopy.16 For each antibody, 10 models were built and the model with the lowest discreet optimized protein energy (DOPE) score43 was chosen. Stereochemical assessment using Ramachandran analysis44 of the non-loop regions showed 40 and 31 outlier residues for Pertuzumab and Trastuzumab IgM models, respectively, representing around 2% of the whole protein, confirming that the models are structurally reasonable. Further structural assessment was based on the analysis of stability during all-atom MD simulations. To build models of the hexameric and pentameric IgM, the protomers were duplicated and aligned within PyMOL45 and two disulfide bridges were constructed between adjacent protomers (Fig. S1). In the pentamers with a J-chain, the J-chain is positioned between two IgM protomers and a disulfide bond was built between the J-chain (residues C37 and C91) and C569 of the protomers. There are currently no crystal structures available for the J-chain and a search in the Protein Databank did not reveal any protein with a significant sequence identity. A previous solution structure of human and mouse IgM (PDB: 2RCJ)14 modelled the J-chain as a constant Ig fold domain as predicted by a secondary structure assignment.32 However, the disulfide bonds in this model were later shown to be incorrect.33 As there is currently no credible structural data, the 137-residue J-chain (UniProt: P01591) was modelled as an unstructured loop. Disulfide bonds within the J-chain were constructed as reported by Frutiger et al.;33 these disulfide bridges maintained the rigidity of the J-chain as shown by the similar RMSF values compared to the Fc domains in the pentamer. To model the binding of HER2 ECDs to Pertuzumab and Trastuzumab IgM hexamers, we aligned the respective HER2–Fab complexes to each of the Fab domains onto the initial hexameric models. Missing residues in HER2 ECDs were built as loops using Modeller. Using CHARMM-GUI Glycan Reader & Modeller,46 we modelled the N-glycosylation of the Cμ1 domain on the Pertuzumab IgM hexamer model bound to HER2 ECDs. Monosialylated complex glycans (Fig. S10C), the most prevalent species found on this Ig domain,29 were added to residue N163 on all of the twelve heavy chains. Steepest descent energy minimization was performed in a vacuum to remove any atomic overlaps. The schematic structure of the glycan was drawn using the DrawGlycan-SNFG webserver.47

Simulation systems and protocols

Atomistic MD simulations were performed to test the structural stability of Pertuzumab and Trastuzumab IgM protomers. The antibodies were parameterized using the CHARMM36 force field48 and solvated in TIP3P water molecules with 0.15 M NaCl added to neutralize the system. Short 100 ps equilibration simulations were performed during which the heavy atoms of the antibodies were positionally restrained using a force constant of 1000 kJ mol−1 nm−2. The temperature was maintained at 310 K using the Nose–Hoover thermostat with a time constant of 1.0 ps,49,50 while the pressure was kept at 1 atm by an isotropic coupling to the Parrinello–Rahman barostat with a time constant of 5.0 ps.51 The smooth particle mesh Ewald method with a real-space cut-off of 1.2 nm was utilized to calculate the electrostatic interactions.52 The van der Waal's interactions were truncated at 1.2 nm with a force switch smoothing function applied from 1.0 to 1.2 nm. All covalent bonds were constrained using the LINCS algorithm,53 which allowed for an integration time step of 2 fs. After the equilibration simulations, the positional restraints were removed and three independent 1 μs production simulations were conducted, each starting with a different distribution of atomic velocities.

To understand the dynamics of Pertuzumab and Trastuzumab IgM hexamers and pentamers with and without HER2 ECDs, we performed coarse-grained (CG) MD simulations. Antibodies and HER2 ECDs were converted to CG representation using the MARTINI 2.2 force field with the elastic network ElNeDyn model applied to maintain the secondary structures.54,55 The systems were solvated using standard MARTINI water molecules and neutralized using 0.15 M NaCl. Equilibration simulations were performed for 10 ns whereby positional restraints with a force constant of 1000 kJ mol−1 nm−2 were applied to the backbone atoms of the proteins. The temperature was maintained at 310 K using a velocity-rescaling thermostat with a time constant of 1 ps,56 while the pressure was kept at 1 atm using an isotropic coupling to a Berendsen barostat with a time constant of 5 ps.57 The electrostatics were computed using the reaction field method, whilst the Lennard-Jones potential was cut off using a potential shift Verlet scheme; the short-range cut-off for both was set to 1.1 nm. Following the equilibration simulations, three independent 10 μs simulations were conducted with different starting velocities. For these production runs, the Parrinello–Rahman barostat was utilized with a coupling constant of 12.0 ps.51 The list of simulations performed for this study is available in Table S1. All simulations were performed using GROMACS package version 2018 (ref. 58) and visualized in PyMOL45 and VMD.59 All analysis was performed using GROMACS package version 2018 (ref. 60) or custom Python scripts utilizing the MDAnalysis toolkit.61

Inhibition of the proliferation of HER2-positive SKBR3

The HER2-positive SKBR3 cell line was maintained in T-25 flasks. Cells were cultured in DMEM (Gibco) and seeded into Corning® Costar® 6-well cell culture plates (Costar) at 2 × 105 mL−1. Respective antibodies previously produced (Pertuzumab IgM, Trastuzumab IgM, Pertuzumab IgG1, Trastuzumab IgG1, see Lua et al.20) were added at 1 μg mL−1. After 7 days, viable cells were counted using a haemocytometer with trypan blue stain 0.4% (Invitrogen). The paired sample t-test was utilized to perform statistical analysis.

Conflicts of interest

There are no conflicts to declare.

Author contributions

Conceptualization, S. K.-E. G., and P. J. B.; computational investigation, F. S.; wet lab investigation, J. Y. Y.; writing – original draft, F. S., S. K.-E. G. and P. J. B.; supervision, S. K.-E. G. and P. J. B.

Acknowledgements

This work was supported by A*STAR BII APD Lab core funds. Simulations were performed on the petascale computer cluster ASPIRE-1 at the National Supercomputing Center of Singapore.

References

  1. M. R. Ehrenstein and C. A. Notley, Nat. Rev. Immunol., 2010, 10, 778–786 CrossRef CAS PubMed .
  2. H. Wang, J. E. Coligan and H. C. Morse, Front. Immunol., 2016, 7, 99 Search PubMed .
  3. N. R. Cooper, Adv. Immunol., 1985, 37, 151–216 CAS .
  4. N. Baumgarth, O. C. Herman, G. C. Jager, L. E. Brown, L. A. Herzenberg and J. Chen, J. Exp. Med., 2000, 192, 271–280 CrossRef CAS PubMed .
  5. P. Seiler, U. Kalinke, T. Rulicke, E. M. Bucher, C. Bose, R. M. Zinkernagel and H. Hengartner, J. Virol., 1991, 65, 3001–3006 CrossRef .
  6. M. S. Diamond, E. M. Sitati, L. D. Friend, S. Higgs, B. Shrestha and M. Engle, J. Exp. Med., 2003, 198, 1853–1862 CrossRef CAS PubMed .
  7. A. Kretschmer, R. Schwanbeck, T. Valerius and T. Rösner, Transfus. Med. Hemotherapy, 2017, 44, 320–326 CrossRef PubMed .
  8. A. C. Davis and M. J. Shulman, Immunol. Today, 1989, 10, 118–122 CrossRef CAS PubMed .
  9. J. W. Brewer, T. D. Randall, R. M. E. Parkhouse and R. B. Corley, Immunol. Today, 1994, 15, 165–168 CrossRef CAS PubMed .
  10. A. Cattaneo and M. S. Neuberger, EMBO J., 1987, 6, 2753–2758 CrossRef CAS PubMed .
  11. M. J. Niles, L. Matsuuchi and M. E. Koshland, Proc. Natl. Acad. Sci. U. S. A., 2006, 92, 2884–2888 CrossRef PubMed .
  12. A. Feinstein and E. A. Munn, Nature, 1969, 224, 1307–1309 CrossRef CAS PubMed .
  13. A. C. Davis, K. H. Roux and M. J. Shulman, Eur. J. Immunol., 1988, 18, 1001–1008 CrossRef CAS PubMed .
  14. S. J. Perkins, A. S. Nealis, B. J. Sutton and A. Feinstein, J. Mol. Biol., 1991, 221, 1345–1366 CrossRef CAS PubMed .
  15. D. M. Czajkowsky and Z. Shao, Proc. Natl. Acad. Sci. U. S. A., 2009, 106, 14960–14965 CrossRef CAS PubMed .
  16. R. Muller, M. A. Grawert, T. Kern, T. Madl, J. Peschek, M. Sattler, M. Groll and J. Buchner, Proc. Natl. Acad. Sci. U. S. A., 2013, 110, 10183–10188 CrossRef PubMed .
  17. R. R. Akhouri, S. Goel, H. Furusho, U. Skoglund and M. Wahlgren, Cell Rep., 2016, 14, 723–736 CrossRef CAS PubMed .
  18. E. Hiramoto, A. Tsutsumi, R. Suzuki, S. Matsuoka, S. Arai, M. Kikkawa and T. Miyazaki, Sci. Adv., 2018, 4, eaau1199 CrossRef CAS PubMed .
  19. J. Baselga, J. Cortés, S.-B. Kim, S.-A. Im, R. Hegg, Y.-H. Im, L. Roman, J. L. Pedrini, T. Pienkowski, A. Knott, E. Clark, M. C. Benyunes, G. Ross and S. M. Swain, N. Engl. J. Med., 2012, 366, 109–119 CrossRef CAS PubMed .
  20. W.-H. Lua, W.-L. Ling, J. Y. Yeo, J.-J. Poh, D. P. Lane and S. K.-E. Gan, Sci. Rep., 2018, 8, 1–9 CrossRef PubMed .
  21. C. T.-T. Su, W.-L. Ling, W.-H. Lua, J.-J. Poh and S. K.-E. Gan, Sci. Rep., 2017, 7, 1–7 CrossRef PubMed .
  22. C. T.-T. Su, W.-H. Lua, W.-L. Ling and S. K.-E. Gan, Antibodies, 2018, 7, 20 CrossRef CAS PubMed .
  23. W.-L. Ling, W.-H. Lua, J.-J. Poh, J. Y. Yeo, D. P. Lane and S. K.-E. Gan, Front. Immunol., 2018, 9, 1–11 CrossRef PubMed .
  24. W. H. Lua, C. T.-T. Su, J. Y. Yeo, J.-J. Poh, W.-L. Ling, S.-X. Phua and S. K.-E. Gan, J. Allergy Clin. Immunol., 2019, 144, 514–523 CrossRef CAS PubMed .
  25. S.-X. Phua, K.-F. Chan, C. T.-T. Su, J.-J. Poh and S. K.-E. Gan, Biosci. Rep., 2019, 39, BSR20181958 CrossRef CAS PubMed .
  26. T. Wan, R. L. Beavil, S. M. Fabiane, A. J. Beavil, M. K. Sohi, M. Keown, R. J. Young, A. J. Henry, R. J. Owens, H. J. Gould and B. J. Sutton, Nat. Immunol., 2002, 3, 681–686 CrossRef CAS PubMed .
  27. H.-S. Cho, K. Mason, K. X. Ramyar, A. M. Stanley, S. B. Gabelli and D. J. Leahy, Nature, 2003, 421, 756–760 CrossRef CAS PubMed .
  28. M. Franklin, K. Carey, F. Vajdos, D. Leahy, A. de Vos and M. Sliwkowski, Cancer Cell, 2004, 5, 317–328 CrossRef CAS PubMed .
  29. J. N. Arnold, M. R. Wormald, D. M. Suter, C. M. Radcliffe, D. J. Harvey, R. A. Dwek, P. M. Rudd and R. B. Sim, J. Biol. Chem., 2005, 280, 29080–29087 CrossRef CAS PubMed .
  30. H. S. Lee, Y. Qi and W. Im, Sci. Rep., 2015, 5, 8926 CrossRef PubMed .
  31. G. M. Cann, A. Zaritsky and M. E. Koshland, Proc. Natl. Acad. Sci. U. S. A., 1982, 79, 6656–6660 CrossRef CAS PubMed .
  32. J. Zikan, J. Novotny, T. L. Trapane, M. E. Koshland, D. W. Urry, J. C. Bennett and J. Mestecky, Proc. Natl. Acad. Sci. U. S. A., 1985, 82, 5905–5909 CrossRef CAS PubMed .
  33. S. Frutiger, G. J. Hughes, N. Paquet, J. C. Jaton and R. Liithy, Biochemistry, 1992, 31, 12643–12647 CrossRef CAS PubMed .
  34. A. Janda, A. Bowen, N. S. Greenspan and A. Casadevall, Front. Microbiol., 2016, 7, 22 Search PubMed .
  35. D. Yang, R. Kroe-Barrett, S. Singh, C. J. Roberts and T. M. Laue, mAbs, 2017, 9, 1231–1252 CrossRef CAS PubMed .
  36. J. Schlessinger, I. Z. Steinberg, D. Givol, J. Hochman and I. Pecht, Proc. Natl. Acad. Sci. U. S. A., 1975, 72, 2775–2779 CrossRef CAS PubMed .
  37. M. Oda, H. Kozono, H. Morii and T. Azuma, Int. Immunol., 2003, 15, 417–426 CrossRef CAS PubMed .
  38. L. Rasche, J. Duell, I. C. Castro, V. Dubljevic, M. Chatterjee, S. Knop, F. Hensel, A. Rosenwald, H. Einsele, M. S. Topp and S. Brändlein, Haematologica, 2015, 100, 377–384 CrossRef CAS PubMed .
  39. L. Rasche, J. Duell, C. Morgner, M. Chatterjee, F. Hensel, A. Rosenwald, H. Einsele, M. S. Topp and S. Brändlein, PLoS One, 2013, 8, e63414 CrossRef CAS PubMed .
  40. S. Brandlein, N. Rauschert, L. Rasche, A. Dreykluft, F. Hensel, E. Conzelmann, H.-K. Muller-Hermelink and H. P. Vollmers, Mol. Cancer Ther., 2007, 6, 326–333 CrossRef PubMed .
  41. M. Liedtke, C. J. Twist, B. C. Medeiros, J. R. Gotlib, C. Berube, M. M. Bieber, N. M. Bhat, N. N. Teng and S. E. Coutre, Haematologica, 2012, 97, 30–37 CrossRef CAS PubMed .
  42. A. Sali and T. Blundell, J. Mol. Biol., 1994, 234, 779–815 CrossRef PubMed .
  43. M. Shen, D. Devos, F. Melo and A. Sali, Protein Sci., 2006, 15, 1653–1666 CrossRef PubMed .
  44. G. N. Ramachandran, C. Ramakrishnan and V. Sasisekharan, J. Mol. Biol., 1963, 7, 95–99 CrossRef CAS PubMed .
  45. W. L. DeLano, CCP4 Newsletter On Protein Crystallography, 2002, vol. 40, pp. 82–92 Search PubMed .
  46. S. J. Park, J. Lee, Y. Qi, N. R. Kern, H. S. Lee, S. Jo, I. Joung, K. Joo, J. Lee and W. Im, Glycobiology, 2019, 29, 320–331 CrossRef CAS PubMed .
  47. K. Cheng, Y. Zhou and S. Neelamegham, Glycobiology, 2017, 27, 200–205 CrossRef CAS PubMed .
  48. J. Huang and A. D. MacKerell, J. Comput. Chem., 2013, 34, 2135–2145 CrossRef CAS PubMed .
  49. S. Nosé, J. Phys. Chem. B, 1984, 81, 511–519 CrossRef .
  50. W. G. Hoover, Phys. Rev. A: At., Mol., Opt. Phys., 1985, 31, 1695–1697 CrossRef PubMed .
  51. M. Parrinello, J. Appl. Phys., 1981, 52, 7182 CrossRef CAS .
  52. U. Essmann, L. Perera, M. L. Berkowitz, T. Darden, H. Lee and L. G. Pedersen, J. Chem. Phys., 1995, 103, 8577–8593 CrossRef CAS .
  53. B. Hess, H. Bekker, H. J. C. Berendsen and J. G. E. M. Fraaije, J. Comput. Chem., 1997, 18, 1463–1472 CrossRef CAS .
  54. L. Monticelli, S. K. Kandasamy, X. Periole, R. G. Larson, D. P. Tieleman and S. J. Marrink, J. Chem. Theory Comput., 2008, 4, 819–834 CrossRef CAS PubMed .
  55. X. Periole, M. Cavalli, S. Marrink and M. A. Ceruso, J. Chem. Theory Comput., 2009, 5, 2531–2543 CrossRef CAS PubMed .
  56. G. Bussi, D. Donadio and M. Parrinello, J. Chem. Phys., 2007, 126, 014101 CrossRef PubMed .
  57. H. J. C. Berendsen, J. P. M. Postma, W. F. van Gunsteren, A. DiNola and J. R. Haak, J. Chem. Phys., 1984, 81, 3684 CrossRef CAS .
  58. M. J. Abraham, T. Murtola, R. Schulz, S. Páall, J. C. Smith, B. Hess and E. Lindah, SoftwareX, 2015, 1–2, 19–25 CrossRef .
  59. W. Humphrey and A. Dalke, J. Mol. Graphics, 1996, 15, 33–38 CrossRef .
  60. M. J. Abraham, T. Murtola, R. Schulz, S. Páll, J. C. Smith, B. Hess and E. Lindahl, SoftwareX, 2015, 1, 19–25 CrossRef .
  61. N. Michaud-agrawal, E. J. Denning, T. B. Woolf and O. Beckstein, J. Comput. Chem., 2011, 32, 2319–2327 CrossRef CAS PubMed .

Footnote

Electronic supplementary information (ESI) available. See DOI: 10.1039/c9sc04722k

This journal is © The Royal Society of Chemistry 2020