BACKGROUND
proteins interacting with other proteins are known to play key roles in almost all cellular and biological processes such as signalling, metabolism, and trafficking.
therefore, understanding protein-protein interactions is fundamental to our ability to comprehend and control cell function.
protein complexes can be either permanent or transient.
polypeptide chains of permanent  complexes remain bound to each other throughout their functional lifetime.
an example hereof is the complex between β and γ subunits of hetero-trimeric g-proteins  <cit> .
in transient complexes, the subunits can either form a complex or detach from each other as a result of specific biological conditions.
the overall structures of these proteins are stable in the unbound as well as in the bound state.
the complex of cyclin and cyclin-dependent protein kinase is an example of such a transient complex  <cit> .
the molecular mechanism of interactions between protein-protein complexes is usually characterized by structure determination methods such as x-ray crystallography, nmr spectroscopy, and cryo-electron microscopy.
however, many complexes are too transient to lend themselves to experimental characterization.
docking simulations have recently gained importance as possible means for predicting the quaternary structures of protein-protein complexes  <cit> .
yet, these algorithms generate a large number of false-positives and false-negatives, mainly due to the inaccuracy of scoring functions  used to evaluate the docked conformations  <cit> .
docking algorithms that reconstruct known complexes using the structures of the proteins in the bound form  generally show reasonable levels of success.
however, when the structures in the unbound form are used as input , the performance of the same algorithms becomes rather poor.
this failure is mainly attributed to the inability of the algorithms and the scoring functions to take account of the conformational changes that  accompany complex formation  <cit> .
in an effort to overcome these deficiencies, one usually resorts to additional biological data, which assist in singling out the native  conformation among putative docked conformations  <cit> .
however, in the absence of prior knowledge regarding the binding site, exhaustive searches of up to  <dig> different positions and orientations are performed for the substrate  <cit> , followed by a refinement and clustering stage.
thus, it is difficult to take into account all the possible changes in internal conformations that can facilitate binding  <cit> , and detect the optimal bound state.
currently, protein-protein docking methods are often successful if the experimentally determined partner proteins undergo little conformational changes upon binding.
however, the realistic and computationally efficient treatment of conformational changes, especially of the protein backbone, during docking remains a challenge.
new promising approaches of flexible refinement, ensemble docking and explicit inclusion of flexibility during the entire docking process, have been developed  <cit> .
various strategies have been adopted to score putative protein-protein docked complexes.
most of them utilize atomic level potentials  <cit> , other employ residue level potentials  <cit> , or a combination of both  <cit> .
many of these potentials are knowledge-based  <cit> , others are force field-based  <cit> , or a combination of force field and knowledge-based potentials  <cit> .
the knowledge-based approach aims at deriving empirical interaction parameters from known structures.
the accuracy of these potentials depends on the ability to obtain reliable statistics and to correctly define the reference state.
the relatively small number of solved transient protein-protein complexes, and the small number of contacts at the interface of these complexes result in poor statistics  <cit> .
therefore, many knowledge-based potentials were derived from a dataset of protein folds   <cit>  or from datasets of transient and non-transient interfaces  <cit> .
several works show differences in amino acid composition between protein folds, transient interfaces and non-transient interfaces  <cit> .
ofran and rost  <cit>  performed a comprehensive survey of proteins, identifying six types of interfaces, each one significantly different in its amino acid composition and pairing preference.
these differences impair the ability to obtain reliable statistics and to define accurate reference states, which may well be type specific.
consequently, potentials solely derived from a dataset of transient protein complexes may be more suitable for docking simulations as only the transient complexes may need computational docking predictions.
further more, the lp technique does not require the definition of a reference state therefore the resulting potentials may not be biases toward a specific type of transient complex as shown by liu et al  <cit> .
our purpose is to develop knowledge-based protein-protein docking potentials  tailored to transient complexes that will efficiently discriminate between the structures of native and non-native complexes.
therefore, the potentials are exclusively derived from a dataset of transient complexes.
for each complex, a large number of misdocked  complexes are generated.
a set of potentials is optimized using an lp technique such that for each complex the native state has a lower energy than any of its decoys.
this approach enables us to include in our training set large amount of false positives  in addition to the true positives , thus alleviating the statistical problem of limited data.
in the present study we develop two sets of potentials, coarse-grained and atom-based, using a functional form of the potentials that is not highly sensitive to small structural changes.
this work is an extension of our previous work  <cit>  on designing coarse-grained protein-docking potentials.
however, in our current study we developed two-step potentials, because we are using a larger dataset that cannot be solved using single-step potentials.
the applicability of the new set of dps is illustrated and compared by performing a jackknife test on the zdock  <dig>  decoys set.
outline of the algorithm
the linear programming approach
we proposed to overcome the sampling problem for deriving dps by applying an lp technique to a set of transient complexes.
the utility of this technique in developing both folding and docking potentials previously has been exemplified, by us  <cit>  and by others  <cit> , after the pioneering work of maiorov and crippen  <cit> .
for each complex we generated a large number of putative docked conformations, a few compliant with the native  conformation, and the remaining 'misdocked'.
the basic requirement from the dp is to yield, for each native complex, a free energy lower than that of any misdocked complex.
each of the non-native decoys defines a constraint  on the energy function and we are looking for a set of potentials satisfying these constraints.
the advantages of this approach over the statistical approaches are twofold:  statistical approaches learn from known native states on other native states, i.e.
the  boltzmann statistics is applied to a set of known protein complexes in order to infer the structure of other native complexes.
in the lp approach, we learn from a set of native states and large sets of non-native states  on how a native state should, and should not, look like.
therefore, we have more varied data to derive the empirical potentials.
the lp approach is not sensitive to over-or under-representation or sequence/structure homologies in the training set.
protein models
the side chain based model was previously used by us to develop docking potentials  <cit> .
in short, each amino acid is represented by three interaction sites: the side chain centroid , the amide nitrogen  and the carbonyl oxygen  on the backbone .
this representation involves six types of interactions, s-s, s-o, s-n, o-n, n-n and o-o, the first three being residue-specific.
this results in a set of  <dig>  +  <dig>  +  <dig>  =  <dig> types of interactions.
the atom based model was first introduced by zhang et al.
<cit> .
each side chain is represented by all its heavy atoms which are divided into eighteen types .
this model results in  <dig> different interaction types.
<dig>  atom types were taken from the work of zhang et al.
<cit> .
<dig>  bkbn designate main chain backbone atoms.
<dig>  atom names are according to the pdb nomenclature.
energy function and parameters optimization
the dp design problem is set as follows: assume x to be the coordinate vector representing a protein structure,  the coordinates of the native complex , and {x^a:xbi}i=1n the set of n decoys generated by alternative rotations and translations of monomer b. the demand upon the energy of the native complex is to be lower than that of any of the misdocked conformations, as expressed by the set of n inequalities
{u-u > εi}i=1n
where u represents the potential energy of a complex and εi is a nonnegative constant, set to  <dig>  in the calculations below.
parameter optimization
the total energy is expressed as a linear combination of a parameter set {pl}l=1l
u=∑l=1lslpl
where {s}l=1l designates the set of l .
the functions {pl}l=1l are linear and, therefore, can be optimized using lp algorithms such as simplex or interior-point.
the functions
sl are expressed in terms of a position  dependent step function
hl={1if|rn−rm|≤rcutoffl0 otherwise
that extracts all interacting pairs m and n  within an interaction range rlcutoff characteristic of the specific type  of interaction, i.e.
sl=∑n∈a,m∈bhlδ
hereδ is the kronecker delta, which equals  <dig> if the type of interaction between sites m and n is of type l, and  <dig> otherwise.
for the side chain-based model we adopted rcutoffl values of  <dig>  Å for b-b,  <dig>  Å for b-s, and  <dig>  Å for s-s interactions, whereas for the atom-based model we adopted an rcutoffl value of  <dig>  Å. the total energy for a given decoy i is reduced to the summation over the effective pl, weighted by the numbers nli of contacts of type l occurring in the examined decoy.
the n inequalities in eq  <dig> are thus replaced by
{∑l=1253pl>εi}i=1n
where n^l is the number of contacts of type l in the native complex.
with one step contact potentials, a single cutoff distance defines the threshold below which interactions are considered to be "on", and above which these are considered to be "off".
such potentials may be too simplified to accurately describe the binding energy of interaction proteins.
alternatively, one can design multi-step contact potentials.
we replaced the step function hl byhl,d, where the superscript d refers to different distance ranges .
the energy of the ith decoy is computed from the summation over the effective docking potentials pl,d, weighted by the numbers nl,di of contacts of type l occurring at the distance interval  are thus replaced by
{∑l=1253∑d=1dpl,d>εi}i=1n
determining the cutoff value for distance ranges in multiple step potentials is not trivial.
to overcome this problem we derived a set of single step potentials using a series of successive cutoff values and measuring the correlation coefficients between potentials with adjacent cutoff values.
if the correlation is relatively low another step is created at this point.
in the side chain based-models the following steps were used: r <dig> ≤  <dig>  and  <dig>  < r <dig> ≤ <dig> Å for b-b, r <dig> ≤  <dig>  and  <dig>  <r <dig> ≤ <dig> Å for b-s, and r <dig> ≤  <dig>  and  <dig>  <r <dig> ≤ <dig>  Å for s-s interactions to get a set of  <dig> ×  <dig> =  <dig> parameters.
in the atom-based models the following steps were used: r <dig> ≤  <dig> and 4<r <dig> ≤ <dig> Å to get a set of  <dig> ×  <dig> =  <dig> parameters.
generating a population of putative docked conformations
decoys for the complexes, in their bound state, were generated in the dataset protein-protein docking benchmark  <dig>  by mintseris et al.
<cit>  .
this dataset includes  <dig> protein data bank  entries.
we adopted the docking algorithm of palma et al.
<cit> , taking into consideration only the surface matching criterion to generate docked conformations.
geometric fit was assessed from the number of overlapping surface nodes between the receptor  and the substrate .
the use of this criterion alone in the docking process enables us to generate an unbiased set of decoys, which form many contacts between the receptor and the ligand yet lacking any physicochemical basis.
for any given rotation angle, all possible docking positions were divided into cells of  <dig> ×  <dig> ×  <dig> nodes according to their translation vector, and the best matching position from each cell was preserved.
this modification of the algorithm enabled us to generate a comprehensive and unbiased set of misdocked conformations for each complex.
to reduce the large numbers of decoys generated for each complex, the conformations were further clustered using a hierarchical clustering algorithm applying a cutoff value of root mean square deviation  ≤  <dig> Å and the central conformation of each cluster was preserved.
on average,  <dig>  non-native decoys  were generated for each complex, which summed up to a total of  <dig> , <dig> non-redundant inequalities or constraints on the potentials' parameters.
another set of complexes , listed in table  <dig>  is based on the work of lo conte et al.
<cit> .
this set contains  <dig> complexes for which previously  <cit>  we had generated decoys excluding those that share more that 60% sequence identity  with any of the complexes in dataset  <dig>  dataset  <dig> was derived from the dataset of mintseris and weng  <cit>  as explained in the results section  and decoys were generated as above.
<dig>  a subset of the lo conte et al.
<cit>  dataset.
the table list pdb codes.
<dig>  a subset of the dataset of  <dig> transient complexes compiled by mintseris and weng  <cit>
numerical solution of linear inequalities
inequalities were solved as explained in our previous work  <cit> .
in short, each of the inequalities divides the parameter space into two regions, one allowed , and one forbidden.
a given inequality may affect the solution in three different ways: it may  further restrict the space allowed for the parameter set ,  have no effect on the allowed space, or  impose an impossible condition .
in the latter case no parameter set could satisfy all inequalities, making the inequalities unsolvable.
linear inequalities were solved using the interior point program bpmpd  <cit> .
test decoys
test decoys were used to compare the performance of resulting potentials in the bound and unbound states.
decoy sets containing near-native complexes were generated for the complexes in dataset  <dig> .
this dataset contains the bound and unbound structures of each of the complexes' subunits.
near-native decoys are defined as complexes having an rmsd ≤  <dig>  Å from the native state.
to each complex, docking simulations close to the native state were performed and near-native decoys were generated.
docking was performed using the algorithm of palma et al.
<cit> , taking into consideration only the surface matching criterion to generate docked conformations.
geometric fit was assessed from the number of overlapping surface nodes between the receptor and the substrate.
the docking was performed close to the native state by rotating the substrate by up to 30° and translating by up to 5Å around the native state.
these decoys were clustered using a cutoff rmsd value of ≤  <dig> Å and the central conformation of each cluster was preserved.
subsequently, the preserved decoys were sorted according to their rmsd values from the native state, and a subset of on average  <dig> decoys, uniformly spanning the entire range of near-native decoys, was selected.
this subset was combined with the above set of non-native decoys generated for each complex.
test decoys in the unbound state were generated by superimposing the unbound structures on the decoys generated in the bound state.
by doing so we generated for each complex a set of decoys containing a representative set of both non-native and near-native decoys that enabled us to accurately compare the performance of the potentials in the bound and unbound states irrespective of the docking sampling problem.
jackknife procedure
the jackknife test was performed each time by randomly selecting ten complexes from dataset  <dig>  thereupon, any complex having more than 35% sequence identity with the selected complexes was removed from datasets  <dig>   <dig>  and  <dig>  and the remaining complexes served as the training set.
new potentials were optimized using the reduced training set and tested on the selected ten complexes.
RESULTS
optimizing docking potentials
when optimizing the single-step potentials using constraints derived from dataset  <dig> we encountered a state of infeasibility both for the atom-based and the side chain-based models.
that is, no set of parameters could satisfy all the constraints.
therefore, two-steps potentials were derived.
with these potentials we were able to solve all the constraints that were derived from both dataset  <dig> and dataset  <dig>  the resulting sdps and initial set of adps  are presented in the additional file  <dig> tables s <dig> and s <dig>  respectively.
next, we scanned the dataset of transient complexes compiled by mintseris and weng  <cit>  and identified a subset of  <dig> complexes , for which not all constraints were satisfied using the present adps-i and added them to our training dataset.
thereupon, we solved all constraints derived from datasets:  <dig>   <dig>  and  <dig> to obtain a new set of atomic contact potentials .
).
the potentials are presented in table s <dig> .
in figure  <dig>  panels  and  display the maps of the first and second steps of the adps-ii, respectively, whereas panel  displays the single step decoys as the reference state   <cit>  potentials.
the average contact energy per atom group is displayed in figure  <dig>  the adps-ii properties' and the qualitative differences between the two potentials will be discussed below.
when attempting to solve these constraints using the side chain-based model a state of infeasibility was encountered, which persisted even when using a different set of cutoff values: r <dig> ≤  <dig>  and  <dig>  < r <dig> ≤ <dig> Å for b-b, r <dig> ≤  <dig>  and  <dig>  <r <dig> ≤ <dig>  Å for b-s, and r <dig> ≤  <dig>  and  <dig>  <r <dig> ≤ <dig>  Å for s-s interactions.
optimal tolerance to structural changes
in order to test the optimal tolerance of the potentials to structural changes between the bound and unbound states of the proteins, the decoys  were ranked order exclusively based on the binding-energies.
this test enabled comparing the performances of the potentials in the bound and the unbound states regardless of the sampling.
ranking success for adps-ii is shown in figure  <dig>  in the bound state a near-native decoy is ranked in the top fifteen places in  <dig> out of  <dig> cases, and in the top fifty in  <dig> out of  <dig> cases.
in the top  <dig> places a near-native was identified in all cases.
these results are expected since these complexes were included in the training dataset.
when ranking the unbound decoys using adps-ii, in the top fifteen a near-native complex was ranked in twenty cases, and in the top fifty in twenty six cases.
in the top  <dig> ranked decoys no near-native complex was detected in sixteen cases.
in order to estimate the tolerance of the potentials to structural changes we measured the cα rmsd and heavy-atom rmsd of the interfacial amino acids  between the bound and unbound structure.
interfacial amino acids are defined as those having at least one atom in close proximity  to any atom of the other subunit.
the average irmsd vs. the ranking success of each ranked group is plotted in figure  <dig>  complexes with a near-native structure ranked in the top  <dig> places do not show significant structural change between the bound and unbound states .
complexes that are ranked in the range of 200- <dig> places are characterized by increased side chains displacements.
complexes with a near native structure not ranked in the top  <dig> places are characterized by large side chain and backbone movements .
rescoring zdock <dig>  decoys
next, we subjected the potentials to a jackknife test using the zdock <dig>  decoys   <cit> .
this test set contains  <dig> predictions for each complex in the unbound state generated by the zdock docking program.
each time ten complexes were selected randomly from dataset  <dig> .
complexes homologous to the test complexes were removed from the training dataset and new sets of adps were developed .
the new potentials were used to rank order the decoys of the test set.
this procedure was repeated five times to incorporate a total of fifty proteins.
the average correlation coefficient between adps-ii and the jackknifed potentials is  <dig>  indicating the stability of the trained parameters.
the adps-ii performance was compared with two other potentials, the common atomic contact energies   <cit>  and the newer dars  <cit> , as well with the zdock <dig>  scoring function.
a hit was defined if the irmsd ≤ 4Å, where irmsd values were taken from the zdock <dig>  decoys' supplied tables.
the results have been summarized in table  <dig>  and depicted graphically in figure  <dig>
a-complex pdb code.
b-rank of the lowest energy near-native decoy with irmsd ≤4Å out of  <dig> predictions.
c-jackknife test
-- no near-native decoy identified among the  <dig> best ranked decoys.
a near-native decoy was ranked in the top  <dig> places for twenty seven complexes using the adps-ii, for twenty two using the zdock <dig>  scoring function, for twenty using dars, and for eleven using ace.
the numbers of times the potentials failed to include a near-native case among the top  <dig> raked decoys is: seven for the adps-ii, twelve for the zdock <dig>  scoring function, seventeen for dars, and twenty two for ace.
thus, the adps-ii outperform the other potentials on this set of decoys generated by a docking program different from the one used to generate the training decoys.
the interfacial surface area of each complex was measured using the program naccess  <cit> .
no correlation was found between the size of the interfacial area and the success rate of the adps-ii.
we also compared the capabilities of the adps-ii and zdock <dig>  scoring functions to rank hits in the top ten places.
in order to be consistent with the zdock  <cit>  studies, a hit is defined if the interfacial rmsd is ≤  <dig> Å. for thirty-nine out of the fifty complexes used in the jackknife test such hits do exist in the zdock <dig>  decoy set.
the results are presented in table  <dig>  a near-native decoy was ranked in the top  <dig> places for eight  complexes using the adps-ii, and for six  using the zdock <dig>  scoring function.
a-rank of the lowest energy near-native decoy with irmsd ≤ <dig> Å out of  <dig> predictions.
b-jackknife test
--no near-native decoy identified among the  <dig> best ranked decoys.
discussion and
CONCLUSIONS
when trying to solve the set of inequalities derived from dataset  <dig> using the one step potentials  we encountered the state of infeasibility.
that is, no set of parameters can satisfy all inequalities.
these results indicate that neither the common functional form of atomic contact potentials  <cit>  , nor the side chain-based potentials we designed  <cit> , were adequate to discriminate at a high level between native and non-native docked complexes.
we note that infeasibility was encountered for atomic contact potentials even when using a rcutoff value of  <dig>  Å or  <dig>  Å  <cit> .
in the current study were able to overcome the problem by developing two-step potentials and doubling the number of parameters for the adps and the sdps.
two types of potentials were developed for protein-protein docking: one heavy atom based , and one side-chain based , using training datasets  <dig> and  <dig> for both cases.
when we added dataset  <dig> to the training database we were still able to solve all derived constraints using the atomic model  but not using the side-chain based model, in spite of the fact that the atomic level model has fewer parameters  than the side-chain based one .
these results are in agreement with two other publications  <cit>  showing that the discriminatory powers of atomic level potentials are superior to those of residue level potentials.
we have tested the tolerance of the resulting adps-ii to structural changes between the bound and unbound structures under optimal conditions.
that is when the complexes were included in the training set and the unbound decoys were generated by superimposing the unbound structures on the decoys generated in the bound state.
this comparison enable us the test the effect of structural changes per se.
the results  indicate that structural changes, measured by irmsd of cα/all heavy atoms of more than  <dig> Å/ <dig> Å, impair the ability of the potentials to rank a near native structure in the top  <dig> places.
here we show that the adps-ii outperforms two other contact potentials: the ace and dars.
in addition, the adps-ii outperforms the zdock <dig>  scoring function, which is a combination of desolvation energy  electrostatics and shape complementarity  <cit> .
the shape complementarity criterion alone was shown to have good ranking ability of near native structures  <cit> .
therefore, the performance of our adps-ii in conjunction with shape complementarity and electrostatics is likely to improve.
the colour coded map of the potentials is shown in figure 1a and 1b for adps-ii and in figure 1c for dars.
analysis of the adps-ii reveals that the first step potentials  are stronger in their magnitude than the second step potentials .
the average absolute value of the first and second step potentials is  <dig>  and  <dig>  respectively.
the standard deviation of the adps-ii is  <dig>  and  <dig>  for the first and second steps respectively and  <dig>  for the dars potentials.
thus, separating the potentials into two steps allows us to generate more varied potentials in the first step than in the second step that contains more information.
while, the variability of the single step dars potentials appear to be closer to that of the adps-ii second step.
the average interaction per atom type for the first and second step is displayed in figure  <dig>  overall, the average interactions of the first step potentials are more repulsive than the second step potentials with hydrogen bond interactions being the most attractive ones in the first step potentials.
on the other hand, in the second step potentials hydrophobic interactions are the most attractive ones.
comparing the adps-ii with the dars reveals some interesting differences.
for example, csγ-csγ interaction is the most repulsive adps-ii potential in both ranges whiles the most attractive dars potential.
we emphasize that our potentials were derived exclusively from transient complexes where the formation of cys-cys interactions are undesired, since these complexes are formed and detach from each other in response to specific biological conditions.
the dars potentials assign repulsive to neutral energies to all backbone interactions , except for the cα-cα interaction, which is attractive.
our adps-ii assign a relatively attractive energy to the c-n interactions at the close range and, therefore, captures more accurately h-bond interactions.
interestingly, and unlike the cα-cα interactions, gcα-gcα interactions also are assigned a large attractive energy.
this may reflect the flexibly and compactness of the gly side chains, allowing those to be in closer proximity with each other relative to other side chains.
most computational approaches for protein-protein docking involve a tradeoff between accuracy and computational power.
therefore, a common scheme for protein docking is an initial on grid exhaustive search  <cit>  using relatively simple scoring functions followed by off grid refinement and reranking stage using more elaborate energy functions.
several types of docking potentials were developed for these purposes.
on one hand are the relatively simple one step potentials such as dars and ace, whereas on the other hand are the multiple step potential such as the itscore-pp   <cit> .
the latter potentials were shown to successfully rank the zdock <dig>  decoy , these potentials rank a near-native  in the top  <dig> places in  <dig> % of the cases compared to a success rate of  <dig> % for our adps-ii .
the itscore-pp scoring function is very appropriate for the reranking stage of the docked structure as it requires pre-minimization of the structures to remove clashes and its resolution  is much finer than the 1Å spacing usually used for the search grid.
here we developed an intermediate level of complexity potentials, closer in their functional form to the dars and ace than to the itscore-pp.
these potentials are likely to improve the initial on grid search step.
the main challenges in designing transient protein-protein docking potentials using the boltzmann statistics are:  compiling a none-redundant and representative dataset of protein complexes and,  setting a reference state.
to develop the ace, zhang et al.
used a set of protein folds as a learning dataset and their shuffled structures as a reference state  <cit> .
chuang et al.
used a set of transient and none-transient protein interfaces as a learning set and a large set of misdocked decoys as the reference state  <cit> .
using the lp approach we were able to train our potentials using a dataset of transient complexes only with no need to define a reference state.
the ability to learn solely from transient complexes may add specificity toward this type of complexes.
the above results show that the lp approach is a good alternative to the boltzmann statistics in designing protein docking potentials.
list of abbreviations
: atom-based docking potentials; : docking potential; : false positive; : linear programming; : root mean square deviation; : interfacial rmsd; : side-chain based docking potentials;   : true positive.
authors' contributions
dt preformed the work and wrote the manuscript.
supplementary material
additional file 1
supplementary material.
supplementary tables s <dig>  s <dig> and s3
click here for file
acknowledgements
this research was supported by the israel science foundation .

 proteins interacting with other proteins are known to play key roles in almost all cellular and biological processes such as signalling, metabolism, and trafficking. therefore, understanding protein-protein interactions is fundamental to our ability to comprehend and control cell function. protein complexes can be either permanent or transient. polypeptide chains of permanent complexes remain bound to each other throughout their functional lifetime. an example hereof is the complex between β and γ subunits of hetero-trimeric g-proteins. in transient complexes, the subunits can either form a complex or detach from each other as a result of specific biological conditions. the overall structures of these proteins are stable in the unbound as well as in the bound state. the complex of cyclin and cyclin-dependent protein kinase is an example of such a transient complex. the molecular mechanism of interactions between protein-protein complexes is usually characterized by structure determination methods such as x-ray crystallography, nmr spectroscopy, and cryo-electron microscopy. however, many complexes are too transient to lend themselves to experimental characterization. docking simulations have recently gained importance as possible means for predicting the quaternary structures of protein-protein complexes. yet, these algorithms generate a large number of false-positives and false-negatives, mainly due to the inaccuracy of scoring functions used to evaluate the docked conformations. docking algorithms that reconstruct known complexes using the structures of the proteins in the bound form generally show reasonable levels of success. however, when the structures in the unbound form are used as input, the performance of the same algorithms becomes rather poor. this failure is mainly attributed to the inability of the algorithms and the scoring functions to take account of the conformational changes that accompany complex formation. in an effort to overcome these deficiencies, one usually resorts to additional biological data, which assist in singling out the native conformation among putative docked conformations. however, in the absence of prior knowledge regarding the binding site, exhaustive searches of up to different positions and orientations are performed for the substrate, followed by a refinement and clustering stage. thus, it is difficult to take into account all the possible changes in internal conformations that can facilitate binding, and detect the optimal bound state. currently, protein-protein docking methods are often successful if the experimentally determined partner proteins undergo little conformational changes upon binding. however, the realistic and computationally efficient treatment of conformational changes, especially of the protein backbone, during docking remains a challenge. new promising approaches of flexible refinement, ensemble docking and explicit inclusion of flexibility during the entire docking process, have been developed. various strategies have been adopted to score putative protein-protein docked complexes. most of them utilize atomic level potentials, other employ residue level potentials, or a combination of both. many of these potentials are knowledge-based, others are force field-based, or a combination of force field and knowledge-based potentials. the knowledge-based approach aims at deriving empirical interaction parameters from known structures. the accuracy of these potentials depends on the ability to obtain reliable statistics and to correctly define the reference state. the relatively small number of solved transient protein-protein complexes, and the small number of contacts at the interface of these complexes result in poor statistics. therefore, many knowledge-based potentials were derived from a dataset of protein folds or from datasets of transient and non-transient interfaces. several works show differences in amino acid composition between protein folds, transient interfaces and non-transient interfaces. ofran and rost performed a comprehensive survey of proteins, identifying six types of interfaces, each one significantly different in its amino acid composition and pairing preference. these differences impair the ability to obtain reliable statistics and to define accurate reference states, which may well be type specific. consequently, potentials solely derived from a dataset of transient protein complexes may be more suitable for docking simulations as only the transient complexes may need computational docking predictions. further more, the lp technique does not require the definition of a reference state therefore the resulting potentials may not be biases toward a specific type of transient complex as shown by liu et al. our purpose is to develop knowledge-based protein-protein docking potentials tailored to transient complexes that will efficiently discriminate between the structures of native and non-native complexes. therefore, the potentials are exclusively derived from a dataset of transient complexes. for each complex, a large number of misdocked complexes are generated. a set of potentials is optimized using an lp technique such that for each complex the native state has a lower energy than any of its decoys. this approach enables us to include in our training set large amount of false positives in addition to the true positives, thus alleviating the statistical problem of limited data. in the present study we develop two sets of potentials, coarse-grained and atom-based, using a functional form of the potentials that is not highly sensitive to small structural changes. this work is an extension of our previous work on designing coarse-grained protein-docking potentials. however, in our current study we developed two-step potentials, because we are using a larger dataset that cannot be solved using single-step potentials. the applicability of the new set of dps is illustrated and compared by performing a jackknife test on the zdock decoys set. outline of the algorithm the linear programming approach we proposed to overcome the sampling problem for deriving dps by applying an lp technique to a set of transient complexes. the utility of this technique in developing both folding and docking potentials previously has been exemplified, by us and by others, after the pioneering work of maiorov and crippen. for each complex we generated a large number of putative docked conformations, a few compliant with the native conformation, and the remaining'misdocked'. the basic requirement from the dp is to yield, for each native complex, a free energy lower than that of any misdocked complex. each of the non-native decoys defines a constraint on the energy function and we are looking for a set of potentials satisfying these constraints. the advantages of this approach over the statistical approaches are twofold: statistical approaches learn from known native states on other native states, i.e. the boltzmann statistics is applied to a set of known protein complexes in order to infer the structure of other native complexes. in the lp approach, we learn from a set of native states and large sets of non-native states on how a native state should, and should not, look like. therefore, we have more varied data to derive the empirical potentials. the lp approach is not sensitive to over-or under-representation or sequence/structure homologies in the training set. protein models the side chain based model was previously used by us to develop docking potentials. in short, each amino acid is represented by three interaction sites: the side chain centroid, the amide nitrogen and the carbonyl oxygen on the backbone. this representation involves six types of interactions, s-s, s-o, s-n, o-n, n-n and o-o, the first three being residue-specific. this results in a set of + + = types of interactions. the atom based model was first introduced by zhang et al.. each side chain is represented by all its heavy atoms which are divided into eighteen types. this model results in different interaction types. atom types were taken from the work of zhang et al.. bkbn designate main chain backbone atoms. atom names are according to the pdb nomenclature. energy function and parameters optimization the dp design problem is set as follows: assume x to be the coordinate vector representing a protein structure, the coordinates of the native complex, and {x^a:xbi}i=1n the set of n decoys generated by alternative rotations and translations of monomer b. the demand upon the energy of the native complex is to be lower than that of any of the misdocked conformations, as expressed by the set of n inequalities {u-u > εi}i=1n where u represents the potential energy of a complex and εi is a nonnegative constant, set to in the calculations below. parameter optimization the total energy is expressed as a linear combination of a parameter set {pl}l=1l u=∑l=1lslpl where {s}l=1l designates the set of l. the functions {pl}l=1l are linear and, therefore, can be optimized using lp algorithms such as simplex or interior-point. the functions sl are expressed in terms of a position dependent step function hl={1if|rn−rm|≤rcutoffl0 otherwise that extracts all interacting pairs m and n within an interaction range rlcutoff characteristic of the specific type of interaction, i.e. sl=∑n∈a,m∈bhlδ hereδ is the kronecker delta, which equals if the type of interaction between sites m and n is of type l, and otherwise. for the side chain-based model we adopted rcutoffl values of Å for b-b, Å for b-s, and Å for s-s interactions, whereas for the atom-based model we adopted an rcutoffl value of Å. the total energy for a given decoy i is reduced to the summation over the effective pl, weighted by the numbers nli of contacts of type l occurring in the examined decoy. the n inequalities in eq are thus replaced by {∑l=1253pl>εi}i=1n where n^l is the number of contacts of type l in the native complex. with one step contact potentials, a single cutoff distance defines the threshold below which interactions are considered to be "on", and above which these are considered to be "off". such potentials may be too simplified to accurately describe the binding energy of interaction proteins. alternatively, one can design multi-step contact potentials. we replaced the step function hl byhl,d, where the superscript d refers to different distance ranges. the energy of the ith decoy is computed from the summation over the effective docking potentials pl,d, weighted by the numbers nl,di of contacts of type l occurring at the distance interval are thus replaced by {∑l=1253∑d=1dpl,d>εi}i=1n determining the cutoff value for distance ranges in multiple step potentials is not trivial. to overcome this problem we derived a set of single step potentials using a series of successive cutoff values and measuring the correlation coefficients between potentials with adjacent cutoff values. if the correlation is relatively low another step is created at this point. in the side chain based-models the following steps were used: r ≤ and < r ≤ Å for b-b, r ≤ and <r ≤ Å for b-s, and r ≤ and <r ≤ Å for s-s interactions to get a set of × = parameters. in the atom-based models the following steps were used: r ≤ and 4<r ≤ Å to get a set of × = parameters. generating a population of putative docked conformations decoys for the complexes, in their bound state, were generated in the dataset protein-protein docking benchmark by mintseris et al.. this dataset includes protein data bank entries. we adopted the docking algorithm of palma et al., taking into consideration only the surface matching criterion to generate docked conformations. geometric fit was assessed from the number of overlapping surface nodes between the receptor and the substrate. the use of this criterion alone in the docking process enables us to generate an unbiased set of decoys, which form many contacts between the receptor and the ligand yet lacking any physicochemical basis. for any given rotation angle, all possible docking positions were divided into cells of × × nodes according to their translation vector, and the best matching position from each cell was preserved. this modification of the algorithm enabled us to generate a comprehensive and unbiased set of misdocked conformations for each complex. to reduce the large numbers of decoys generated for each complex, the conformations were further clustered using a hierarchical clustering algorithm applying a cutoff value of root mean square deviation ≤ Å and the central conformation of each cluster was preserved. on average, non-native decoys were generated for each complex, which summed up to a total of, non-redundant inequalities or constraints on the potentials' parameters. another set of complexes, listed in table is based on the work of lo conte et al.. this set contains complexes for which previously we had generated decoys excluding those that share more that 60% sequence identity with any of the complexes in dataset dataset was derived from the dataset of mintseris and weng as explained in the results section and decoys were generated as above. a subset of the lo conte et al. dataset. the table list pdb codes. a subset of the dataset of transient complexes compiled by mintseris and weng numerical solution of linear inequalities inequalities were solved as explained in our previous work. in short, each of the inequalities divides the parameter space into two regions, one allowed, and one forbidden. a given inequality may affect the solution in three different ways: it may further restrict the space allowed for the parameter set, have no effect on the allowed space, or impose an impossible condition. in the latter case no parameter set could satisfy all inequalities, making the inequalities unsolvable. linear inequalities were solved using the interior point program bpmpd. test decoys test decoys were used to compare the performance of resulting potentials in the bound and unbound states. decoy sets containing near-native complexes were generated for the complexes in dataset. this dataset contains the bound and unbound structures of each of the complexes' subunits. near-native decoys are defined as complexes having an rmsd ≤ Å from the native state. to each complex, docking simulations close to the native state were performed and near-native decoys were generated. docking was performed using the algorithm of palma et al., taking into consideration only the surface matching criterion to generate docked conformations. geometric fit was assessed from the number of overlapping surface nodes between the receptor and the substrate. the docking was performed close to the native state by rotating the substrate by up to 30° and translating by up to 5Å around the native state. these decoys were clustered using a cutoff rmsd value of ≤ Å and the central conformation of each cluster was preserved. subsequently, the preserved decoys were sorted according to their rmsd values from the native state, and a subset of on average decoys, uniformly spanning the entire range of near-native decoys, was selected. this subset was combined with the above set of non-native decoys generated for each complex. test decoys in the unbound state were generated by superimposing the unbound structures on the decoys generated in the bound state. by doing so we generated for each complex a set of decoys containing a representative set of both non-native and near-native decoys that enabled us to accurately compare the performance of the potentials in the bound and unbound states irrespective of the docking sampling problem. jackknife procedure the jackknife test was performed each time by randomly selecting ten complexes from dataset thereupon, any complex having more than 35% sequence identity with the selected complexes was removed from datasets and and the remaining complexes served as the training set. new potentials were optimized using the reduced training set and tested on the selected ten complexes.  optimizing docking potentials when optimizing the single-step potentials using constraints derived from dataset we encountered a state of infeasibility both for the atom-based and the side chain-based models. that is, no set of parameters could satisfy all the constraints. therefore, two-steps potentials were derived. with these potentials we were able to solve all the constraints that were derived from both dataset and dataset the resulting sdps and initial set of adps are presented in the additional file tables s and s respectively. next, we scanned the dataset of transient complexes compiled by mintseris and weng and identified a subset of complexes, for which not all constraints were satisfied using the present adps-i and added them to our training dataset. thereupon, we solved all constraints derived from datasets: and to obtain a new set of atomic contact potentials. ). the potentials are presented in table s. in figure panels and display the maps of the first and second steps of the adps-ii, respectively, whereas panel displays the single step decoys as the reference state potentials. the average contact energy per atom group is displayed in figure the adps-ii properties' and the qualitative differences between the two potentials will be discussed below. when attempting to solve these constraints using the side chain-based model a state of infeasibility was encountered, which persisted even when using a different set of cutoff values: r ≤ and < r ≤ Å for b-b, r ≤ and <r ≤ Å for b-s, and r ≤ and <r ≤ Å for s-s interactions. optimal tolerance to structural changes in order to test the optimal tolerance of the potentials to structural changes between the bound and unbound states of the proteins, the decoys were ranked order exclusively based on the binding-energies. this test enabled comparing the performances of the potentials in the bound and the unbound states regardless of the sampling. ranking success for adps-ii is shown in figure in the bound state a near-native decoy is ranked in the top fifteen places in out of cases, and in the top fifty in out of cases. in the top places a near-native was identified in all cases. these results are expected since these complexes were included in the training dataset. when ranking the unbound decoys using adps-ii, in the top fifteen a near-native complex was ranked in twenty cases, and in the top fifty in twenty six cases. in the top ranked decoys no near-native complex was detected in sixteen cases. in order to estimate the tolerance of the potentials to structural changes we measured the cα rmsd and heavy-atom rmsd of the interfacial amino acids between the bound and unbound structure. interfacial amino acids are defined as those having at least one atom in close proximity to any atom of the other subunit. the average irmsd vs. the ranking success of each ranked group is plotted in figure complexes with a near-native structure ranked in the top places do not show significant structural change between the bound and unbound states. complexes that are ranked in the range of 200- places are characterized by increased side chains displacements. complexes with a near native structure not ranked in the top places are characterized by large side chain and backbone movements. rescoring zdock decoys next, we subjected the potentials to a jackknife test using the zdock decoys. this test set contains predictions for each complex in the unbound state generated by the zdock docking program. each time ten complexes were selected randomly from dataset. complexes homologous to the test complexes were removed from the training dataset and new sets of adps were developed. the new potentials were used to rank order the decoys of the test set. this procedure was repeated five times to incorporate a total of fifty proteins. the average correlation coefficient between adps-ii and the jackknifed potentials is indicating the stability of the trained parameters. the adps-ii performance was compared with two other potentials, the common atomic contact energies and the newer dars, as well with the zdock scoring function. a hit was defined if the irmsd ≤ 4Å, where irmsd values were taken from the zdock decoys' supplied tables. the results have been summarized in table and depicted graphically in figure a-complex pdb code. b-rank of the lowest energy near-native decoy with irmsd ≤4Å out of predictions. c-jackknife test -- no near-native decoy identified among the best ranked decoys. a near-native decoy was ranked in the top places for twenty seven complexes using the adps-ii, for twenty two using the zdock scoring function, for twenty using dars, and for eleven using ace. the numbers of times the potentials failed to include a near-native case among the top raked decoys is: seven for the adps-ii, twelve for the zdock scoring function, seventeen for dars, and twenty two for ace. thus, the adps-ii outperform the other potentials on this set of decoys generated by a docking program different from the one used to generate the training decoys. the interfacial surface area of each complex was measured using the program naccess. no correlation was found between the size of the interfacial area and the success rate of the adps-ii. we also compared the capabilities of the adps-ii and zdock scoring functions to rank hits in the top ten places. in order to be consistent with the zdock studies, a hit is defined if the interfacial rmsd is ≤ Å. for thirty-nine out of the fifty complexes used in the jackknife test such hits do exist in the zdock decoy set. the results are presented in table a near-native decoy was ranked in the top places for eight complexes using the adps-ii, and for six using the zdock scoring function. a-rank of the lowest energy near-native decoy with irmsd ≤ Å out of predictions. b-jackknife test --no near-native decoy identified among the best ranked decoys. discussion and  when trying to solve the set of inequalities derived from dataset using the one step potentials we encountered the state of infeasibility. that is, no set of parameters can satisfy all inequalities. these results indicate that neither the common functional form of atomic contact potentials, nor the side chain-based potentials we designed, were adequate to discriminate at a high level between native and non-native docked complexes. we note that infeasibility was encountered for atomic contact potentials even when using a rcutoff value of Å or Å. in the current study were able to overcome the problem by developing two-step potentials and doubling the number of parameters for the adps and the sdps. two types of potentials were developed for protein-protein docking: one heavy atom based, and one side-chain based, using training datasets and for both cases. when we added dataset to the training database we were still able to solve all derived constraints using the atomic model but not using the side-chain based model, in spite of the fact that the atomic level model has fewer parameters than the side-chain based one. these results are in agreement with two other publications showing that the discriminatory powers of atomic level potentials are superior to those of residue level potentials. we have tested the tolerance of the resulting adps-ii to structural changes between the bound and unbound structures under optimal conditions. that is when the complexes were included in the training set and the unbound decoys were generated by superimposing the unbound structures on the decoys generated in the bound state. this comparison enable us the test the effect of structural changes per se. the results indicate that structural changes, measured by irmsd of cα/all heavy atoms of more than Å/ Å, impair the ability of the potentials to rank a near native structure in the top places. here we show that the adps-ii outperforms two other contact potentials: the ace and dars. in addition, the adps-ii outperforms the zdock scoring function, which is a combination of desolvation energy electrostatics and shape complementarity. the shape complementarity criterion alone was shown to have good ranking ability of near native structures. therefore, the performance of our adps-ii in conjunction with shape complementarity and electrostatics is likely to improve. the colour coded map of the potentials is shown in figure 1a and 1b for adps-ii and in figure 1c for dars. analysis of the adps-ii reveals that the first step potentials are stronger in their magnitude than the second step potentials. the average absolute value of the first and second step potentials is and respectively. the standard deviation of the adps-ii is and for the first and second steps respectively and for the dars potentials. thus, separating the potentials into two steps allows us to generate more varied potentials in the first step than in the second step that contains more information. while, the variability of the single step dars potentials appear to be closer to that of the adps-ii second step. the average interaction per atom type for the first and second step is displayed in figure overall, the average interactions of the first step potentials are more repulsive than the second step potentials with hydrogen bond interactions being the most attractive ones in the first step potentials. on the other hand, in the second step potentials hydrophobic interactions are the most attractive ones. comparing the adps-ii with the dars reveals some interesting differences. for example, csγ-csγ interaction is the most repulsive adps-ii potential in both ranges whiles the most attractive dars potential. we emphasize that our potentials were derived exclusively from transient complexes where the formation of cys-cys interactions are undesired, since these complexes are formed and detach from each other in response to specific biological conditions. the dars potentials assign repulsive to neutral energies to all backbone interactions, except for the cα-cα interaction, which is attractive. our adps-ii assign a relatively attractive energy to the c-n interactions at the close range and, therefore, captures more accurately h-bond interactions. interestingly, and unlike the cα-cα interactions, gcα-gcα interactions also are assigned a large attractive energy. this may reflect the flexibly and compactness of the gly side chains, allowing those to be in closer proximity with each other relative to other side chains. most computational approaches for protein-protein docking involve a tradeoff between accuracy and computational power. therefore, a common scheme for protein docking is an initial on grid exhaustive search using relatively simple scoring functions followed by off grid refinement and reranking stage using more elaborate energy functions. several types of docking potentials were developed for these purposes. on one hand are the relatively simple one step potentials such as dars and ace, whereas on the other hand are the multiple step potential such as the itscore-pp. the latter potentials were shown to successfully rank the zdock decoy, these potentials rank a near-native in the top places in % of the cases compared to a success rate of % for our adps-ii. the itscore-pp scoring function is very appropriate for the reranking stage of the docked structure as it requires pre-minimization of the structures to remove clashes and its resolution is much finer than the 1Å spacing usually used for the search grid. here we developed an intermediate level of complexity potentials, closer in their functional form to the dars and ace than to the itscore-pp. these potentials are likely to improve the initial on grid search step. the main challenges in designing transient protein-protein docking potentials using the boltzmann statistics are: compiling a none-redundant and representative dataset of protein complexes and, setting a reference state. to develop the ace, zhang et al. used a set of protein folds as a learning dataset and their shuffled structures as a reference state. chuang et al. used a set of transient and none-transient protein interfaces as a learning set and a large set of misdocked decoys as the reference state. using the lp approach we were able to train our potentials using a dataset of transient complexes only with no need to define a reference state. the ability to learn solely from transient complexes may add specificity toward this type of complexes. the above results show that the lp approach is a good alternative to the boltzmann statistics in designing protein docking potentials. list of abbreviations: atom-based docking potentials;: docking potential;: false positive;: linear programming;: root mean square deviation;: interfacial rmsd;: side-chain based docking potentials;: true positive. authors' contributions dt preformed the work and wrote the manuscript. supplementary material additional file 1 supplementary material. supplementary tables s s and s3 click here for file acknowledgements this research was supported by the israel science foundation.
our results show that, using the lp approach, we were able to train our potentials using a dataset of transient complexes only the newly developed potentials outperform three other known potentials in this test. to achieve this we considered a set of transient complexes and generated a large set of putative docked structures, based on a shape complementarity criterion, for each complex. the demand on the potentials was to yield, for the native structure, a potential energy lower than those of any of the non-native structures. the performance of the atom-based potentials was evaluated by a jackknife test on a set of complexes taken from the zdock decoys set.
 protein-protein docking is a challenging computational problem in functional genomics, particularly when one or both proteins undergo conformational change upon binding. the major challenge is to define a scoring function soft enough to tolerate these changes and specific enough to distinguish between near-native and "misdocked" conformations.  using a linear programming technique, we developed two types of potentials: side chain-based and heavy atom-based. to achieve this we considered a set of transient complexes and generated a large set of putative docked structures, based on a shape complementarity criterion, for each complex. the demand on the potentials was to yield, for the native structure, a potential energy lower than those of any of the non-native structures. we show that the heavy atom-based potentials were able to comply with this requirement but not the side chain-based one. thus, despite the smaller number of parameters, the capability of heavy atom-based potentials to discriminate between native and "misdocked" conformations is improved relative to those of the side chain-based potentials. the performance of the atom-based potentials was evaluated by a jackknife test on a set of complexes taken from the zdock decoys set.  our results show that, using the lp approach, we were able to train our potentials using a dataset of transient complexes only the newly developed potentials outperform three other known potentials in this test. 