Introduction

Human immunodeficiency virus type 1 (HIV-1) is the most prevalent strain of HIV responsible for the development of acquired immunodeficiency syndrome (AIDS). The HIV-1 envelop (Env) consists of a host cell-derived lipid membrane and viral glycoproteins that play a crucial role in mediating viral entry into host cells. The Env glycoprotein is initially synthesized in the endoplasmic reticulum (ER) as a precursor gp160 and cleaved by furin into two subunits, gp120 and gp41. The non-covalently associated gp120–gp41 complex is transported to the cell surface in the form of a trimer, where it is subsequently incorporated into the envelope of nascent virions during viral assembly. The exposure of Env protein is essential for binding to the primary receptor CD4 and the co-receptors CCR5 or CXCR4, triggering membrane fusion and viral entry. However, this exposure also renders the virus susceptible to immune attack. In response to host immune pressure, Env is densely coated with N-linked glycans added during post-translational modification in the ER and Golgi apparatus, which effectively shield vulnerable epitopes from immune recognition.

Since HIV-1 was identified as the cause of AIDS in the early 1980s, extensive research has been conducted to elucidate the mechanisms of viral infection and immune response. A major breakthrough was the design of soluble gp140 trimers, composing gp120 and part of the gp41 subunit. Native-like Env trimers, such as SOSIP,¹ native flexibly linked (NFL),² and uncleaved prefusion-optimized (UFO) constructs,³ mimic the structure of viral spike and serve as valuable antigen targets for developing small molecule inhibitors and broadly neutralizing antibodies (bNAbs). With advances in structure determination techniques, such as X-ray crystallography and cryo-electron microscopy (cryo-EM), numerous high-resolution structures of gp140, both unliganded and antibody-bound, are now available in the Protein Data Bank (PDB). Even with such progresses, however, the remaining portions of gp41, including the membrane-proximal external region (MPER), transmembrane domain (TMD), and cytoplasmic tail (CT), are relatively understudied. The MPER is a highly conserved region targeted by several bNAbs, including 10E8, 2F5, 4E10, and Z13e1.^4–7 However, in most studies, MPER is examined as a monomeric peptide bound to antibodies or within membrane mimetics such as bicelles and nanodiscs, which do not capture the structure and dynamics of unliganded MPER in the context of the gp120–gp41 trimer embedded in a membrane bilayer. Despite the critical roles of the MPER and TMD in mediating the fusion of viral and host cell membranes,^8,9 they are often excluded from structural studies due to difficulties in crystallization of hydrophobic TM segments while maintaining their native conformations. Therefore, nuclear magnetic resonance (NMR) spectroscopy remains one of the few viable methods for obtaining structural information of MPER and TMD. However, the NMR studies from different research groups have yielded conflicting conclusions regarding the oligomeric state of the TMD. Reported TMD structures include trimeric coiled coils,^10–13 monomeric helices,¹⁴ and trimers that are not tightly bundled.¹⁵ Similarly, there are different conclusions about the conformation and orientation of MPER. Some studies suggest that gp41 bends at the MPER-TMD boundary (around residue 673), with the entire MPER adopting a helical conformation that lies nearly parallel to the membrane,^11,16 while others show that the C-terminal residues of MPER forms a continuous helix with the TMD.^10,12,13,17 In addition, the exceptionally long CT plays an important role in facilitating the incorporation of the Env glycoprotein into virions,¹⁸ While it is known to contain three conserved amphipathic α-helical segments, referred to as lentiviral lytic peptides (LLPs), proposed models differ in the arrangement of LLPs,^13,19 and the complete structure of CT and its location in the membrane remains inconclusive.

In this work, we built a model of full-length gp120–gp41 trimer embedded in a lipid bilayer mimicking the lipid composition of the mammalian plasma membrane²⁰ (Figure 1). We prepared simulation systems varying in the presence of the cleavage site and CT, as well as the initial position of protein in the membrane. Multiple microsecond-long all-atom molecular dynamics (MD) simulations were performed for each system to explore the motions of individual protein domains and the membrane, and to examine how their conformational variability is affected by the difference in the initial configurations. Our results show that the ectodomain undergoes substantial tilting relative to the membrane plane while maintaining a rigid internal structure. In contrast, the MPER and TMD display highly diverse conformations, and their structural variations are influenced by the presence of the CT and the initial TMD position in the membrane. Moreover, we selected several bNAbs targeting the epitopes across different regions of the Env protein and demonstrate that the simulation trajectories can be used to assess the epitope accessibility.

Figure 1.

Results

The ectodomain maintains a rigid internal structure and tilts independently of the TMD

As described in the Methods section, multiple simulation systems were constructed, varying in cleavage, TMD positioning, and CT truncation. Simulations trajectories are denoted as CH^CT1, UL^CT2, CL^ΔCT3, UH^ΔCT3, etc., where the first letter (C/U) indicates cleaved or uncleaved, the second letter (H/L) indicates the high or low initial TMD position, CT/ΔCT indicates the presence or absence of the CT, and the numeric suffix specifies the trajectory index among three replicas (Table S1). In all simulations, both the ectodomain and TMD adopted variable orientations relative to the bilayer plane. To quantify the tilt of two domains, we defined the tilt angles for the ectodomain (θ_EC) and TMD (θ_TM) (Figure 2A). θ_EC is measured between the bilayer normal and the vector from the center of mass (COM) of G594 (on all three protomers) to the COM of D167, while θ_TM is measured between the bilayer normal and the vector from the COM of I684 to the COM of V705. Across trajectories, θ_EC typically ranges from 0° to 40°, with occasional tilting up to 50°. Representative structures for different θ_EC values are shown in Figure 2A. In contrast to the wide range of θ_EC, θ_TM generally remained within 0° to 20° with rare excursions to 30° observed in a small fraction of trajectories. We split each 1-µs trajectory into four 0.25-µs intervals, and the variations of θ_EC and θ_TM over time show that large conformational changes occurred during the first 0.5 µs, followed by convergence of θ_EC and θ_TM distributions during the second 0.5 µs of most trajectories (Figures S2-S5). For CH^ΔCT, CL^ΔCT, UH^ΔCT and UL^ΔCT systems, the θ_EC and θ_TM distributions calculated from the second 0.5 µs of three independent trajectories are shown in Figure 2B, and the data for CH^CT, CL^CT, UH^CT and UL^CT systems are shown in Figure S6. Although the combined dataset spans a wide range of θ_EC (0°–50°) and θ_TM (0°–30°), each individual trajectory explores only a portion of this space. For instance, CH^ΔCT1 samples large θ_EC with small θ_TM, whereas CH^ΔCT3 samples small θ_EC with medium-to-large θ_TM. Thus, multiple independent simulations are essential to capture a wide region of the conformational space. Most importantly, no consistent correlation was observed between θ_EC and θ_TM, either within single trajectories or across all trajectories combined. Pearson correlation coefficients of θ_EC and θ_TM in single trajectories varies between -0.5 and 0.5, with examples such as CH^ΔCT1 and CH^ΔCT2 showing similar θ_TM but distinct θ_EC values, and CL^ΔCT1 and CL^ΔCT3 showing similar θ_EC but distinct θ_TM values.

Figure 2.

Despite the considerable tilting of the ectodomain relative to the membrane, its internal structure remains rigid and well-preserved throughout the simulations. The root-mean-square fluctuation (RMSF) and root-mean-square deviation (RMSD) of the ectodomain were calculated after the snapshots from each trajectory were aligned to the initial structure by maximizing the overlap in the ectodomain. The majority of the ectodomain displays low RMSF (< 2 Å), and the RMSD stabilized around 4 Å after an initial rise during the first 250 ns (Figures 3A and S7). Higher RMSF values were observed in the residues missing from the cryo-EM structure highlighted in red in Figures 1A and B, and in part of HR1 (Q551-H564) that forms a flexible loop at the interface between two neighboring protomers. In addition, the RMSF of the MPER was also calculated with the trajectories aligned by the ectodomain. The entire MPER (K665–R683) and the adjacent HR2 segment (L660–D664) at the C-terminus of the ectodomain exhibited elevated fluctuations. The gp120–gp41 model was built based on the NMR structure in which the MPER adopts a bent conformation consisting of two helices joined by a sharp turn. The resulting trimer widens from the HR2 helix to the midpoint of the MPER (F673) and narrows from F673 to the TMD. However, this specific conformation was not maintained throughout the simulations. We measured the inter-chain distances between the Cα atoms of corresponding residues (G644, E654, D664, and F673) on neighboring protomers to characterize the structural variation along the HR2 helix and the MPER (Figure 3B). The inter-chain distances of G644 and E654 maintained narrow distributions centered on their initial values, while those of D664 and F673 exhibited broader distributions, reflecting inward shifts of three protomers. This effect was more pronounced in uncleaved systems. In cleaved systems, the HR2 helix interacted with residues M530-N543 and L619-N625 of neighboring protomers, which helps to stabilize the bent MPER conformation (Figures 3C and D). However, these interactions were not consistently observed across all protomers and all simulations, and therefore we still observed considerable conformational variability in this region. In uncleaved systems, the closed loop at the cleavage site occupied the space between the HR2 helix and the neighboring protomer, disrupting their interactions and facilitating inward shifts of the HR2 helix (Figures 3E and F).

Figure 3.

The energetically unfavorable R696 in the hydrophobic core results in asymmetric, kinked TMD conformations and disrupts membrane integrity

Unlike the predominantly hydrophobic TMDs commonly found in many viral envelop proteins, the gp41 TMD contains multiple charged residues: R683 at the N-terminal boundary, R707 and R709 at the C-terminal boundary, and a central arginine, R696 (Figure 4A). In simulations, R683 consistently interacted with lipid headgroups in the exoplasmic leaflet, while R707 and R709 interacted with those in the cytoplasmic leaflet, together acting as anchors that secure the TMD within the bilayer. In the cryo-EM structure (PDB ID: 7LOI), the side chain of R696 forms cation-π interaction with the side chain of F699 and a hydrogen bond with the backbone carbonyl group of L692, but it is oriented outward from the helical bundle. When embedded in the bilayer, however, this configuration is energetically unfavorable, as the positively charged side chain directly contacts hydrophobic lipid tails.²² Early in the simulations, the TMD rapidly rearranged to allow R696 to interact with more favorable partners, including negatively charged lipid headgroups from either leaflet, ions and water molecules diffusing into the bilayer center, as well as the polar and positively charged groups in the CT when it is present. Because the limited space at the TMD core can only accommodate at most two inward-facing arginine residues, at least one R696 is forced outward to interact with lipid headgroups or CT residues. The differences in arginine orientation and interacting partners give rise to asymmetric protomer conformations and distinct TMD tilts (Figures 4A–F). When R696 points outward, its interactions with lipid head groups or CT residues can be strong enough to destabilize the local helix, introducing a kink into the TMD. Representative snapshots from different trajectories illustrate these asymmetric kinked conformations (Figure S8–S15). Beyond local deformation of the TMD, R696–lipid interactions perturb bilayer organization, inducing the translocation of lipid headgroups and water molecules toward the bilayer center (Figure 4G). We calculated the interaction frequencies of each TMD residue with lipid headgroups, lipid tails, and water molecules. The results show that the membrane disruption is persistent since many non-terminal TMD residues maintain frequent contacts with water and lipid headgroups throughout the simulation (Figure 4H).

Figure 4.

To explore whether R696 exhibits a preference for the exoplasmic versus cytoplasmic leaflet, or interacts with either randomly, we generated two initial structures (high and low) with the TMD positioned at two distinct positions, separated by 4 Å along the membrane normal (see Methods). In simulations initiated from the “high” TMD configuration, R696 residues in three protomers interacted with the lipid headgroups in either leaflet (Figures 4B, C, S8 and S10) In contrast, in simulations initiated from the “low” TMD configuration, R696 residues interacted exclusively with the headgroups in the cytoplasmic leaflet (Figures 4E, F, S9 and S11). In the full-length systems, the plate-shaped CT occupies substantial space in the cytoplasmic leaflet, displacing lipids during model construction. Because the CT is not thick enough to fully span the cytoplasmic leaflet, an empty gap remained between the CT and the exoplasmic leaflet (Figure S1D). Over time, lipids in the exoplasmic leaflet shifted downward while the CT residues moved upward to fill this space, causing local bilayer thinning. Under these conditions, upward-oriented R696 can still interact with the headgroups in the exoplasmic leaflet, but downward-oriented R696 primarily contacted CT residues, with rare cases of lipids migrating upward from the cytoplasmic leaflet and approaching the protein. A comparison of all full-length systems (Figures S12–S15) shows that R696 preferentially adopted downward orientations in the simulations initiated from the “low” TMD configuration (Figures S13 and S15). In these cases, upward-shifted CT residues interact with the C-terminal half of the TMD, rather than solely with R696, resulting in deeper burial of the TMD in the membrane.

MPER adopts diverse conformations, and its exposure depends on both MPER and TMD conformations

Starting from the initial helix-turn-helix conformation consisting of two separate helical segments, the N-terminal half (MPER-N) and the C-terminal half (MPER-C), MPER underwent rapid rearrangements, and a wide variety of conformations were sampled across all trajectories. In the initial structure, the trimeric MPER was positioned perpendicular to the membrane, with the helical MPER-N tilted inward and MPER-C tilted outward (Figure 4A). Such conformation and orientation were maintained in some trajectories such as CL^ΔCT3 (Figure S9C), but in others the helix-turn-helix MPER shifted into a horizontal orientation parallel to the membrane surface (Figure S13A) or a more vertical arrangement with both MPER-N and MPER-C tilted outward (Figures 4E and S111A). It was also observed that the HR2 helix in the ectodomain, MPER, and TMD merged into a continuous long helix (Figures 4C, F, and S10C). In addition, it was common for part of the MPER, particularly the MPER-C, to lose its helical structure and become a random coil. Therefore, the distinct MPER conformations reported in various experimental studies can be all valid as each captures a possible state within the highly flexible conformational landscape of the MPER. It is noteworthy that the MPER in three protomers can adopt different conformations and orientations, resulting in asymmetric local structures.

Because the MPER is a target of multiple bNAbs, we next examined how its exposure is affected by TMD conformation. As described above, interactions between R696 and lipid headgroups affect the burial depth of the TMD, thereby influencing the positioning of the adjacent MPER. To quantify MPER exposure, we measured the vertical distance from the Cɑ of F673, approximately the midpoint of MPER, to the highest point of the neighboring lipid headgroups, denoted by d_F673. Positive values indicate that F673 lies above the bilayer surface, and negative values indicate its membrane burial. In the initial ‘low’ and ‘high’ TMD configurations, d_F673 was 6.1 Å and 9.1 Å, respectively, but across simulations it spanned a wide range from -15 Å to 20 Å (Figures 5A and B). Two examples illustrate this variability. In the first example, all three R696 residues interacted with the cytoplasmic leaflet, drawing the TMD deeper into the membrane. Consequently, the entire MPER-C and most of MPER-N were buried in the membrane, with one F673 positioned 11.3 Å below the membrane surface, thus bringing the ectodomain in close proximity to the membrane (Figures 5C and E). In the second example, two R696 residues interacted with the exoplasmic leaflet, while one interacted with the cytoplasmic leaflet. In the protomer with upward-oriented R696, the MPER-N, MPER-C, and TMD formed a continuous helix. As a result, most of the MPER extended outside the bilayer, with d_F673 reaching 18.0 Å and the ectodomain displaced farther from the membrane (Figures 5D and F). The d_F673 distributions across all cleaved CT-truncated systems suggest that d_F673 tends to be smaller when the simulations started from the “low” TMD configuration (Figure 5B). In both “high” and “low” configurations, the mean d_F673 calculated from the simulation trajectories is smaller than its value calculated from the initial structure, due to both protein conformational change and lipid diffusion. In the “low” configuration, d_F673 decreased by 5.0 Å (from 6.1 Å to 1.1 Å), a slightly larger reduction than in the “high” configuration where d_F673 dropped by 3.5 Å (from 9.1 Å to 5.6 Å), indicating that deeper MPER burial arises from a larger-scale protein motion rather than merely from a lower initial placement.

Figure 5.

Ectodomain epitopes are conditionally accessible, whereas MPER epitopes are virtually inaccessible in the closed prefusion state

To access the accessibility of epitopes on different parts of the trimeric Env protein in the prefusion state, we quantified how frequent their epitopes were exposed without steric clashes from neighboring protein residues, glycans or membrane lipids to six selected antibodies: PGT128 targeting the V3 loop and the N332 glycan (PDB ID: 5JSA); PG9 targeting the V1/V2 loop (PDB ID: 3U2S); VRC01 targeting the CD4 binding site (PDB ID: 4LST); 35O22 targeting the gp120–gp41 interface (PDB ID: 4TVP); and 10E8 (PDB ID: 6VPX) and 4E10 (PDB ID: 1TZG) targeting the MPER (Figure 6A, Table S2-S5).

Figure 6.

The epitope of PGT128 consists of the V3 loop and N332 glycan, and it is widely recognized that PGT128 binding is mediated by the protein-protein interaction with V3 loop and the protein-glycan interaction with N332 glycan, facilitated by its extra-long HCDR3 loop penetrating the glycan shield (Figure S16). Due to the absence of defined secondary structures, glycans exhibit greater conformational flexibility than the protein, although the crowded surface of gp120–gp41 may partially restrict the motion of glycans. N332 glycan must adopt specific conformations to enable specific interactions with PGT128, whereas most conformations block antibody approach. Additional glycans on N137, N156, and N301 can further occlude the site (Figure 6B). Epitope accessibility to PGT128 varied widely across protomers and trajectories, with many cases exceeding 35%. In about half trajectories, at least one protomer exhibited >35% accessibility (Table S2). For instance, in system CH^ΔCT1, two protomers showed accessibility >40%, whereas in system CL^CT2, all three protomers are <5%. The second antibody, PG9, targets the V1/V2 apex, where binding can be hindered by six glycans, N156, N160, and N185E on the same protomer, and N160, N185E, and N185H on the neighboring protomer (Figure S17). In most trajectories, the epitopes on all three protomers were occluded in >95% of the snapshots, with the exception of a few systems (CH^ΔCT1, CH^ΔCT2, and CL^ΔCT3) where accessibility was non-negligible. (Table S3). The third antibody, VRC01, targets the CD4 binding site where six glycans, N185H, N197, N276, N363, and N462 on the same protomer and N301 on the neighboring protomer, can interfere with antibody binding (Figure S18). Similar to PGT128, the VRC01 epitope is moderately to highly accessible to at least one protomer in many trajectories, whereas in some trajectories, the epitopes on all three protomers are nearly completely shielded (Table S4).

35O22 targets the gp120–gp41 interface that is not intrinsically membrane-proximal. However, its binding orientation combined with ectodomain tilting can lead to clashes of 35O22 with the membrane (Figure S19). Comparing accessibility with and without including the steric effects of membrane lipids shows that the membrane can hinder 35O22 binding, particularly when it approaches from the direction of ectodomain tilt (Figure 6B). For example, in system UH^ΔCT2, the accessibility frequencies of the three protomers were 49%, 5%, and 25%, when only the shielding of three glycans (N88, N234, and N618) was considered. When the steric effects of both glycans and membrane lipids were included, the first protomer’s frequency dropped from 49% to 1%, while the other two remain essentially unchanged (Table S5).

The last two antibodies, 10E8 and 4E10, target the MPER. In the PDB structures of 10E8, the epitope is a 17-residue peptide corresponding to residues 671-687 in this work. In the PDB structure of 4E10, the epitope is a 12-residue peptide corresponding to residues 669-680. Superposition of these antibodies onto the initial simulation structure revealed extensive clashes with the protein and glycans linked to N611, N618, and N637 on the neighboring protomer, as well as slight steric clashes with the membrane (Figures S20 and S21). Using the strict criteria (≤10 heavy-atom clashes with protein/glycan and ≤20 with lipids), no snapshot was identified in which the epitope of either antibody was accessible. Even with relaxed criteria (≤20 heavy-atom clashes with protein/glycan or ≤40 with lipids), accessible cases remain rare, which require either large ectodomain tilts to create space on the opposite side or the MPER transitioning into continuous helices to displace the ectodomain from the membrane (Figures 6C, S22–S25). Given the extremely low frequency of such specific conformations, MPER epitopes are effectively inaccessible in the prefusion trimer, suggesting that MPER-targeting antibodies such as 10E8 and 4E10 act at later stages of viral entry.

Discussion

The hydrophobic and flexible nature of the MPER has made it difficult to crystallize, and most structural insights have come from NMR studies. However, structural information on the MPER and TMD has been inconsistent in the literature. For example, one study of the MPER in DPC micelles suggested that the MPER adopted a distorted helical structure lying parallel to the membrane surface,¹⁶ while NMR studies using gp41 fragments containing part of the MPER and TMD embedded in bicelles or bilayers reached conflicting conclusions. Some proposed a kink¹⁴ or turn¹¹ between the MPER and TMD, whereas others suggested that the MPER and TMD formed a continuous helix aligned perpendicular to the membrane.^10,17 In the NMR structure (PDB ID: 7LOI) used to build our model structure, a kink appears between the MPER and TMD, along with a turn in the middle of the MPER.^13,23 Crystal structures of MPER peptides bound to various antibodies show that in addition to the helix-turn-helix conformation, MPER can also adopt alternative structures in which the MPER-C remains helical while part of the MPER-N becomes unstructured.²⁴ Such variability underscores the intrinsic structural flexibility of the MPER, with each experimental method capturing only a subset of possible states under specific conditions. In this study, we employed MD simulations to sample a broader conformational landscape and gain a more comprehensive view of the structural heterogeneity of the MPER. The results reveal that the MPER can adopt all of the conformations reported experimentally. While MPER plasticity has been linked to its role in virus-host membrane fusion because it enables the ectodomain and TMD to adopt distinct orientations during large-scale structural rearrangements, our results show that this flexibility is already present in the prefusion state. At this stage, the MPER functions as a flexible hinge that facilitates ectodomain tilting, which can be critical for the spatial alignment of the CD4-binding site with the host receptor for efficient receptor engagement.

To probe antibody recognition, we assessed epitope accessibility using snapshots extracted from simulation trajectories. Unlike static experimental structures, this approach incorporates protein dynamics as well as the influence of glycans and membrane lipids, yielding a frequency-based measure indicating how often each epitope is accessible to its corresponding antibody, rather than a binary accessible/occluded classification. It is worth noting that viral glycans are highly variable in both site occupancy and the specific glycoform present at each glycosylation site, as indicated by mass spectrometry data.^25,26 For structure modeling, however, we assumed full occupancy at all glycosylation sites and selected a single representative glycoform per site. As a result, the estimated glycan shield may not fully reflect the in vivo situation, and the calculated frequency of antibody accessibility should be regarded as an approximation. Our results indicate that epitopes on the ectodomain, though heavily shielded by glycans, can still become transiently accessible in the closed prefusion state. The frequency provides a quantitative measurement of how vulnerable each epitope is to antibody binding. When estimating steric clashes between the membrane and antibodies, we adopted relaxed cutoffs to account for lipid fluidity, but the membrane due to fast lipid dynamics may in reality accommodate even greater steric overlap. Although the MPER can transiently extend out of the membrane when gp120–gp41 remains in the prefusion state, MPER-targeting antibodies still encounter substantial steric hindrance from the bulky gp120 subunit, surrounding glycans, and the membrane throughout nearly all simulations. Therefore, these antibodies are unlikely to adopt the proper orientations required to bind their epitopes. This is consistent with experiment studies indicating that MPER-targeting antibodies bind effectively only after the gp120–gp41 trimer undergoes major conformational rearrangements toward a fusion-intermediate or post-fusion state.^27–29

Overall, the data presented here demonstrate that structural modeling integrated with molecular dynamic simulations can be applied to complex biomolecular systems, such as the gp120–gp41 trimer embedded in a membrane bilayer, which enables detailed characterization of protein dynamics at the molecular level and facilitates the investigation of antigen-antibody interactions, thereby offering potential guidance for rational vaccine design.

Methods

System building

Modeling of full-length gp120–gp41 trimer

The full-length gp120–gp41 trimeric model was generated by combining the crystal structure of the ectodomain with the NMR structure containing the MPER, TMD, and CT (Figure 1A and S1A). Considering structure resolution and the number of missing residues, we selected the crystal structure of a soluble Env trimer exhibiting a closed-from, native-like prefusion conformation (PDB ID: 6B0N),³⁰ which contains the entire ectodomain including the V1-V5 loops, fusion peptide (FP), heptad repeat 1 (HR1), and heptad repeat 2 (HR2) (Figure 1B). In contrast to the ectodomain well represented in the PDB with numerous crystal and cryo-EM structures, structural data for the TMD are scarce. For this region, we chose the NMR structure (PDB ID: 7LOI)¹³ containing the MPER, TMD, and CT. Both 6B0N and 7LOI are trimeric structures and the C-terminus of 6B0N includes five residues (L660–D664) overlapping with the N-terminus of 7LOI. Measurement of the Cα–Cα distances for D664 in the three protomers revealed that the inter-protomer spacing in 7LOI (16 Å) was much shorter than in 6B0N (39 Å) (Figure S26A). To reconcile this difference, we performed short MD simulations on 6B0N with residues E32–Q640 fixed and a distance restraint applied to D664 in order to reduce the Cα– Cα distance to 33 Å. Similarly, for 7LOI, residues I675–L856 were fixed, and a distance restraint was applied to D664 to increase the distance to 33 Å (Figure S26B). The two adjusted structures were then combined by aligning the overlapping residues. The missing loops in the ectodomain (T63, D149, E185A–N185I, S401–G409) and the unresolved region between the TM and CT (F717–G738) were grafted from a modeled structure generated by I-TASSER³¹ (Figure S26C). 6B0N is a native flexibly linked (NFL) trimer in which the furin cleavage site ⁵⁰⁸REKR⁵¹¹ is substituted with a 10-residue linker ⁵⁰⁸GGGGSGGGGS⁵¹¹, but the underlined residues are unresolved in the structure. To construct the cleaved gp120–gp41 model, we removed the initial “GG” and final “S” from 6B0N. For the uncleaved gp160 model, we closed the cleaved loop by adding ⁵⁰⁸REKR⁵¹¹ back to the structure. If we simply connect G507 and A512 without adjusting the conformation of flanking loops on both sides of the cleavage site, two neighboring protomers become entangled, yielding a knot-like fold (Figures S27A and B). In 6B0N, the helical HR2 region lies at the C-terminus of each protomer, allowing it to thread through a loop formed by the neighboring protomer (Figure S27C). However, in the context of the full-length gp120–gp41, it is implausible that the three protomers are entangled as it does not conform to a realistic folding pathway. Therefore, we performed short MD simulations with restraints to adjust the conformations of the flanking loops and the HR2 helix before closing the loop with REKR fragment (Figures S27D and E). Two mutations (S764C and S837C) were applied to CT to reintroduce the palmitoylation sites, and lipid tails oriented towards the hydrophobic core of the bilayer were added to the palmitoylation sites (Figure S1D). In addition to the full-length model, we built the CT-truncated model by removing the modeled residues F717–G738 that are originally unresolved in the PDB structure 7LOI, and the residues E739– L856 that form the large CT plate.

Glycosylation

N-linked glycans were modeled using Glycan Reader & Modeler^32–34 in CHARMM-GUI.³⁵ The reported mass spectrometry (MS) experiments have revealed the probability of different types of N-linked glycans at each glycosylation site.^25,26 For each site, we selected the one with the highest probability. As the MS data provide only the glycan type (high-mannose, hybrid, and complex) and composition (number of each monosaccharide unit), a representative isomer was chosen when multiple isomers corresponded to the same composition. For example, HexNAc(2)Hex(9), i.e., the Man₉ glycan, has a single isomer, whereas HexNAc(2)Hex(8), i.e., the Man₈ glycan, has multiple isomers depending on which non-reducing terminal mannose is trimmed from Man₉. For complex N-linked glycans, one composition can correspond to multiple isomers differing in the number of branches on the α1-3 and α1-6 arms and in the placement of non-reducing terminal neuraminic acid (Neu5Ac). In the absence of linkage-specific information, one isomer was arbitrarily selected in cases of multiple possibilities. The selected glycan sequences for the 27 glycosylation sites in each promoter are summarized in Table S6.

Membrane

The full-length and CT-truncated gp120–gp41 models were embedded into an asymmetric lipid bilayer with the lipid composition corresponding to a mammalian plasma membrane,²⁰ which consists of phosphatidylcholine (PC), phosphatidylethanolamine (PE), phosphatidylinositol (PI), phosphatidylserine (PS), phosphatidic acid (PA), sphingomyelin (SM), cholesterol (CHOL), and glucosylceramide (GlcCer). In the CT-truncated system, the exoplasmic and cytoplasmic leaflets contain similar numbers of lipids, although their compositions differ. In the full-length system, the cytoplasmic leaflet contains approximately 100 fewer lipids than the exoplasmic leaflet due to the space occupied by the CT. The lipid composition is summarized in Table S7. Although it is recognized that the charged residues R683, R707, and R709, flanking the N- and C-termini of the TMD, are anchored within the lipid headgroups, there remains some flexibility to shift the TMD slightly up or down within the membrane. Assuming the lipid bilayer is aligned parallel to the xy-plane, the protein was positioned at two distinct depths along the z-axis, differing by approximately 4 Å, referred to as the “high” and “low” TMD configurations. To examine whether this subtle difference in the initial configurations leads to distinct conformational changes during simulations, we initiated independent simulations from each configuration. The glycoprotein-membrane system was solvated in a box of approximately 210 × 210 × 260 ų using the TIP3P water model,³⁶ and KCl was added at a concentration of 0.15 M to neutralize the system. Final simulation input files were generated using Membrane Builder^37–41 in CHARMM-GUI.

Simulation Details

The combination of cleaved vs. uncleaved, full-length vs. CT-truncated, and high vs. low TMD positions in the membrane results in eight distinct configurations, and we performed three independent 1-μs all-atom MD simulations for each configuration. The CHARMM36(m) force field^42–46 was used for proteins, carbohydrates, and lipids. The total number of atoms is approximately 1,100,000 (including ∼300,000 water molecules and ∼1,700 ions), with slight variations between systems. van der Waals interactions were smoothly switched off over 10–12 Å using a force-based switching function,⁴⁷ and long-range electrostatic interactions were calculated using the particle-mesh Ewald method⁴⁸ with a mesh size of ∼1 Å. All simulations were performed using GROMACS.⁴⁹ Bond lengths and angles involving hydrogens were constrained using the LINCS algorithm.⁵⁰ The system was equilibrated in the canonical (NVT) ensemble at 310.15 K for 2 × 1.25 ns with a 1-fs time step using the Berendsen thermostat⁵¹ (coupling constant τ_t = 1 ps), and then switched to the isothermal–isobaric (NPT) ensemble for 1.25 ns with a 1-fs time step, followed by 3 × 5 ns with a 2-fs time step using the Berendsen thermostat (τ_t = 1 ps) and barostat (reference pressure = 1 bar, coupling constant τ_p = 5 ps, compressibility = 4.5×10⁻⁵ bar⁻¹). Positional and dihedral restraints were applied to proteins, glycans, and lipids, with force constants progressively reduced over successive intervals. In the production run, temperature was maintained using the Nosé–Hoover thermostat^52,53 and the pressure coupling was applied using the semi-isotropic Parrinello–Rahman barostat.^54,55 A 4-fs time step was used with the hydrogen mass repartitioning technique⁵⁶ and all restraint potentials were removed. The Python library MDTraj⁵⁷ was used to analyze the simulation trajectories.

Assessment of Antibody Epitope Accessibility

We selected six antibodies based on epitope location and the availability of PDB structures containing both antibody and epitope: PGT128 targeting the V3 loop and the N332 glycan (PDB ID: 5JSA);³ PG9 targeting the V1/V2 loop (PDB ID: 3U2S);⁵⁸ VRC01 targeting the CD4 binding site (PDB ID: 4LST);⁵⁹ 35O22 targeting the gp120–gp41 interface (PDB ID: 4TVP);⁶⁰ 10E8 (PDB ID: 6VPX)⁶¹ and 4E10 (PDB ID: 1TZG)⁶ targeting the MPER. For each antibody– epitope complex, TM-align⁶² was used to identify the optimal alignment between the epitope from the PDB structure and each protomer of the trimeric protein, and the resulting rotation matrix was applied to place the antibody relative to each protomer. Steric clashes were then assessed by counting any heavy atom of protein, glycan, or lipid within 2 Å of a heavy atom of the antibody. An epitope was considered occluded by protein and glycans if more than 10 heavy-atom clashes occurred with these components. To account for membrane flexibility, we adopted a more permissive cutoff for lipids: the epitope was classified as occluded by the membrane if more than 20 lipid heavy-atom clashes were detected. To calculate the frequency of epitope accessibility, snapshots were extracted from simulation trajectories, and structural alignment and assessment of steric clashes were performed on each snapshot.

Data availability

The modeled structures, force field parameters, and GROMACS input files are available on Zenodo [https://doi.org/10.5281/zenodo.17259909]. The code used in this article can be found at https://doi.org/10.5281/zenodo.17259951.

Acknowledgements

This study was supported by NIH R35 GM153458. We thank Dr. Thomas Walz and Giorgos Hiotis for invaluable discussions on our work.

Additional information

Contributions

Y.C. and W.I. conceived the concepts. Y.C. performed the simulations. Y.C. conducted the data analysis. Y.C. and W.I. wrote this paper. Y.C. and W.I. discussed the results and commented on the manuscript at all stages.

Funding

HHS | National Institutes of Health (NIH) (R35-GM153458)

• Wonpil Im

Additional files

Supporting Tables and Figures