MainProteins commonly function in the form of complexes, such as transport1, signal transduction2 and gene regulation3. Structural characterization of these assemblies is crucial for understanding molecular mechanisms and related drug discovery. However, experimental determination of protein complex structures is expensive or even prohibitive owing to some technical problems4,5,6. Therefore, computational methods are necessary for the prediction of protein complex structures7.Traditional computational methods for protein complex structure prediction include template-based and docking-based approaches. Template-based methods8,9,10 leverage homologous structures to address interface flexibility but are limited by template availability and quality. Docking-based approaches11,12,13, initially developed for dimers, incur prohibitive computational costs for multimers owing to the combinatorial explosion of multiple-chain arrangement. Although symmetry constraints14,15,16 and advanced sampling strategies17,18 can reduce computational costs, the performance of docking methods remain limited in handling multi-chain systems.Recent advances in deep learning have revolutionized protein-structure prediction19,20,21,22,23 and related applications24,25,26,27,28,29,30. AlphaFold2 (AF2)19-derived variants enable complex modeling without retraining31,32,33,34, whereas the retrained AlphaFold-Multimer (AFM)35 exhibits superior performance in predicting complex structures. However, these tools face major challenges in large assemblies, including degraded accuracy with increasing number of chains36, insufficient conformational diversity and memory limitation of graphics processing units (GPUs). AFsample237 enhances conformational diversity but does not overcome memory limitations. For example, AFM v2 can handle only ~3,300 residues on a 40-GB NVIDIA A100 GPU. AlphaFold3 (AF3)38 extends the input limit to 5,000 residues on its server and ~5,120 residues for local execution on an 80-GB GPU. However, it is still a challenge for structure prediction of larger proteins.To address the challenge in modeling large protein complexes, assembly-based methods have been proposed, which first predict smaller subcomponents using AlphaFold and subsequently assemble them into larger complexes. MoLPC36, for example, uses Monte Carlo tree search to assemble complexes from dimers or trimers. CombFold39 applies a hierarchical combinatorial assembly strategy, enabling broader exploration of inter-subunit interactions and yielding higher-quality predictions. Nevertheless, these methods remain inherently dependent on the quality of subcomponents predicted by AFM. Symmetry-aware approaches, such as EvoDock40 and SymProFold41, improve modeling of symmetric complexes but lack generalizability to heteromeric and asymmetric assemblies. Several methods have leveraged deep learning to guide assembly42,43, yet their improvement remains limited in realistic applications.To overcome the limitations of current methods, we present HDOCK-Multimer (HDM), an automated framework that integrates docking and assembly strategies for accurate and efficient modeling of large protein complexes. HDM incorporates a genetic algorithm (GA) into combinatorial assembly for better leveraging the correct pairwise interactions of predicted subcomplexes. It further introduces efficient multi-body and symmetric docking approaches to enhance sampling and reduce dependence on the accuracy of predicted subcomponents. Across three benchmark datasets, HDM achieves top-1 success rates of up to 86% with template modeling score (TM-score) > 0.7, outperforming state-of-the-art assembly approaches and deep learning-based methods.ResultsThe workflow of HDMHDM takes subunit sequences as input and outputs predicted complex structures. HDM comprises three stages: (1) generation of monomeric and subcomponent structures, (2) determination of modeling strategies, and (3) structure modeling followed by scoring and clustering (Fig. 1a). In the first stage, monomer and subcomponent structures are generated using AF2 and AFM. Here, subcomponents are defined as all possible trimeric subunit combinations, which can capture more cooperative interactions than pairwise combinations while avoiding expensive computational cost or GPU memory overflow of tetrameric combinations. AF2 outputs five models for each prediction, and only the top-ranked model is retained. Monomers are ranked by the average predicted local distance difference test (pLDDT) score, whereas subcomponents are selected by model confidence scores (equation (1)). In the second stage, HDM selects appropriate modeling strategies, including asymmetric docking, symmetric docking and combinatorial assembly, based on the complex homology and pairwise docking results. Symmetric docking is restricted to homo-oligomers and supports multiple symmetry types, including cyclic (Cn), dihedral (Dn), tetrahedral (T) and octahedral (O) symmetries. Asymmetric docking is applicable to homo-oligomers as well as heteromeric complexes containing fewer than five unique chains. For complexes with five or more unique chains, only the assembly strategy is adopted as the docking accuracy decreases substantially with increasing number of chains. In the third stage, HDM constructs complex structures using selected strategies. Models generated by different approaches are evaluated using our iterative knowledge-based scoring function (IT-score) for protein–protein interactions44 and ranked accordingly. To reduce redundancy, all models are clustered using a Cα root mean square deviation (RMSD) cut-off of 5 Å, with only the best-scored model retained within each cluster. By default, HDM outputs the top-ten models.Fig. 1: The framework of HDM.Full size imagea, An overview of HDM. The input is the sequences of subunits in the complex and the output is the predicted complex structure. HDM integrates multiple modeling strategies, including asymmetric docking, symmetric docking and assembly methods, to comprehensively explore the conformational space of protein complexes. b, Asymmetric docking workflow. The full complex is iteratively constructed through pairwise docking and assembly of subunits subject to score and additional interaction screening. c, Tetrahedral docking workflow. C2 and C3 symmetric candidates are screened for compatibility. Valid combinations undergo sliding optimization; otherwise, candidates are generated from standard tetrahedral axes and rotational sampling before optimization. d, Octahedral docking workflow. C2–C3 and C4 symmetric candidates are screened for compatibility. Valid combinations undergo sliding optimization; otherwise, candidates are generated from standard octahedral axes and rotational sampling before optimization. e, Combinatorial assembly workflow. Representative subunits and pairwise transformations extracted from predicted subcomponents are used for random assembly sampling, followed by GA optimization and IT-score re-ranking. Here nA and nB denote the representative structures of subunits A and B, TA→A′ and TB→B′ stand for the transformations from the representative structures of subunits A and B to their corresponding positions A′ and B′ in the subcomponent, and eAB represents a pair of transformations that defines a possible interaction between subunits A and B.Performance on Benchmark 1 (35 large heteromers)HDM is first validated on Benchmark 1 of 35 large heteromeric complexes with 5 to 20 chains per complex. When evaluated by TM-score45,46, HDM substantially outperforms competing methods, including CombFold, MoLPC and AFM version 2.2 (AFM v2). Specifically, HDM obtains top-1 success rates of 83% and 66% for acceptable accuracy and high-quality predictions, respectively, compared with 74% and 54% for CombFold, 29% and 11% for MoLPC, and 29% and 26% for AFM v2 (Fig. 2a and Supplementary Fig. 1a). Similar improvements are also observed when considering the top-ten predictions (Fig. 2a). In addition, HDM tends to perform better with the higher-quality-weighted pairwise connectivity (QPC), giving a correlation coefficient of 0.7 (Supplementary Fig. 2a). The performance of HDM also shows a negative correlation with the increasing complex size (Supplementary Fig. 2b).Fig. 2: Performance of HDM and other methods on Benchmark 1.Full size imagea, Top-N (N = 1, 10) success rates of HDM, CombFold, MoLPC and AFM v2, evaluated by TM-score. MoLPC produces only one prediction per target. b, TM-score comparison of the models predicted by CombFold and HDM. c, TM-score comparison of the models predicted by MoLPC and HDM. d, TM-score comparison of the models predicted by AFM v2 and HDM. e, Top-N (N = 1, 10) success rates of HDM, CombFold, MoLPC and AFM v2, evaluated by DockQ. f, DockQ comparison of the models predicted by CombFold and HDM. g, DockQ comparison of the models predicted by MoLPC and HDM. h, DockQ comparison of the models predicted by AFM v2 and HDM. i, Polytomella complex I (PDB 7ARC): high-quality HDM model (left), inaccurate model for CombFold (middle) and MoLPC (right). The predicted models (colored by chain) are superimposed on the native structure (gray). The dashed circles indicate that those chains assembled into incorrect positions. The dashed boxes show enlarged views of wrongly assembled chains. j, EIF2B–eIF2 complex (PDB 6I3M): acceptable-quality HDM model (left), high-quality CombFold model (middle) and inaccurate MoLPC model (right). The predicted models (colored by chain) are superimposed onto the native structure (gray). Dashed circles indicate those chains that assembled into incorrect positions.Source dataThe head-to-head TM-score comparisons indicate that HDM improves the performance over other methods in different ways (Fig. 2b–d). Compared with CombFold, HDM mainly improves predictions for targets where CombFold fails with TM-score 0.8 for many cases with interface similarities below 0.1, indicating that HDM is not primarily driven by previously observed interfaces, supporting its ability to assemble novel complexes.Computational timeAF2 and AFM predictions are performed on NVIDIA A100 GPUs with 40 GB of memory. The runtime of AFM varies depending on the size of structures. For example, predicting a subcomponent comprising 2,600 amino acids typically takes about 35,000 s. Docking and assembly procedures are executed on 2.6-GHz Intel central processing units (CPUs) using multi-process parallelization, with 20 CPU cores used by default. The average running times for the assembly method on two benchmarks are 5,635 s and 5,587 s, respectively. The docking stage generally consumes more time than the assembly process, particularly for heteromeric complexes. On Benchmark 2, the average running time for asymmetric docking is 10,252 s, whereas symmetric docking takes an average of 4,543 s. Overall, the running time for building a full complex structure can be normally finished within hours.DiscussionWe have proposed a hybrid modeling framework for structure prediction of large protein complexes by integrating ab initio docking and combinatorial assembly from predicted subcomponent structures. Extensive benchmark validation demonstrates the superior performance and scalability of HDM over existing state-of-the-art methods. Nevertheless, there is still room for further improvement. First, our symmetric docking method covers only a subset of symmetries. Extending the framework to additional symmetries, such as helical symmetry, will broaden its applicability. Second, current model ranking relies on simply summing up pairwise interaction scores, which may not fully capture higher-order cooperativity and collective effects in large assemblies. Therefore, there is a need to develop a true multi-body scoring function. Third, incorporating AlphaFold-predicted interface information during early docking stages may improve multi-body docking, particularly for heteromeric complexes. Fourth, HDM currently requires known stoichiometric information, although we have demonstrated its ability in predicting or narrowing down the stoichiometry space. Enhancing this capability will much reduce computational costs and enable broader applications. Fifth, subcomponent selection currently relies on AFM-derived confidence scores, whereas alternative metrics, such as local interaction score51 and predicted DockQ score33, may provide alternative effective assessments. Sixth, for convenience, we use RMSD to cluster sampled models in the current HDM protocol. Other metrics such as fraction of common contacts or orientational RMSD could also be explored to remove the redundant complexes.Although HDM shows the capacity to predict stoichiometry, it should be noted that none of these tools (HDM, CombFold and so on) provide truly absolute stoichiometry predictions, as this ultimately relies on the quality of predicted subcomponents. Relative stoichiometries may be inferred from complementary experimental techniques such as mass spectrometry. Integrating experimental data, such as crosslinking mass spectrometry52, fluorescence resonance energy transfer53 and cryo-electron microscopy (cryo-EM)54, into the prediction process may further improve the modeling accuracy. Currently, HDM supports crosslinking restraints, as demonstrated on target PDB 7EGF (Supplementary Fig. 7). The future study will focus on incorporating additional experimental data. In addition, HDM also shows promise for modeling antibody–antigen complexes. Among four antibody–antigen targets in Benchmarks 1 and 2, HDM achieves successful predictions on two targets (PDB 6SSI and 4RDQ), but fails on the other two (PDB 2VYR and 6X04) (Supplementary Table 4). A more extensive evaluation on antibody–antigen complexes is needed to demonstrate the robustness of HDM on such targets.MethodsBenchmark datasetsThree benchmark datasets were used to evaluate HDM. Benchmark 1, taken from the study of CombFold39, consists of 35 heteromeric complexes with 5 to 20 chains and 1,300 to 8,000 amino acids per complex (Supplementary Fig. 8a). All complexes in Benchmark 1 were released after April 2018 and not included in the training set of AFM v235. HDM was compared with CombFold, MoLPC36 and AFM v2 on this benchmark. Benchmark 2 is derived from the study of MoLPC after excluding overlapped complexes with Benchmark 1. The resulting set contains 172 complexes with size ranging from 500 to 10,000 amino acids and 10 to 30 chains per complex (Supplementary Fig. 8b), including 114 homomers and 58 heteromers. Benchmark 2 comprises 48 cyclic (Cn), 70 dihedral (Dn), 21 tetrahedral (T), 6 octahedral (O) and 3 helical (H) symmetric complexes, as well as 24 asymmetric complexes, which were all released before January 2022. On benchmark 2, HDM was compared with CombFold, MoLPC and the latest end-to-end deep learning method, AF338. Benchmark 3 consists of 7 CASP1555 targets with more than 3,000 amino acids, including 2 homomers, 3 heteromers and 2 single-chain proteins. In addition to the methods mentioned previously, HDM was further compared with the top-ranking competitors from CASP15 on this benchmark.In this study, the smallest component used for docking or assembly is called a subunit, which can be a full chain or part of a long chain. Across three benchmark datasets, full-length chains are used as subunits for nearly all targets except a few cases as follows. For the nucleosome-free yeast Swi/Snf complex (PDB 7C4J) in Benchmark 2, the two long chains are divided into two subunits, with the first 1,000 residues forming one subunit and the remaining residues forming the other. For two long single-chain targets in Benchmark 3 (T1165 and T1169), each chain was partitioned into domains according to IUPred348. Predicted disordered linkers between domains were removed, and each domain was split into equally sized segments as input subunits. Notably, during the docking stage, the full-length chains were still used as subunits. Therefore, our HDM method can predict the structure of large protein complexes of any sizes.AlphaFold structure predictionHere, AlphaFold is used to predict the structures of monomers and trimeric subcomponents. Here, trimers are selected as a balance between accuracy and computational cost, as dimers will result in reduced modeling accuracy (Supplementary Fig. 9), whereas exhaustive tetramer prediction is often prohibitively expensive. For monomers, we first search each sequence against the AlphaFold Protein Structure Database (AFDB)56. If a sequence-matching hit is retrieved, we retain the structure with the highest average pLDDT value. Otherwise, the monomer structure is predicted using AF2, and the model with the highest pLDDT is selected for docking. For trimers, all possible trimeric subunit combinations are generated using AF2. It is noted that replacing the software used to predict subcomponent structures, for example, from AF2 to AF3, does not substantially affect HDM performance (Supplementary Fig. 10). Only the highest-confidence subcomponent structures are retained for the subsequent assembly process. The model confidence of AF2-predicted multimer is defined as:$$\mathrm{Model}\,\mathrm{confidence}=0.8\times \mathrm{ipTM}+0.2\times \mathrm{pTM}$$(1)where ipTM denotes the interface predicted TM-score, which assesses the confidence of the predicted interfaces between chains, and pTM denotes the predicted TM-score, assessing the confidence of the overall structure. For comparison purpose, AF2 and AF3 are also used to predict the structures of the complexes in Benchmarks 1 and 2, respectively.As the disordered terminal regions of AF2-predicted monomers may interfere with the docking process by introducing steric clashes and unfavorable scores, they are therefore trimmed before docking. For chains with an average pLDDT ≥ 80, we remove terminal residues with pLDDT < 70. For chains with an average pLDDT between 50 and 80, secondary structure is predicted by STRIDE57, and terminal residues with pLDDT < 50 or outside helices and sheets are removed. Chains with an average pLDDT < 50 are retained without any trimming. After model generation, removed terminal regions can be reconstructed by homology modeling, MD simulations or structure alignment (Supplementary Fig. 11).Modeling strategiesHDM uses three modeling strategies, including asymmetric docking, symmetric docking and combinatorial assembly, to generate full complex structures. For each target, one or more strategies are selected based on its sequence composition, chain number and structural symmetry. For homomers and heteromers with fewer than five unique chains, asymmetric docking and assembly strategies are selected by default. For complexes comprising five or more unique chains, only the assembly strategy is used. Symmetric docking is also applied to homomeric complexes with certain symmetry (Fig. 1a).As the exact symmetry type may be unknown in real applications, HDM does not require explicit symmetry input and instead infers candidate symmetries, including cyclic, dihedral and cubic symmetries, from chain number and pairwise docking results. Specifically, the total chain number N is factorized as N = n1 × n2, where each factor pair (n1, n2) corresponds to two global cyclic symmetries, \({C}_{{n}_{1}}\) and \({C}_{{n}_{2}}\). The pair (1, N) corresponds to CN symmetry. If n1 = 2, both \({C}_{{n}_{2}}\) and \({D}_{{n}_{2}}\) symmetries are considered. For N = 12 or 24, tetrahedral (T) or octahedral (O) symmetries are also included.Subunits are then docked pairwise using HDOCK, with rotational symmetry inferred from docking results to further constrain the symmetric search space. Each of the top-10 transformations from docking solutions is iteratively applied m times (m = 1, 2, …, N − 1) to generate assemblies containing n chains (n = 2, 3, …, N). For each n that divides N, steric clashes are quantified as the number of backbone atom pairs within 2 Å across all chain pairs. Structures with clashes below the following threshold are retained:$$\mathrm{clash}\_\mathrm{count}\_\mathrm{threshold}=m\times \mathrm{clash}\_\mathrm{count}\_2\mathrm{chains}+50$$(2)where clash_count_2chains represents the number of clashes in the dimer (m = 1).For retained candidates, the potential Cn symmetry axis is estimated from chain centroids. Each chain is rotated by 2π/n around the axis, and the Cα RMSD between rotated and adjacent chains is calculated. Structures with RMSD ≤ 2 Å are considered Cn symmetric. If multiple Cn symmetries are detected, the highest-order one is retained, denoted as Cn(pair). The detected symmetry is subsequently validated by Cn symmetric docking using CHDOCK. Among the top-10 CHDOCK models, those that can be superimposed onto Cn(pair) within 2-Å Cα RMSD are denoted as Cn(chdock), with the lowest-energy Cn(chdock) model retained for asymmetric docking. The factor pair (n, N/n) determines the symmetric configurations. If no valid Cn symmetry is identified, the monomer is used for asymmetric docking and all compatible symmetry types are explored. Taking N = 12 as an example, detection of C2 symmetry leads to the sampling of C2, C6, D6 and T symmetries, whereas C3 detection limits sampling to C3, C4 and T symmetries. No symmetry detection requires exploring all seven symmetries (C2, C3, C4, C6, C12, D6, T).Asymmetric dockingThe asymmetric docking method is based on our rigid-body protein–protein pairwise docking method, HDOCK13, which uses a fast Fourier transform-based search strategy. During HDOCK sampling, the receptor is fixed and the ligand is rotated by a 15° Euler angle interval, producing 4,392 orientations. At each orientation, both proteins are projected onto a three-dimensional grid with 1.2 Å spacing for translational search. The calculation of shape complementarity is accelerated using three-dimensional fast Fourier transform58, and the top-ten translations for each rotation are refined using our iterative scoring function, IT-score44. The best-scored docking pose is retained for each rotation, yielding a total of 4,392 binding modes.For multi-body docking without symmetry constraints (Fig. 1b), each chain is treated as a node ni (i ∈ [1, N]), and pairwise docking interactions are considered as edges eij (i, j ∈ [1, N]). Thus, complex modeling becomes a minimum spanning tree construction problem. The method integrates principles from Prim’s59 and Kruskal’s60 algorithms and iteratively updates edge scores based on docking predictions to guide assembly.HDOCK first performs pairwise docking for all input chains, with only the best-scored transformation retained as the initial edge for each pair. These edges are sorted in ascending order of IT-score, and the minimum spanning tree construction starts from the best-scored edge. When adding a new node nj to a partial tree Tk through edge eij, both score and interface screening are applied. The IT-score of eij cannot exceed half of the lowest edge score in Tk (denoted as min_scorek), that is, IT_score(eij) ≤ min_scorek/2; otherwise, the expansion of Tk is terminated. The expansion is also discarded when more than 30 backbone atoms of nj are within 10 Å of any chain in Tk except ni, thereby avoiding undesired interfaces. If the full complex is not obtained in the first iteration, the resulting subcomponents are treated as new nodes for subsequent pairwise docking and interaction graph updates. The docking and screening cycle is iterated until the full assembly is achieved. For homomers or complexes containing homologous chains, the asymmetric docking process can be substantially simplified because identical chains require only one docking run, enabling simultaneous incorporation of multiple chains and reducing computational cost.Symmetric dockingFor symmetric docking, we used our previously developed rigid-body docking methods, CHDOCK61 and DHDOCK16, for handling Cn and Dn symmetries, respectively. CHDOCK sets the z axis as the rotational symmetry axis and samples reduced rotational space (ϕ = 0, θ ∈ [0, π/2], ψ ∈ (0, 2π]) with an interval of 10°, generating 360 orientations. The translational sampling and scoring procedures are similar to those in HDOCK, ultimately producing 360 Cn symmetric binding modes. On the basis of Cn symmetric models generated by CHDOCK, DHDOCK performs additional C2 symmetric sampling along axes perpendicular to the corresponding Cn axes.C n/D n symmetric dockingTo construct a \({C}_{{n}_{1}}\) symmetric structure with N chains, we adopt a stepwise cyclic symmetry sampling strategy, yielding a \({C}_{{n}_{2}}\)–\({C}_{{n}_{1}}\) symmetric mode, where N = n1 × n2. Here, \({C}_{{n}_{2}}\) and \({C}_{{n}_{1}}\) represent local and global symmetries, respectively. CHDOCK first performs \({C}_{{n}_{2}}\) symmetric docking on the monomer, and the top-10 models are retained for subsequent \({C}_{{n}_{1}}\) symmetric docking, yielding 100 candidate models. Candidates are re-ranked by IT-score and clustered, with the 10 best-scored structures selected as the final output. For CN symmetry, only one cyclic sampling round is required. \({D}_{{n}_{1}}\) symmetric complexes are modeled using DHDOCK, and the ten best-scored models are retained as output.T symmetric dockingTetrahedral (T) symmetric structures, characterized by both C2 and C3 symmetries, consist of 12 chains for the simplest homomeric case. Therefore, for homomers with 12 chains, additional T symmetric docking is performed. C2 and C3 symmetric docking are first applied to generate candidate local symmetries. Compatible C2–C3 combinations are identified by: (1) superimposing a dimer chain onto a trimer, (2) expanding a new trimer from the other chain in the dimer, and (3) verifying whether the expanded chains form C2 symmetry. Valid combinations undergo sliding optimization, where each trimer is translated along its cyclic axis to minimize steric clashes. For each combination, the best-scored structure is retained. The top-ten structures are selected as the final output. If no valid combination is found, each candidate trimer is replicated three times, with its C3 axis aligned to a standard tetrahedral axis. These trimers are then rotated around their axes in steps of 15° until C2 symmetries are formed between different chain pairs, followed by sliding optimization to obtain the final output (Fig. 1c).O symmetric dockingFor homomers comprising 24 chains, we perform O symmetric docking. C2–C3 and C4 symmetric docking are first conducted to generate candidate local symmetries, with compatible C2–C3–C4 combinations identified by: (1) superimposing two chains of a tetramer onto a hexamer, (2) expanding a new hexamer with the remaining two chains as anchors, and (3) verifying C4 symmetries of the non-anchor chains. Valid hexamer sets also undergo sliding optimization to obtain the top-ten structures. If there is no valid C2–C3–C4 combination, O symmetric systems are constructed using tetramer structures and standard octahedral axes (Fig. 1d).GA-based combinatorial assemblyThe assembly workflow of HDM is illustrated in Fig. 1e, which consists of three stages: (1) extracting representative subunits and pairwise transformations, (2) initial random sampling, and (3) refinement via GA iteration. The third stage is only applied when the initial sampled structures exhibit suboptimal confidence. Final models can be relaxed optionally to reduce steric clashes.Extracting representative subunits and pairwise transformationsSubunits are first extracted from predicted subcomponents and ranked by their average pLDDT scores, with the highest-ranking structures selected as representatives. Pairwise interactions are identified between subunit pairs if any Cα–Cα distance is within 8 Å. For interacting pairs, transformations from the representative subunits to the corresponding structures are calculated. Considering possible structural variations in identical subunits from different models, alignment is performed by MMalign62,63 to obtain more precise transformations. To improve alignment, residues with pLDDT < 80 are removed from representative subunits with median pLDDT > 90. Each transformation is scored by the average pLDDT values of interface residues. Local symmetry often exists in the complexes with multiple identical chains. Therefore, for complexes with five or more identical subunits, transformation scores extracted from homomeric subcomponents are rescaled using S = S + S × (1 − S/100), where S is the original score.Initial sampling stageThis stage includes 500 independent random sampling rounds. In each round, all transformations are randomly shuffled and traversed sequentially. Three assembly operations are applied to construct new structures:Create: assemble two unpaired subunitsAppend: incorporate an unpaired subunit into an existing subcomplexMerge: combine two different subcomplexes into a larger complexTo avoid steric clashes from inaccurate transformations, clash examination is performed for the backbone atoms of residues with pLDDT > 80 in each operation. Clashes are defined as backbone atoms whose distances are within 2 Å across different subunits. If the clash ratio exceeds a predefined threshold, the operation is discarded. Thresholds vary by operations to account for the effects of subcomplex size: 2% of the shorter subunit/subcomplex length for Create and Merge, and 5% of the unpaired subunit length for Append. For subunits from the same chain, additional sequence connectivity is checked. The operation is discarded if the distance between consecutive residues from two subunits exceeds the number of linker residues multiplied by 4 Å.For fully assembled complexes, the confidence score is calculated as a weighted transformation score39, and defined by the following expression:$$\mathrm{Confidence}\,\mathrm{score}=\sum {w}_{i}{S}_{i}\,/\,\sum {w}_{i}$$(3)where Si represents the transformation score and the weight wi is the number of residues in the shorter subunit or subcomplex in each assembly operation.Each failed sampling yields a set of partial complexes, which is assigned a confidence score based on its largest subcomplex. After 500 iterations, fully assembled complexes are ranked by their confidence scores and clustered. The highest-confidence structure in each cluster is retained as the representative. The top-100 representatives are selected to form the initial population. The population quality is represented by the lowest confidence score of the top-20 structures, denoted as SPQ, and classified into three categories: high quality (SPQ ≥ 85.0), medium quality (80.0 ≤ SPQ < 85.0) and low quality (SPQ < 80.0). Medium- and low-quality populations proceed to the GA optimization stage. High-quality populations or those with fewer than ten clusters are re-ranked by IT-score and clustered, with the top representatives selected for output. If all 500 iterations fail to build the full complex, the set of structures containing the largest and highest-confidence subcomplex is retained. Each chain in the set is filtered using the same criteria as in monomer preprocessing, and the full complex is subsequently obtained using the multi-body docking algorithm.GA optimization stageWe adopt a GA to optimize medium- and low-quality populations. In this work, the population size M is set to 100, and the mutation rate is set to 0.1. Two GA operations, crossover and mutation, are used to generate new structures. Crossover randomly selects two parent structures, removes up to half of the subunits from each, and reconstructs the complex by reapplying transformations from parents. Each crossover operation can yield up to two new candidates. Mutation randomly removes a subunit from a newly generated structure and rebuilds the complex using the transformations from the full set. Both operations follow the same assembly process as that for the initial sampling stage. The clash threshold for the Append operation is adjusted to 2% in mutation to reduce the risk of erroneous operation.Each GA generation performs 2M (that is, 200 here) operations. All structures are ranked by confidence score and then clustered. The top-M representatives are retained for the next generation. We implement a dynamic termination strategy that stops the iteration once the population quality improves by one level (that is, medium-quality populations continue until SPQ ≥ 85.0; low-quality until SPQ ≥ 80.0). In addition, an early-stop criterion is applied to avoid redundant computation under low-quality transformations. Specifically, the iteration is terminated if no new structures are generated or if SPQ fails to improve over three consecutive generations. Upon termination, the final population is re-ranked by IT-score and clustered. The top-ten models are selected for output.Structure relaxationStructural variations between representative subunits and their corresponding structures in the subcomponents may introduce steric clashes during assembly. Therefore, an additional relaxation step is recommended to alleviate steric clashes and improve the reliability of IT-score evaluation. Structural relaxation is also advised during the sliding optimization step in the construction of cubic symmetric systems. In this study, all models generated by assembly process and O symmetric docking approach are refined via 500 steps of MD simulations using the ff14SB force field64 of AMBER (version 14).Model scoringFinal complexes are ranked by IT-score44, an iterative knowledge-based scoring function for protein-–protein interactions. Initially developed for evaluating pairwise interactions, IT-score is derived by iteratively optimizing interatomic pair potentials to maximize the discrimination between native and predicted structures based on the distributions of atom pairs. For a dimer, the score is computed as a sum of the interatomic pair potentials within a 10 Å cut-off between chains. Here, we extend the IT-score to multimeric systems by decomposing the complex into chain pairs and summing the pairwise scores. Models are ranked according to IT-score, with lower scores indicating higher structural quality.Model clusteringTo balance the result diversity and computational efficiency, different RMSD thresholds are applied at different clustering stages. A backbone RMSD of 5 Å is adopted for local symmetry sampling in Cn and Dn symmetric docking stage, whereas a 5 Å Cα RMSD threshold is applied for global sampling, assembly, and final model clustering.For complexes containing homologous subunits, the optimal subunit correspondences are identified using a heuristic search strategy39 to avoid exhaustive enumeration. Specifically, the centroids of each chain are calculated for both complexes. An initial superposition is performed by aligning the centroids in their original order. Then, for each pair of homologous subunits in the second complex, swap centroids and recalculate the centroid RMSD. If a swap leads to a lower RMSD, the new correspondence is retained. This process is iterated until no further RMSD reduction is observed. The final chain mapping is used for calculating the Cα RMSD between two complexes, thereby enhancing structural diversity in the final output.Integration of crosslinking experimental dataWe implement an extension to the assembly algorithm of HDM for allowing the integration of distance restraints derived from crosslinking experimental data. Specifically, after each assembly operation, we calculate the satisfaction ratio of crosslink-derived distance restraints involving the newly formed subcomplex. A restraint is considered satisfied if the Cα–Cα distance between the corresponding residues is below the specified threshold. For restraints involving repeated subunits, the restraint is considered satisfied if at least one corresponding pair meets the distance threshold. The satisfaction ratio of a subcomplex is defined as follows.$$\mathrm{Satisfaction}\,\mathrm{ratio}=\sum _{\mathrm{satisfied}}{w}_{1}{w}_{2}\,/\,\sum _{\mathrm{all}}{w}_{1}{w}_{2}$$(4)where w1 represents the experimental confidence of the restraint, such as the false discovery rate52, capturing experimental uncertainty. w2 represents the average pLDDT score of the residues involved in the satisfied restraint and is used to account for uncertainty in predicted subunit structures. If the satisfaction ratio of a newly formed subcomplex falls below the predefined threshold, the corresponding assembly operation is discarded. During the initial sampling and GA optimization stage, the satisfaction ratio threshold is set to 0.7 and 0.8, respectively. For the final full-complex models, both the confidence score and IT-score are multiplied by the satisfaction ratio to prioritize models fulfilling more distance restraints.Reintegrating disordered terminal regionsTo avoid interference in docking, the automated HDM pipeline removes disordered terminal regions during preprocessing. The trimmed regions can be reconstructed by homology modeling, MD simulation or structural alignment. For example, for calcium/calmodulin-dependent protein kinase (PDB 1HKX), low-confidence disordered terminal residues are removed before docking (Supplementary Fig. 11). HDM then generates a high-quality complex model (TM-score = 0.97), with the full complex obtained by aligning the full-length monomer onto each chain of the predicted model using MMalign62,63.Quality-weighted pairwise connectivityPairwise connectivity39 measures the maximum fraction of a complex that can be assembled from the predicted subcomponents. To account for the influence of monomer accuracy on assembly, we introduce QPC. Specifically, each chain is treated as a node, with low-quality nodes removed based on monomer structural similarity (TM-score45 < 0.65 between the extracted and native monomers). Edges are then added between remaining nodes if they are in contact within the native complex and a similar transformation (DockQ47 score ≥ 0.23) exists among the extracted transformations. After building the graph, we identify its connected components and calculate the QPC as follows.$$\mathrm{QPC}=\frac{{\sum }_{i\in \mathrm{LCC}}{L}_{i}{T}_{i}}{{L}_{\mathrm{total}}}$$(5)where Ltotal denotes the total residue number in the full complex, LCC represents the largest connected component, and Li and Ti are the length and TM-score of the monomer i in the LCC.Overall, QPC estimates the maximum potential to accurately assemble the full complex. Higher values indicate that a larger fraction of the complex could be accurately modeled, with a value of 1.0 indicating that the complete complex can potentially be assembled from the predicted subcomponents.Accuracy assessmentTwo metrics are used to evaluate the accuracy of predicted models, including TM-score45,46 for overall structural similarity and DockQ47 for interface quality. TM-scores are calculated using the MMalign program62,63. Consistent with previous studies36,39, models with TM-score > 0.7 and 0.8 are classified as acceptable and high quality, respectively. DockQ v2 program65 is used for calculating DockQ values. Following the CAPRI criteria66, models are considered as acceptable-, medium- and high-quality for 0.23 ≤ DockQ < 0.49, 0.49 ≤ DockQ < 0.80 and DockQ ≥ 0.80, respectively. The success rate is defined as the fraction of targets yielding at least one acceptable-quality model for top-one predictions, unless otherwise specified.Comparison with other methodsHDM was compared with CombFold39 and MoLPC36 using the same subcomponent inputs on Benchmark 1 and Benchmark 2. CombFold and MoLPC were downloaded and run locally using the publicly available implementations at https://github.com/dina-lab3D/CombFold and https://github.com/patrickbryant1/MoLPC/. The results of AFM v235 on Benchmark 1 were generated locally. The AF338 predictions on Benchmark 2 and Benchmark 3 were obtained using the AF3 server at https://alphafoldserver.com/. Restricting AF3 evaluation to only high-confidence regions slightly reduces TM-scores, probably owing to fewer residues available for alignment (Supplementary Fig. 12). Therefore, full-length AF3 predictions are used during comparative evaluations. The results of other methods on Benchmark 3 were collected from the published paper39.Model visualizationThe images for all protein structures in this study were prepared using ChimeraX67.Reporting summaryFurther information on research design is available in the Nature Portfolio Reporting Summary linked to this article.Data availabilitySource Data for Figs. 2–6 and Supplementary Figs. 1–12 are provided with this paper. All published datasets used in this paper were taken from the PDB (accession codes specified in the figure captions and in Source Data). The prediction results are available via Zenodo at https://doi.org/10.5281/zenodo.21807856 (ref. 68).Code availabilityThe HDOCK-Multimer package is freely available at https://github.com/huang-laboratory/HDOCK-Multimer or via Zenodo at https://doi.org/10.5281/zenodo.21805320 (ref. 69).ReferencesSim, S. I., Chen, Y., Lynch, D. L., Gumbart, J. C. & Park, E. Structural basis of mitochondrial protein import by the TIM23 complex. Nature 621, 620–626 (2023).Article Google Scholar Banerjee, T. et al. A dynamic partitioning mechanism polarizes membrane protein distribution. Nat. Commun. 14, 7909 (2023).Article Google Scholar Wang, H., Bienz, M., Yan, X. X. & Xu, W. Structural basis of the interaction between BCL9–Pygo and LDB–SSBP complexes in assembling the Wnt enhanceosome. Nat. Commun. 14, 3702 (2023).Article Google Scholar Robinson, C. V., Sali, A. & Baumeister, W. The molecular sociology of the cell. Nature 450, 973–082 (2007).Article Google Scholar Ward, A. B., Sali, A. & Wilson, I. A. Biochemistry. Integrative structural biology. Science 339, 913–915 (2013).Article Google Scholar Lyumkis, D. Challenges and opportunities in cryo-EM single-particle analysis. J. Biol. Chem. 294, 5181–5197 (2019).Article Google Scholar Li, H., Nithin, C., Kmiecik, S. & Huang, S. Y. Computational methods for modeling protein–protein interactions in the AI era: current status and future directions. Drug Discov. Today 30, 104382 (2025).Article Google Scholar Mukherjee, S. & Zhang, Y. Protein–protein complex structure predictions by multimeric threading and template recombination. Structure 19, 955–966 (2011).Article Google Scholar Zhang, Q. C. et al. Structure-based prediction of protein–protein interactions on a genome-wide scale. Nature 490, 556–560 (2012).Article Google Scholar Yan, Y., Wen, Z., Wang, X. & Huang, S. Y. Addressing recent docking challenges: a hybrid strategy to integrate template-based and free protein-protein docking. Proteins 85, 497–512 (2017).Article Google Scholar Pierce, B. G. et al. ZDOCK server: interactive docking prediction of protein–protein complexes and symmetric multimers. Bioinformatics 30, 1771–1773 (2014).Article Google Scholar Kozakov, D. et al. The ClusPro web server for protein–protein docking. Nat. Protoc. 12, 255–278 (2017).Article Google Scholar Yan, Y., Tao, H., He, J. & Huang, S. Y. The HDOCK server for integrated protein–protein docking. Nat. Protoc. 15, 1829–1852 (2020).Article Google Scholar Pierce, B., Tong, W. & Weng, Z. M-ZDOCK: a grid-based approach for Cn symmetric multimer docking. Bioinformatics 21, 1472–1478 (2005).Article Google Scholar Schneidman-Duhovny, D., Inbar, Y., Nussinov, R. & Wolfson, H. J. PatchDock and SymmDock: servers for rigid and symmetric docking. Nucleic Acids Res. 33, W363–W367 (2005).Article Google Scholar Yan, Y., Tao, H. & Huang, S. Y. HSYMDOCK: a docking web server for predicting the structure of protein homo-oligomers with Cn or Dn symmetry. Nucleic Acids Res. 46, W423–W431 (2018).Article Google Scholar Inbar, Y., Benyamini, H., Nussinov, R. & Wolfson, H. J. Prediction of multimolecular assemblies by multiple docking. J. Mol. Biol. 349, 435–447 (2005).Article Google Scholar Esquivel-Rodríguez, J., Yang, Y. D. & Kihara, D. Multi-LZerD: multiple protein docking for asymmetric complexes. Proteins 80, 1818–1833 (2012).Article Google Scholar Jumper, J. et al. Highly accurate protein structure prediction with AlphaFold. Nature 596, 583–589 (2021).Article Google Scholar Baek, M. et al. Accurate prediction of protein structures and interactions using a three-track neural network. Science 373, 871–876 (2021).Article Google Scholar Lin, Z. et al. Evolutionary-scale prediction of atomic-level protein structure with a language model. Science 379, 1123–1130 (2023).Article MathSciNet Google Scholar Zheng, W. et al. Deep-learning-based single-domain and multidomain protein structure prediction with D-I-TASSER. Nat. Biotechnol. 44, 641–653 (2026).Article Google Scholar Peng, Z., Wang, W., Han, R., Zhang, F. & Yang, J. Protein structure prediction in the deep learning era. Curr. Opin. Struct. Biol. 77, 102495 (2022).Article Google Scholar Li, T., Chen, J., Li, H., Cao, H. & Huang, S. Y. EMProt improves structure determination from cryo-EM maps. Nat. Struct. Mol. Biol. 33, 341–350 (2026).Article Google Scholar Chen, J., Li, T., He, J. & Huang, S. Y. Protein model building for intermediate-resolution cryo-EM maps by integrating evolutionary and experimental information. Structure 34, 375–384.e5 (2026).He, J., Lin, P., Chen, J., Cao, H. & Huang, S. Y. Model building of protein complexes from intermediate-resolution cryo-EM maps with deep learning-guided automatic assembly. Nat. Commun. 13, 4066 (2022).Article Google Scholar Lyu, J. et al. AlphaFold2 structures guide prospective ligand discovery. Science 384, eadn6354 (2024).Article Google Scholar Terwilliger, T. C. et al. Improved AlphaFold modeling with implicit experimental information. Nat. Methods 19, 1376–1382 (2022).Article Google Scholar Scardino, V., Di Filippo, J. I. & Cavasotto, C. N. How good are AlphaFold models for docking-based virtual screening? iScience 26, 105920 (2022).Article Google Scholar Wong, F. et al. Benchmarking AlphaFold-enabled molecular docking predictions for antibiotic discovery. Mol. Syst. Biol. 18, e11081 (2022).Article Google Scholar Ko, J. & Lee, J. Can AlphaFold2 predict protein–peptide complex structures accurately? Preprint at bioRxiv https://doi.org/10.1101/2021.07.27.453972 (2021).Tsaban, T. et al. Harnessing protein folding neural networks for peptide–protein docking. Nat. Commun. 13, 176 (2022).Article Google Scholar Bryant, P., Pozzati, G. & Elofsson, A. Improved prediction of protein–protein interactions using AlphaFold2. Nat. Commun. 13, 1265 (2022).Article Google Scholar Gao, M., Nakajima, A. D., Parks, J. M. & Skolnick, J. AF2Complex predicts direct physical interactions in multimeric proteins with deep learning. Nat. Commun. 13, 1744 (2022).Article Google Scholar Evans, R. et al. Protein complex prediction with AlphaFold-Multimer. Preprint at bioRxiv https://doi.org/10.1101/2021.10.04.463034 (2021).Bryant, P. et al. Predicting the structure of large protein complexes using AlphaFold and Monte Carlo tree search. Nat. Commun. 13, 6028 (2022).Article Google Scholar Kalakoti, Y. & Wallner, B. AFsample2 predicts multiple conformations and ensembles with AlphaFold2. Commun. Biol. 8, 373 (2025).Article Google Scholar Abramson, J. et al. Accurate structure prediction of biomolecular interactions with AlphaFold 3. Nature 630, 493–500 (2024).Article Google Scholar Shor, B. & Schneidman-Duhovny, D. CombFold: predicting structures of large protein assemblies using a combinatorial assembly algorithm and AlphaFold2. Nat. Methods 21, 477–487 (2024).Article Google Scholar Jeppesen, M. & André, I. Accurate prediction of protein assembly structure by combining AlphaFold and symmetrical docking. Nat. Commun. 14, 8283 (2023).Article Google Scholar Buhlheller, C. et al. SymProFold: structural prediction of symmetrical biological assemblies. Nat. Commun. 15, 8152 (2024).Article Google Scholar Gao, Z et al. Deep reinforcement learning for modelling protein complexes. Preprint at https://doi.org/10.48550/arXiv.2405.02299 (2024).Gao, Z et al. Protein multimer structure prediction via prompt learning. Preprint at https://doi.org/10.48550/arXiv.2402.18813 (2024).Huang, S. Y. & Zou, X. An iterative knowledge-based scoring function for protein–protein recognition. Proteins 72, 557–79 (2008).Article Google Scholar Zhang, Y. & Skolnick, J. Scoring function for automated assessment of protein structure template quality. Proteins 57, 702–10 (2004).Article Google Scholar Xu, J. & Zhang, Y. How significant is a protein structure similarity with TM-score = 0.5? Bioinformatics 26, 889–95 (2010).Article Google Scholar Basu, S. & Wallner, B. DockQ: a quality measure for protein–protein docking models. PLoS ONE 11, e0161879 (2016).Article Google Scholar Erdös, G., Pajkos, M. & Dosztányi, Z. IUPred3: prediction of protein disorder enhanced with unambiguous experimental annotation and visualization of evolutionary conservation. Nucleic Acids Res. 49, W297–W303 (2021).Article Google Scholar Berman, H. M. et al. The Protein Data Bank. Nucleic Acids Res. 28, 235–42 (2000).Article Google Scholar Kim, W. et al. Rapid and sensitive protein complex alignment with Foldseek-Multimer. Nat. Methods 22, 469–472 (2025).Article Google Scholar Kim, A.R. et al. Enhanced protein–protein interaction discovery via AlphaFold-Multimer. Preprint at bioRxiv https://doi.org/10.1101/2024.02.19.580970 (2024).Leitner, A., Faini, M., Stengel, F. & Aebersold, R. Crosslinking and mass spectrometry: an integrated technology to understand the structure and function of molecular machines. Trends Biochem. Sci. 41, 20–32 (2016).Article Google Scholar Dimura, M. et al. Quantitative FRET studies and integrative modeling unravel the structure and dynamics of biomolecular systems. Curr. Opin. Struct. Biol. 40, 163–185 (2016).Article Google Scholar Schmidt, C. & Urlaub, H. Combining cryo-electron microscopy (cryo-EM) and cross-linking mass spectrometry (CX-MS) for structural elucidation of large protein assemblies. Curr. Opin. Struct. Biol. 46, 157–168 (2017).Article Google Scholar Kryshtafovych, A., Schwede, T., Topf, M., Fidelis, K. & Moult, J. Critical assessment of methods of protein structure prediction (CASP)-Round XV. Proteins 91, 1539–1549 (2023).Article Google Scholar Varadi, M. et al. AlphaFold Protein Structure Database: massively expanding the structural coverage of protein-sequence space with high-accuracy models. Nucleic Acids Res. 50, D439–D444 (2022).Article Google Scholar Frishman, D. & Argos, P. Knowledge-based protein secondary structure assignment. Proteins 23, 566–79 (1995).Article Google Scholar Yan, Y. & Huang, S. Y. Pushing the accuracy limit of shape complementarity for protein–protein docking. BMC Bioinform. 20, 696 (2019).Article Google Scholar Prim, R. C. Shortest connection networks and some generalizations. Bell Syst. Tech. J. 36, 1389–401 (1957).Article Google Scholar Kruskal, J. B. On the shortest spanning subtree of a graph and the traveling salesman problem. Proc. Am. Math. Soc. 7, 48–50 (1956).Article MathSciNet Google Scholar Yan, Y. & Huang, S. Y. CHDOCK: a hierarchical docking approach for modeling Cn symmetric homo-oligomeric complexes. Biophys. Rep. 5, 65–72 (2019).Article Google Scholar Mukherjee, S. & Zhang, Y. MM-align: a quick algorithm for aligning multiple-chain protein complex structures using iterative dynamic programming. Nucleic Acids Res. 37, e83 (2009).Article Google Scholar Zhang, C., Shine, M., Pyle, A. M. & Zhang, Y. US-align: universal structure alignments of proteins, nucleic acids, and macromolecular complexes. Nat. Methods 19, 1109–1115 (2022).Article Google Scholar Maier, J. A. et al. ff14SB: improving the accuracy of protein side chain and backbone parameters from ff99SB. J. Chem. Theory Comput. 11, 3696–713 (2015).Article Google Scholar Mirabello, C. & Wallner, B. DockQ v2: improved automatic quality measure for protein multimers, nucleic acids, and small molecules. Bioinformatics 40, btae586 (2024).Article Google Scholar Collins, K. W. et al. CAPRI-Q: the CAPRI resource evaluating the quality of predicted structures of protein complexes. J. Mol. Biol. 436, 168540 (2024).Article Google Scholar Pettersen, E. F. et al. UCSF ChimeraX: structure visualization for researchers, educators, and developers. Protein Sci. 30, 70–82 (2021).Article Google Scholar Yao, X., Ya, Y., Li, H. & Huang, S. Y. HDOCK-Multimer: integrating docking and combinatorial assembly for structure prediction of large protein complexes. Zenodo https://doi.org/10.5281/zenodo.21807856 (2026).Yao, X., Ya, Y., Li, H. & Huang, S. Y. xyao7/HDOCK-Multimer: HDOCK-Multimer program (version v1.0.0) [computer software]. Zenodo https://doi.org/10.5281/zenodo.21805320 (2026).Download referencesAcknowledgementsWe thank J. He, P. Lin, H. Tao and T. Li for helpful discussions, and the computing support from the HPC Platform of Huazhong University of Science and Technology.FundingThis work was supported by the National Key Research and Development Program of China (grant numbers 2025YFC3410200 and 2025YFC3410203), the National Natural Science Foundation of China (grant numbers 32430020, 32161133002 and 62072199), the Major Project of Guangzhou National Laboratory (grant number GZNL2023A03007), and the start-up grant of Huazhong University of Science and Technology.Author informationAuthors and AffiliationsSchool of Physics and Key Laboratory of Molecular Biophysics of MOE, Huazhong University of Science and Technology, Wuhan, ChinaXuan Yao, Yifan Ya, Hao Li & Sheng-You HuangHubei Fundamental Research Center for Physics, Wuhan, ChinaSheng-You HuangAuthorsXuan YaoView author publicationsSearch author on:PubMed Google ScholarYifan YaView author publicationsSearch author on:PubMed Google ScholarHao LiView author publicationsSearch author on:PubMed Google ScholarSheng-You HuangView author publicationsSearch author on:PubMed Google ScholarContributionsS.-Y.H. conceived of and supervised the project. X.Y. and Y.Y. designed and performed the experiments. X.Y., Y.Y., H.L. and S.-Y.H. analyzed the data. X.Y. and S.-Y.H. wrote the paper. All authors reviewed and approved the final version of the paper.Corresponding authorCorrespondence to Sheng-You Huang.Ethics declarationsCompeting interestsThe authors declare no competing interests.Peer reviewPeer review informationNature Computational Science thanks Raphaelle Versini and the other, anonymous, reviewer(s) for their contribution to the peer review of this work. Peer reviewer reports are available. Primary Handling Editor: Jie Pan, in collaboration with the Nature Computational Science team.Additional informationPublisher’s note Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.Supplementary informationSource dataRights and permissionsOpen Access This article is licensed under a Creative Commons Attribution-NonCommercial-NoDerivatives 4.0 International License, which permits any non-commercial use, sharing, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if you modified the licensed material. You do not have permission under this licence to share adapted material derived from this article or parts of it. The images or other third party material in this article are included in the article’s Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by-nc-nd/4.0/.Reprints and permissionsAbout this article