Introduction

Biological electron transport chains (ETC) comprise a series of membrane-bound enzyme complexes (CI-CIV) that transfer electrons towards oxygen and protons across a biological membrane, creating a proton motive force (PMF) that powers the synthesis of ATP and active transport (1, 2). Biological membranes are often envisaged in the light of the classical fluid mosaic model (3), where the membrane proteins and lipids independently diffuse as a 2-dimensional solution. However, the biological membranes are highly crowded (4-6), particularly the inner mitochondrial membrane (IMM), with a protein content of around 45% (7). In this regard, the majority of the respiratory complexes do not diffuse independently within the membrane but operate rather as larger supercomplex (SC) assemblies, as first revealed by blue native gels (8) and confirmed by structural analyses (9-11). Moreover, several lipid molecules, particularly cardiolipin, a central anionic lipid of the IMM, are often tightly bound to the SCs (12-15). Yet, despite these structural insights, the functional role of the SC assemblies and the physical principles leading to their formation remains elusive and much debated.

SCs have been suggested to provide kinetic advantages (12, 16), favor substrate channeling ((10, 17), but cf. (18)), decrease formation of reactive oxygen species (ROS) (19, 20), stabilize individual proteins (21-24), and/or reduce non-specific protein-protein interactions (18). It has also been recognized that the high protein density in membranes may play an important role in the SC formation (18), and influence, e.g., the lipid curvature (7, 25). On a physiological level, decreased SC formation has been linked to various patho-physiological conditions such as diabetes (26), heart failure (27), and apoptosis (28), although recent experiments (29) on mice unable to form SC showed no significant differences relative to WT mice under the studied conditions.

In recent years, cryo-electron microscopic (cryo-EM) and cryo-electron tomographic (cryo-TEM) studies revealed the structure of several SCs, including the mammalian 1.5 MDa SCI/III₂ complex and the respirasome 1.7 MDa SCI/III₂/IV complex (30-34), as well as various exotic megacomplexes isolated from single-celled eukaryotes (35, 36). SCs are also essential for some bacteria, such as the SCIII₂/IV₂ of actinobacteria (15, 37-39) that catalyze quinol oxidation coupled to oxygen reduction in a tight obligate protein assembly. Inhibition of such bacterial SC provides potential avenues for treating pathogenic infections, such as tuberculosis, with emerging multi-drug resistant strains.

To probe the physical principles leading to the SC formation and its molecular consequences, we study here the structure and dynamics of the mammalian mitochondrial SCI/III₂ by combining large-scale atomistic and coarse-grained molecular dynamics simulations with analysis of cryo-EM data and mathematical models that provide insight into crowding effects, “phase transitions” underlying the SC formation, as well as the dynamic changes in the SCs.

The studied SCI/III₂ catalyzes an NADH-driven quinone reduction and offers a unique system to probe how the quinone (Q) substrate is transported within the SC assembly, and how the individual proteins interact within the membrane. In this regard, the quinol (QH₂) produced by the CI module of the SC is re-oxidized by the CIII module, whilst the proton transport activity (CI: 2H⁺/e⁻; CIII₂: 1H⁺/e⁻) powers the synthesis of ATP (Fig. 1). Interestingly, prior cryo-EM studies of this SC (40) revealed an asymmetric protein assembly with a trapped quinone in the proximal/distal Q_o site, suggesting possible functional consequences. However, the dynamics and interaction of the SC and its substrates within the surrounding membrane environment and the individual active site remain poorly understood.

Figure 1.

Results

Protein-membrane interactions modulate the lipid dynamics and membrane strain

To probe how the SC formation affects the dynamics of the OXPHOS proteins, the lipid membrane, and the quinone pool, we studied the ovine SCI/III₂ by both atomistic (aMD) and coarse-grained molecular dynamics (cgMD) simulations. While the aMD simulations give insight into the microsecond dynamics (∼11 μs in total) of the system with atomic details, our cgMD simulations allowed us to probe the SC dynamics and its interactions in the membrane on much longer timescales (∼0.3 ms) at a more approximate level (see Methods). The SC models were further compared against simulations of the individual CI and CIII₂ complexes, as well as simulations of a lipid membrane, with a POPC:POPE:CDL ratio (2:2:1) as in the IMM.

We find that the individual OXPHOS complexes, CI and CIII₂, induce prominent molecular strain in the lipid membrane that decreases the local membrane thickness at the protein-lipid interface (Fig. 2G, Fig S2A,D,E). Remarkably, the SC assembly partially releases this strain energy as a result of a smaller solvation area established at the SC-lipid interface (Fig. 2G, Fig. S5, Fig. S7, Table 4). This effect is likely to arise from the thinner hydrophobic belt region of the OXPHOS protein (ca. 30 Å, Fig. S1A) relative to the lipid membrane (40.5 Å, Fig. S1B), and leads to a ∼30% accumulation of cardiolipin (CDL) around the SC (Fig. 2H, Fig. S2B) that thermodynamically prefers thinner membranes (∼38 Å) (Fig S1B). The membrane perturbation is supported by both our aMD (Fig. S2A) and cgMD simulations (Fig. 2G), suggesting that the overall results are robust. We find that an increase in the lipid tail length decreases the relative stability of the SC (Fig. S5C), further supporting that the hydrophobic mismatch between the OXPHOS proteins and the lipid membrane modulates the strain effect (Fig. S5A-C). In addition to the perturbations in the membrane thickness, the OXPHOS proteins also distort the local membrane plane (Fig. 2A-F, Fig. S4A-C, Fig. S7) that could modulate the accessibility of the Q/QH₂ substrate into the active sites of CI and CIII₂ (see below, section Discussion).

Figure 2.

Although the present cryo-EM maps do not resolve the molecular details of the surrounding lipid membrane (but see below), we performed a spatial integration of the experimental cryo-EM density map by training a neural network model (see Methods) that allowed us to identify and calculate the membrane thickness of the lipid belt around Complex I (PDB ID: 6RFR, EMD-4873) (41). In this regard, we find that the experimental data also shows a statistically significant membrane thinning around CI relative to the SC at the protein-membrane interface (Fig. 2J-L), strongly supporting our simulation results. We note that during the finalization of this work, a membrane distortion effect was also reported for in situ high-resolution cryo-EM structures of respiratory SCs (42), thus providing additional support for our findings.

We observe that the SC assembly influences the lipid and water dynamics at the protein-lipid interface and perturbs the local dielectric at the membrane plane (Fig. S1C). The dielectric constant (ε⊥) shows a local increase up to 1.5 nm from the membrane before it decays to the aqueous bulk dielectric (Fig. S1C). We note that this effect is strongly affected by the lipid type (Fig. S1C), particularly by the CDL that accumulates around the SC (Fig. S1D) and creates a micro-environment, which could enhance a local proton gradient.

Taken together, our combined findings suggest that the SC formation is driven by entropic effects that reduce the molecular strain in the lipid membrane, whilst the perturbed micro-environment also affects the lipid and Q dynamics, as well as the dynamics of the OXPHOS proteins (see below).

The SC assembly alters the quinone dynamics but does not support substrate channeling

In addition to the phospholipids, the IMM contains 3-5% ubiquinone (Q), which functions as the electron carrier between CI and CIII₂ in the form of ubiquinol (QH₂). Our simulations, performed with ca. 5% Q/QH₂ (1:1) suggest that the long isoprenoid tail of ubiquinone (Q₁₀) has similar physico-chemical properties as the lipid tails, whilst the Q/QH₂ headgroup is polar and localizes at the membrane surface (Fig. 2), with a flip-flop rate between the leaflets of ca. 100-150 ns (Fig. S3A,B). The fast flip-flop motion could support the quinol diffusion from CI to CIII₂, in which the quinol exits the negatively charged (N-side / matrix side) membrane surface of CI and enters the Q_o site of CIII₂ near the positively charged (P-side / intermembrane space) of the membrane (see Discussion). Interestingly, the Q pool does not accumulate uniformly around the SC. Instead, we observe a ca. 3% increase in the local Q concentration near the substrate binding sites of both CI and CIII₂ relative to the bulk membrane, whilst the immediate OXPHOS surroundings show an overall Q/QH₂ depletion (Fig. 2I, Fig S2C, Fig S6E). This local Q pool could arise from specific interactions between Q and the OXPHOS proteins, e.g., in subunits ND1 and NDUFA9 of CI that contain several non-polar residues and surface charges that interact with the quinones (Fig. 2M, Fig. S13). In contrast, on the proximal side of CIII₂, we observe a subtle decrease in the local Q concentration relative to CIII₂ alone (Fig 2I, Fig S2C, Fig S6D,E), an effect that could arise from a shift in the membrane plane near the Q_o site and affect the substrate access to the proximal Q_o site (Fig S4, S7). These observations are consistent with the occupied distal Q_o site and empty proximal Q_o site observed in the cryo-EM structure of the SCI/III₂ (40).

Overall, while our simulations indicate that quinones accumulate around the OXPHOS proteins, our data do not support substrate channeling (10) (see Discussion), as we observe neither a contiguous region of an elevated Q pool between the complexes, nor a directed diffusion pathway of the quinol from CI to CIII₂ (Fig. 2I, Fig. S3D-I).

Moreover, the locally hindered access to the proximal Q_o site due to the membrane shift could further hamper the substrate binding into the active site.

The SC assembly alters the conformational dynamics of CI and CIII₂

We next probed how the SC assembly affects the dynamics of the individual OXPHOS proteins based on our atomistic MD simulations of CI, CIII₂, and the SCI/III₂. Despite the large (1.5 MDa) SC, we note the SC assembly is stabilized by only a few specific hydrogen-bonding and charged interactions at the interface of NDUFB9/NDUFB4 and UQCRC1, particularly the Arg29^FB4-Asp260^C1/Glu259^C1 ion-paired network that interacts by ca. -10 kcal mol^-1 during the MD simulations (Fig. 3A,C). Indeed deletion of this Glu258-Asp260 region in UQCRC1 led to a drastic reduction in the stability of SC in mice (29) (see Fig. S11), thus supporting the importance of these interactions in the SC.

Figure 3.

Based on essential dynamics analysis (EDA, see Methods, Movies S1-S5) that projects out dominant protein motion, we find that the CI and CIII₂ modules of the SC undergo a back-and-forth rocking motion around the interface region, leading to a breathing motion between the hydrophilic domains of the OXPHOS complexes on the N-side of the membrane (matrix side) (Fig. 3B, Mode 1, Movie S1), while the second dominant motion (Mode 2, Movie S1) couples the opening/closing motion of CI (cf. also ref (43)) with a rocking motion of CIII₂. The third mode arises from a combination of the twisting motion of the hydrophilic domain of CI and a minor rocking motion of CIII₂ (Mode 3, Movie S1), which could influence the membrane accessibility of the Q sites. The EDA thus suggests that the dynamics of CI and CIII₂ are allosterically coupled.

We note that the SC assembly induces a subtle conformational change in CI that increases the angle between the hydrophilic and membrane domains (Fig. 3F). This motion is linked to the active/deactive-transition (A/D) (44), which regulates CI activity, particularly under hypoxic / anoxic conditions and reverse electron transfer (RET) conditions (but cf. also ref. (45)). These changes in CI lead to a subtle shift in the conformation of key transmembrane helices (TM) and surrounding conserved loops. In the SCI/III₂, the TM3 of subunit ND6 undergoes a twist towards the deactive state (Fig. 3H, Fig. S12A-D), while the β1-β2 loop of subunit NDUFS2 exhibits lower flexibility in both the apo and the QH₂ bound state relative to CI alone (Fig. S12E-H). These regions are central in modulating the coupling between the proton transport and the electron transfer activities ((46), cf. also Ref. (47)), and could thus have functional consequences. CI also shows a dominant vibrational motion that twists the hydrophilic domain relative to the membrane domain (Fig. 3B,D), and samples a wider angle between hydrophilic and membrane domains in the SC, (Fig. 3F) that result in a larger twisting angle relative to the isolated CI (Fig. 3F). These alternations suggest that the opening/closing and twisting motions could be sterically hindered by the adjacent CIII₂ (see Movies S1-S5), and that the active/deactive transition is coupled to the global motions of the SC. Consistent with our previous findings (43), CI shows also another dominant bending motion, where the membrane domain rotates relative to the hydrophilic domain within the membrane plane.

CIII₂ also undergoes dynamical changes upon the SC formation, particularly by enhancing the normal modes connected to the opening/closing transitions around the iron-sulfur (FeS) Rieske center (Fig. 3B,E) that could affect the electron transfer activity of CIII₂ (48). In the SCI/III₂, QH₂ binds further away from the FeS center in the proximal Q_o site, while in the distal site, the FeS center moves closer to the heme c₁ (Fig 3G,I,J). This asymmetry suggests that the global motion of the SC could regulate the ‘preferred’ Q_o site for the electron bifurcation process and possibly favor the electron transfer onwards to the Complex IV (CIV) that resides on the distal side of the CIII₂ module in the respirasome (SCI/III₂/IV) (30-32). These findings further support a complex allosteric crosstalk within the SC, as also suggested by the recent in situ cryo-EM study of mitochondrial SCs (40). In this regard, we find that the ligand state of CI (apo or QH₂) affects the conformational dynamics and the interaction interface of the SCI/CIII₂ (Fig. S11A-E). This surprising long-range effect is likely to result from the increased flexibility of the apo state (Fig. S11G) that, in turn, modulates the interaction at the interface of the SC. Interestingly, similar conformational changes affecting both the CI and CIII₂ domains of the SC are also supported by the recent in situ cryo-EM structure of mitochondrial SCs (42).

Enthalpy-entropy compensation drives SC formation

Our molecular simulations suggest that while the membrane strain provides an entropic driving force for the SC formation, the molecular interactions at the interface of the assembly are essential for enthalpically stabilizing the SC over non-specific protein assemblies (Fig. S11). To assess how the SC formation is affected by the strain effects, protein-protein interactions, and the protein/lipid ratio, we developed a simple statistical-mechanical lattice model of the CI and CIII₂ diffusion in the IMM (Fig. 4A). To this end, the protein-protein interactions were described by tunable interaction energies (E_specific and E_non-specific), modeling the specific hydrogen-bond / salt-bridges interactions and the non-specific interaction, while the membrane-protein interactions were tuned by the “entropic” strain energy (E_strain). Our model suggests that the ratio between the strain energy and specific interactions indeed provides the key driving force for the SC formation, with a small (30%) decrease in E_strain leading to a ca. 60% decrease in the SC population (Fig. 4B). Our model further predicts that the SC population drastically increases at a specific interaction energy threshold (E_specific<-2 k_BT) relative to the strain contribution (Fig. 4B), whilst the number of non-specific assemblies (adjacent proteins in non-SC orientations) relative to SCs is determined by the ratio of E_specific to E_non-specific. With a high membrane strain, we find that an increase in temperature also favors the formation of SC assemblies, as it is entropically favored to reduce the overall number of ordered (“strained”) lipid molecules around the OXPHOS proteins. Similarly, a high protein-to-lipid ratio significantly increases the population of SCs, suggesting that a high protein packing in IMMs favors the SC formation, whilst the strength of specific interactions determines the relative population of non-specific assemblies in crowded environments. Taken together, our lattice model, despite its simplicity, captures and validates key features observed in our molecular simulations and supports that the SC formation is strongly affected by entropic forces.

Figure 4.

Discussion

We have shown here that the respiratory chain complexes perturb the IMM by affecting the local membrane dynamics. In this regard, we suggest that the SC assemblies form by condensation of local high-entropic membrane regions around the OXPHOS proteins (Fig. 5A). However, as the entropic effects also favor the formation of non-specific protein assemblies, unique interactions, such as the ion-paired network around UQCRC1 and NDUFB9/NDUFB4 (Fig. 3A), must enthalpically stabilize the SC assemblies. The suggested principles show similarities to the hydrophobic effect driving protein folding by condensation of locally ordered water clusters around unfolded protein patches (49).

Figure 5.

We further probed the thermodynamic effects underlying the SC formation by developing a 2D lattice model. Despite its simplicity, the model supports that the SC stability is determined by a delicate balance between (entropic) membrane strain effects and specific (enthalpic) protein-protein interactions, but also strongly affected by the protein concentration and temperature effects (Fig. 4B). Moreno-Loshuertos et al. (50) recently suggested that elevated temperatures (>43^°C) may indeed disrupt SCs, although the stability of the individual OXPHOS proteins was also decreased in the studied conditions. As the enthalpic effects are of electrostatic origin, the SC stability could also be sensitive to the ionic strength and the PMF, which in turn depends on the metabolic state of the mitochondria.

At the molecular level, the SC formation leads to an accumulation of CDL at the protein-membrane interface, as well as a local Q / QH₂ pool near the substrate channel of CI and CIII₂ (Fig. 2M). We find that CDL prefers thinner membranes relative to the neutral phospholipids (PE/PC), and could thus partially compensate for the hydrophobic mismatch between the OXPHOS proteins and the membrane (Fig. S1). The Q diffusion is affected by both specific interactions with OXPHOS proteins, as well as the local membrane thickness (Fig. 2I). The CDL around the SC is thus likely to have both structural and functional consequences, consistent with its effect on the activity and dynamics of several membrane proteins ((51-54), cf. also Ref. (43)). In this regard, CDL was suggested to enhance the substrate dynamics within the Q-tunnel of CI that requires a twisting-bending motion around the membrane and hydrophilic domains (43), whilst destabilization of CDL-binding sites has indeed been shown to disrupt SCs, e.g., the III₂IV₁ SC in yeast (12). Moreover, defects in the CDL synthesis, e.g., in the Barth syndrome, (55), result in the disassembly of SC, indirectly supporting the involvement of CDL as a “SC glue”.

In addition to the changes in the membrane properties, we observed that the SC formation modulates the conformational dynamics of the individual OXPHOS proteins, especially the large-scale bending-twisting motion of CI, the conformation of individual TM helices and conserved loops around proton channels in CI, as well as the motion of the Rieske domain in CIII₂ (Fig. 3, Fig. S12). The dynamics of these regions are likely to module the activity of the OXPHOS proteins, e.g., the active/deactive transition of CI (56, 57) that regulate the ΔpH-driven quinol oxidation and reverse electron transfer (58, 59). Indeed, blocking loop motions surrounding these regions inhibits the proton pumping activity of CI (60), whilst the perturbed motion of the Rieske domain could modulate the electron bifurcation in CIII₂ and subsequent electron transfer to CIV. In this regard, we suggest the preferred QH₂ binding in the distal Q_o site has functional implications for the respirasome (SCI/III₂/IV), where the CIV module is located on the distal side of the CIII₂ protomer. The changes in the conformational dynamics upon SC formation may thus affect ROS production (17, 61) via the active/deactive transition of CI (56), although it should be noted that no differences in ROS generation were observed for mice unable to form SCs relative to WT mice (with ca. 75% SCs) under normal conditions (29). It is possible that differences occur only under more strained conditions, e.g., in hypoxia or together with disease-related mutations in the OXPHOS proteins, where the SCs could become functionally more important (see below).

To understand how SCs could influence the charge currents in the IMMs, we note that the local increase of the quinol concentration near the active site of CI, and the local decreased quinol pool around the proximal Q_o site of CIII₂ could create a substrate gradient (∇c) and affect the quinol flux between CI and CIII₂, J_CI→CIII = -D_mem∇c (D_mem – Q/QH₂ diffusion constant) – if the quinol concentration is rate-limiting for the function of CI or CIII₂ (but see below). The 2D-diffusion time for the quinol, τ = <r²>/4D_mem, between the OXPHOS complexes depends on the effective protein-protein distance (r), which can be estimated from the protein packing density (cf. (7, 34, 62-66)). Using the experimental protein copy numbers in IMM (cf. Refs. (7, 34, 62-66) and Extended SI Text, Table S5), we obtain average edge-to-edge distance between CI and CIII₂ of around 12 nm (Fig. 4E, cf. also (7)), which can be compared to the edge-to-edge distance of ca. 7 nm / 12 nm between the Q tunnel of CI and the proximal / distal Q_o sites of CIII₂ within the SC. The reduced diffusion distance could thus provide a subtle rate enhancement (τ_SC(r=7 nm)/τ_nonSC(r=12 nm)∼0.3) in favor of the SCs under substrate-limited conditions (Fig. 4E-F). We expect that the small variations in the lateral quinol diffusion (0.4 nm ns^-1 near CI, 0.8 nm ns^-1 in bulk) around the SC (Fig S3E-I) could favor the diffusion along the membrane patches, although this diffusion is also affected by non-specific collisions. However, due to the overall 3-5% Q / QH₂ concentration in IMMs, we estimate that each OXPHOS protein is surrounded by around 6 quinone / quinol species depending on the respiratory conditions (Fig. 5, Extended SI Text), and leading to a small (0.3 nm) nearest neighbor distance between Q and CI/CIII₂ (Fig. 4E, see SI Appendix). This implies that each OXPHOS protein has a saturated substrate Q pool within its “reaction sphere” upon high charge flux conditions (here 50% reduced Q pool). In other words, the CIII₂ node of the SC does not need to rely on the quinol generated at the CI node, but can instead oxidize a quinol molecule from the surrounding Q pool. In this regard, Hirst and co-workers (67) found that the alternative oxidase (AOX) co-reconstituted into liposomes outcompetes the quinol re-oxidation rate of the SCI/III₂, suggesting that the diffusion of the Q-species between the OXPHOS complexes is not rate-limiting under the studied conditions. However, at low QH₂/Q ratios (<5% in the model), the minimal nearest neighbor distance between the Q pool and CI/CIII₂ drastically increases (Fig. 4D-F) that in turn leads to a kinetic preference for the SC (Fig. 4D).

At the microscopic level, the SC influences the lipid and water dynamics at the protein-lipid interface and the physico-chemical properties of the IMM, including its local dielectric properties of the membrane. It is possible that such local micro-environments have functional implications that could affect the proton conduction along the membrane, with important bioenergetic consequences. Taken together, our combined findings suggest that SC forms as a result of a complex interplay between molecular interactions and membrane strain effects that control the functional dynamics of the OXPHOS proteins.

Materials and Methods

Atomistic MD Simulations

Atomistic molecular dynamics (MD) simulations of the ovine Complex I (CI), Complex III₂ (CIII₂), and the I-III₂ supercomplex (SCI/III₂) were performed in a POPC:POPE:cardiolipin (2:2:1) membrane containing 5 mol% QH₂ / Q (1:1 ratio). The system was solvated with TIP3P water molecules and 150 mM NaCl. The simulations were performed in different ligand states (apo, Q/QH₂ bound, Table S1) at T=310 K and p=1 bar using NAMD2.14 (68) with a 2 fs integration timestep and long-range electrostatic interactions treated using the Particle Mesh Ewald approach. The systems comprised 0.8-1.65 million atoms and were modeled using the CHARMM36 force field (69) in combination with in-house DFT-based parameters (70) of the co-factors. Protonation states were established based on electrostatic calculations with Monte Carlo sampling techniques using APBS/Karlsberg+ (71-73).

Molecular dynamics (MD) simulations of the ovine CI were conducted based on a cryo-EM structure (PDB ID: 6ZKC (44)), with simulations performed for both the quinol-bound and apo forms. Unresolved regions in the supernumerary subunits NDUFA7 and NDUFB6 were modeled using ColabFold (74). The system was equilibrated for 100 ns in the QH₂-bound state, which was then used to propagate both states for 2×0.5 μs each.

The ovine Complex III₂ was modeled based on the cryo-EM structure (PDB ID: 6Q9E (40)), with simulations also performed in different ligand states (Table S1). After equilibration of the system for 100 ns with the distal and proximal Q sites occupied, the system was propagated in each state for 2×0.5 μs.

A fully atomistic model of the SC was constructed based on the high-resolution structures of CI and CIII₂, that were merged into the experimental structure of the supercomplex (PDB ID: 6QBX (40)). After equilibration of the system for 100 ns with all Q sites occupied in CI and CIII₂, the complete system was propagated for 2×0.5 μs in each state (see Table S1).

Membrane systems with 1:1 POPC/POPE or CLD were constructed using CharmmGUI with a membrane area of 80×80 Ų, a hydration layer of 45 Å, and NaCl or KCl concentrations of 150 mM. See Table S2 for a detailed description of the membrane compositions. See the SI Appendix, for a detailed description of the simulations and analysis.

Coarse-grained MD simulations

Coarse-grained MD (cgMD) simulation models of the ovine CI, CIII₂, and the SCI₁/III₂ were created based on the atomistic models using the MARTINI3 force field (75) and Gromacs (76). All simulations were constructed with identical amounts of lipid molecules, a 2:2:1 POPC:POPE:cardiolipin ratio, and a 3 mol% mixture of quinone/quinol embedded in a simulation system with dimensions 47×47×31 nm³ (for the SC, see Table S1). All simulations were carried out at 310 K in an NPT ensemble using the velocity rescaling thermostat (77, 78), the Parrinello-Rahman barostat (79) and a 20 fs timestep using GROMACS (76). The protein structures were stabilized with elastic networks on backbone beads using a cutoff distance of 0.9 nm with a force constant of 500 kJ mol⁻¹ nm⁻². The elastic network was also applied between residues of different subunits, but not between CI and CIII₂, whilst cgMD parameters for all cofactors were also developed. Two replicas were run for CI and CIII₂ (2×50 μs each), and for the SCI/III₂ (2×50 μs), as well as 23 μs cgMD simulations of the membrane with the same number of lipid and Q/QH₂ molecules as in the protein simulations. Additional 5 μs cgMD simulations of the membrane systems, as well as for all protein models, were performed with longer lipids (0.44 nm, 0.50 nm, and 0.53 nm instead of 0.47 nm).

Lattice model of SC formation

A lattice model of the CI and CIII₂ was constructed (Fig. 4A,B) by modeling the OXPHOS proteins in unique grid positions on a 2D N×N lattice, with N=[6,15]. The protein-protein interaction was described by specific interactions (E_specific) and non-specific interactions (E_non-specific), whereas the membrane-protein interaction determined the strain energy of the membrane (E_strain), based on the number of neighboring “lipid” occupied grids that are in contact with proteins (Fig. 4A). The interaction between the lipids was indirectly accounted for by the background energy of the model. The proteins could occupy four unique orientations on a grid (g_i=[North, East, South, West]). The total Hamiltonian of the lattice model can be written as,

with

and leading to the partition function,

where g_i is the degeneracy of the state, and β=1/k_BT (k_B, Boltzmann’s constant; T, temperature). The probability of a given state i was calculated as,

with the free energy (G) defined as,

The conformational landscape was sampled by Monte Carlo (MC) using 10⁷ MC iterations with 100 replicas (see SI Appendix). Temperature effects were modeled by varying β, and the effect of different protein-to-lipid ratios by increasing the grid area.

Statistical model of the membrane distribution in the IMM

Based on the experimental protein copy numbers (cf. Ref (7)) and average protein areas, the proteins in the IMM were modeled as randomly positioned circles on a membrane square with a side length of 163 nm (see Table S5). The system was minimized according to the Hamiltonian,

where d_ij is the distance between two proteins i and j with respective radii, r_i and r_j, while the force constant k, were introduced to avoid steric clashes and set to 1 kcal mol^-1 nm^-1.

Acknowledgements

This work was supported by grants from the European Research Council (ERC) under the European Union’s Horizon 2020 research and innovation program/grant agreement 715311, the Swedish Research Council (VR), and the Knut and Alice Wallenberg foundation. (2019.0251, 2019.0043). We are thankful for computing time provided by the Partnership for Advanced Computing in Europe (PRACE project: pr127) to access Piz Daint hosted by the Swiss National Supercomputing Center (CSCS). This work was also supported by the National Academic Infrastructure for Supercomputing in Sweden (NAISS 2024/1-28, NAISS 2023/1-31, 2023/6-128) and the Swedish National Infrastructure for Computing (SNIC 2022/1-29, 2022/6-190).