David W Ritchie
email: dave.ritchie@inria.fr

Calculating and scoring high quality multiple flexible protein structure alignments

Motivation: Calculating multiple protein structure alignments (MSAs) is important for understanding functional and evolutionary relationships between protein families, and for modeling protein structures by homology. While incorporating backbone flexibility promises to circumvent many of the limitations of rigid MSA algorithms, very few flexible MSA algorithms exist today. This article describes several novel improvements to the Kpax algorithm which allow high quality flexible MSAs to be calculated. This article also introduces a new Gaussian-based MSA quality measure called "M-score", which circumvents the pitfalls of RMSD-based quality measures. Results: As well as calculating flexible MSAs, the new version of Kpax can also score MSAs from other aligners and from previously aligned reference datasets. Results are presented for a large-scale evaluation of the Homstrad, SABmark, and SISY benchmark sets using Kpax and Matt as examples of state-of-the-art flexible aligners and 3DCOMB as an example of a state-of-the-art rigid aligner. These results demonstrate the utility of the M-score as a measure of MSA quality and show that high quality MSAs may be achieved when structural flexibility is properly taken into account.

Introduction

Calculating multiple structure alignments (MSAs) of proteins is an important step towards understanding the evolutionary and functional relationships between different proteins and protein families [START_REF] Sierk | Déjà vu all over again: Finding and analyzing protein structure similarities[END_REF][START_REF] Hasegawa | Advances and pitfalls of protein structure alignment[END_REF]. Calculating high quality structural alignments can be crucial for the accurate modelling of threedimensional (3D) protein structures by homology [START_REF] Madhusudhan | Alignment of multiple protein structures based on sequence and structure features[END_REF][START_REF] Braberg | SALIGN: a web server for alignment of multiple proteins and structures[END_REF], and can improve the generation of MSAbased sequence profiles for protein fold prediction [START_REF] Pei | PROMALS3D: a tool for multiple protein sequence and structure alignments[END_REF][START_REF] Ghouzam | Improving protein fold recognition with hybrid profiles combining sequence and structure evolution[END_REF]. Many algorithms for calculating protein structure alignments (SAs) have been described and reviewed in recent years [START_REF] Sierk | Déjà vu all over again: Finding and analyzing protein structure similarities[END_REF][START_REF] Kolodny | Comprehensive evaluation of protein structure alignment methods: Scoring by geometric measures[END_REF][START_REF] Hasegawa | Advances and pitfalls of protein structure alignment[END_REF][START_REF] Ma | Algorithms, applications, and challenges of protein structure alignment[END_REF], and several of these can calculate MSAs. Some examples of multiple structural aligners are SSAP [START_REF] Taylor | Multiple protein structure alignment[END_REF], MUSTA [START_REF] Leibowitz | MUSTA -a general, efficient, automated method for multiple structure alignment and detection of common motifs: Application to proteins[END_REF], CE-MC [START_REF] Guda | CE-MC: a multiple protein structure alignment server[END_REF], MultiProt (Shatsky et al., 2004b), POSA [START_REF] Ye | Multiple flexible structure alignment using partial order graphs[END_REF], MAMMOTH-multi [START_REF] Lupyan | A new progressive-iterative algorithm for multiple structure alignment[END_REF], Vorolign [START_REF] Birzele | Vorolign -fast structural alignment using Voronoi contacts[END_REF] MUSTANG [START_REF] Konagurthu | MUSTANG: a multiple structural alignment algorithm[END_REF], PROMALS3D [START_REF] Pei | PROMALS3D: a tool for multiple protein sequence and structure alignments[END_REF], BLOMAPS [START_REF] Wang | Fast multiple alignment of protein structures using conformational letter blocks[END_REF], SALIGN [START_REF] Madhusudhan | Alignment of multiple protein structures based on sequence and structure features[END_REF], MISTRAL [START_REF] Micheletti | MISTRAL: a tool for energy-mbased multiple structure alignment of proteins[END_REF], 3DCOMB [START_REF] Wang | Alignment of distantly related protein structures: algorithm, bound and implications to homology modeling[END_REF], MAPSCI [START_REF] Ilinkin | Multiple structure alignment and consensus identification for proteins[END_REF], Smolign [START_REF] Sun | Smolign: a spatial motifs based protein multiple structural alignment method[END_REF], msTALI [START_REF] Shealy | Multiple structure alignment with msTALI[END_REF], and mulPBA [START_REF] Joseph | Progressive structure-based alignment of homologous proteins: Adopting sequence comparison strategies[END_REF].

Most of the above approaches first use dynamic programming (DP) to align pairs of structures, and they then progressively merge the pair-wise alignments to make a MSA. A few approaches subsequently optimise the resulting MSA using distance-based (e.g. CE-MC and MAMMOTH-multi) or energy-based (e.g. MISTRAL) minimisation techniques. However, when aligning many structures, the progressive alignment approaches can become computationally expensive because they require an initial set of all-against-all pair-wise alignments to build a guide tree, which means that the computational cost scales as O(N 2 ) in the number of structures, N . In order to avoid such computational costs, the MAPSCI and mulPBA algorithms build a multiple alignment with respect to a given "pivot" structure, which leads to a computational cost of "kpax_2016" -2016/10/10 -11:17 -page 2 -#2

i i i i i i i i 2 David W. Ritchie *
O(N ). But this requires prior knowledge of which structure should serve as the pivot. An effective alternative to conventional "row-first" multiple alignment algorithms is 3DCOMB [START_REF] Wang | Alignment of distantly related protein structures: algorithm, bound and implications to homology modeling[END_REF], which uses a "column-first" approach to identify "highly similar fragment blocks" to guide the calculation of MSAs. While 3DCOMB has been shown to give excellent MSAs [START_REF] Wang | Alignment of distantly related protein structures: algorithm, bound and implications to homology modeling[END_REF], the complexity of its column-first algorithm also appears to scale as O(N 2 ) in the number of structures. Because it is difficult to optimise simultaneously the number of aligned residues and their Cα root mean squared deviation (RMSD) [START_REF] Zemla | LGA a method for finding 3D similarities in protein structures[END_REF], different structure aligners make different number/RMSD trade-offs which can lead to producing locally different alignments, even within proteins that share common motifs. Thus, current SA algorithms can introduce errors or inconsistencies even in pair-wise alignments [START_REF] Sadowski | Evolutionary innacuracy of pairwise structural alignments[END_REF], and these errors can accumulate considerably in MSA algorithms that progressively merge multiple pair-wise alignments. Furthermore, such effects can confound efforts to compare different MSA algorithms against reference SAs, because the programs used to calculate the reference SAs are themselves subject to the same number/RMSD trade-off problem.

Incorporating backbone flexibility offers one way to tackle such problems. Although some flexible SA algorithms have been described, such as FlexProt [START_REF] Shatsky | Flexible protein alignment and hinge detection[END_REF](Shatsky et al., , 2004a)), FATCAT [START_REF] Ye | Flexible structure alignment by chained aligned fragment pairs allowing twists[END_REF], POSA [START_REF] Ye | Multiple flexible structure alignment using partial order graphs[END_REF][START_REF] Li | POSA: a user-driven, interactive multiple protein structure alignment server[END_REF], FlexSnap [START_REF] Salem | FlexSnap: flexible non-sequential protein structure alignment[END_REF], Matt [START_REF] Menke | Matt: local flexibility aids protein multiple structure alignment[END_REF], and RAPIDO [START_REF] Mosca | Alignment of protein structures in the presence of domain motions[END_REF], still only a very few algorithms can flexibly align multiple structures. Furthermore, other "flexible" aligners such as Vorolign [START_REF] Birzele | Vorolign -fast structural alignment using Voronoi contacts[END_REF] and GR-Align (Malod-Dognin and Pržulj, 2014) can align pairs of structures that cannot be completely superposed rigidly, but they do not perform superpositions.

To my knowledge, the only aligners that can calculate multiple superpositions flexibly are POSA [START_REF] Ye | Multiple flexible structure alignment using partial order graphs[END_REF], Matt [START_REF] Menke | Matt: local flexibility aids protein multiple structure alignment[END_REF], and Smoalign [START_REF] Sun | Smolign: a spatial motifs based protein multiple structural alignment method[END_REF]. POSA essentially extends the pair-wise FATCAT approach to MSAs. Both POSA and Smolign calculate partial-order alignment graphs, i.e. lists of locally aligned segments which are common to all structures but which are not necessarily sequential in the primary protein sequences. In contrast, both Matt and Kpax calculate strictly sequential MSAs which do not allow non-sequential permutations of residues. Unfortunately, POSA is publicly available only as a web server [START_REF] Li | POSA: a user-driven, interactive multiple protein structure alignment server[END_REF], which makes it infeasible for large-scale studies, and Smolign is no longer available for download.

While the original Kpax pair-wise alignment algorithm is sufficiently fast to be able to search large 3D structural domain databases in just a few seconds, the first version was limited to performing only rigid-body superpositions of entire protein structures. This article describes several extensions to the earlier Kpax algorithm [START_REF] Ritchie | Fast protein structure alignment using Gaussian overlap scoring of backbone peptide fragment similarity[END_REF] which allow pairs of proteins of arbitrary size (including multi-domain proteins) to be aligned flexibly, and which allow multiple proteins to be assembled into a MSA using a pivot-based approach.

It seems obvious that a pre-requisite for calculating high quality MSAs is to begin with an algorithm that can calculate high quality pairwise SAs. However, there does not yet exist an accepted standard for assessing the quality of pair-wise alignments or MSAs. 3DCOMB [START_REF] Wang | Alignment of distantly related protein structures: algorithm, bound and implications to homology modeling[END_REF], calculates its MSA quality measure from the pair-wise TMscore [START_REF] Zhang | TM-align: a protein structure alignment algorithm based on TM-score[END_REF] taken over all fully aligned columns. On the other hand, for columns with an occupancy of 70% or more, mulPBA uses a quality measure based on a weighted average of the numbers of residues whose pair-wise distances fall within certain given distance ranges [START_REF] Joseph | Progressive structure-based alignment of homologous proteins: Adopting sequence comparison strategies[END_REF]. This may be considered as an extension of the Local-Global Alignment (LGA) scoring method which calculates a weighted average over several number/RMSD bins of a pairwise alignment [START_REF] Zemla | LGA a method for finding 3D similarities in protein structures[END_REF]. More recently, the STOVCA program [START_REF] Slater | Towards the development of standardized methods for comparison, ranking and evaluation of structure alignments[END_REF] has been proposed to compare the quality of structural superpositions calculated by other pair-wise aligners, but this program is limited to evaluating pair-wise superpositions using its own internal pair-wise alignment algorithm. [START_REF] Collier | A new statistical framework to assess structural alignment quality using information compression[END_REF] proposed a method for measuring SA quality based on information compression, but this approach is also limited to pair-wise alignments. Thus, there exists a need for a simple and general way to compare the quality of both the alignments and the corresponding superpositions of different MSA algorithms. To try to address the problem of comparing different SA algorithms, this article proposes a novel atomic Gaussian based MSA scoring function, called "M-score", which circumvents the number/RMSD trade-off problem.

In previous studies, [START_REF] Wang | Alignment of distantly related protein structures: algorithm, bound and implications to homology modeling[END_REF] showed that 3DCOMB gives highly competitive rigid MSAs compared to MAPSCI, MAMOTH, Matt, BLOMAPS, MUSTANG, and MultiProt when applied to the large Homstrad [START_REF] Mizuguchi | HOMSTRAD: a database of protein structure alignments for homologous protein families[END_REF][START_REF] Stebbings | HOMSTRAD: recent developments of the homologous protein structure alignment database[END_REF] and SABmark [START_REF] Van Walle | SABmark-a benchmark for sequence alignment that covers the entire known fold space[END_REF] benchmark sets, and [START_REF] Berbalk | Accuracy analysis of multiple structure alignments[END_REF] showed that Matt gave the best overall performance on a sub-set of their challenging SISY benchmark. Therefore, a similar large-scale comparison using these datasets is presented here to compare Kpax with 3DCOMB as an example of a state-of-the-art rigid multiple aligner, and with Matt as the best available flexible multiple aligner.

Methods

Normalised Atomic Gaussian Overlap SAs

Kpax pair-wise SAs are calculated using two Cα-based scoring functions, both of which are based on Gaussian density distributions. In the main scoring function, each Cα atom is represented as a normalised 3D atomic Gaussian shape-density distribution

φ(x) = 1 πσ 2 3/4 e -r 2 /2σ 2 , (1)
where r = |x| is the distance from the nucleus and σ is the atomic radius. I put σ = 1.4 Å, which is slightly less than the van der Waals radius of 1.5 Å often used in molecular mechanics force fields. A similar 3D Gaussian representation has been shown to provided an effective model of molecular shape for small-molecule shape comparison and virtual drug screening [START_REF] Grant | A fast method of molecular shape comparison: A simple application of a Gaussian description of molecular shape[END_REF].

It is well known that Gaussian functions have special mathematical properties. For example, it has long been known that the overlap between a product of two Gaussian functions is also a Gaussian (see Equation 10in [START_REF] Boys | Electronic wave functions I. A general method of calculation for the stationary states of any molecular system[END_REF]). In particular, it can be shown that the 3D overlap between two of the above normalised 3D Gaussian functions on atoms i and j located at x i and x j and separated by a distance R ij = |x i -x j | is given by

G ij = φ(x i )φ(x j )dx = e -R 2 ij /4σ 2 . ( 2
)
This overlap expression has the nice property that it gives a value of one unit when two Cα atoms are perfectly superposed because e 0 = 1. Furthermore, when the distance between atoms i and j exceeds the sum, S, of their Cα radii, we have S = 2σ and e -R 2 /S 2 tends to zero very rapidly when R > S. Thus, while other Gaussian scoring functions have been used to provide probabilistic or "elastic" scoring functions in previous SA algorithms, notably STAMP [START_REF] Russell | Multiple sequence alignment from tertiary structure comparison: Assignment of global and residue confidence levels[END_REF] and DALI [START_REF] Holm | Protein structure comparison by alignment of distance matrices[END_REF], it can be seen that the form of Equation 2 is rather special, and is certainly not arbitrary. Using the above physical model of protein shape, a SA between two proteins may be scored by summing pairs of Gaussian overlaps to give a global alignment score, or "G-score"

G = i,j µ ij G ij , (3)
where the coefficient µ ij = 1 when residue i of the first protein is aligned to residue j of the second protein, and zero otherwise.

"kpax_2016" -2016/10/10 -11:17 -page 3 -#3

i i i i i i i i
Kpax structure alignment 3 However, the above scoring function requires that a pair-wise alignment and 3D superposition are already known, which is initially not the case. Therefore, in order to identify candidate aligned fragment pairs (AFPs) for least-squares fitting and superposition, Kpax calculates an initial pose-invariant similarity score for all possible residue pairs using a scoring function based on a sliding window of 7 residues on each chain (i.e. 6 Cα pairs, around the pair of interest), and using products of normalised Gaussian overlaps at each position, which I call "K-scores", as described previously [START_REF] Ritchie | Fast protein structure alignment using Gaussian overlap scoring of backbone peptide fragment similarity[END_REF]. Products, rather than sums, of local Gaussian overlap scores are used to calculate an initial alignment in order that a pair-wise K-score tends rapidly to zero whenever any of its component distances are large. Thus, only those residue pairs that share highly similar local environments will have large K-scores. Applying dynamic programming (DP) to a matrix of such K-scores, produces an initial seed alignment, or "K-alignment", from which runs of 4 or more pairs (AFPs) having large K-scores may be identified and used as 3D superposition seeds.

The K-alignment of a pair of similar protein domains typically contains up to about 10 AFPs, whereas a pair of large domains with many similar secondary structure elements may give up to around 50 AFPs. For each candidate 3D superposition seed proposed by a K-alignment AFP, an allagainst all matrix of Cα G-scores (the "G-matrix") is calculated, from which a global alignment may be calculated by DP with no gap penalties. The residue pairs extracted from the resulting "G-alignment" are then used to define a new 3D least-squares fitting transformation, and the cycle is iterated until convergence (typically just 3 or 4 iterations are required). The G-alignment that gives the best global G-score is then retained as the the "optimal" global alignment. From my experiments with the examples presented here, using any value of σ that differs much from 1.4 Å generally gives poorer alignments (details not shown). Hence, putting σ = 1.4Å appears to be numerically "optimal" in practice, as well as theory.

Detecting Conserved Motifs using Tiled DP

Although Kpax uses the very fast and robust quaternion minimisation least squares fitting method of [START_REF] Liu | Fast determination of the optimal rotational matrix for macromolecular superpositions[END_REF], and only a relatively small number of least-squares fits are necessary, the overall computational cost of the above procedure scales as O(N 2 ) in the number of residues to be aligned. This makes it is relatively expensive to align large multi-domain proteins by iterating on multiple candidate G-matrices. Furthermore, it is possible that the optimal DP trace-back from the K-matrix might not contain any 3D fitting seeds that correspond to an optimal 3D G-alignment. Therefore, when one or both proteins has more that 120 residues, a series of small "tiled" SAs is calculated, in which each protein is divided into overlapping strips of 72 residues, and a local SA is calculated between each strip of one protein and the whole of the other protein using DP on the tiles of the K and G matrices. For example, when aligning two proteins each having 144 residues, a total of 6 tiled alignments are calculated using the residue ranges (1:71,1:144), (37:108,1:144), (73:144, 1:144), (1:144,1:71), (1:144, 37:108), (1:144,73:144). The 3D transformation from the tile that gives the greatest G-score is retained and used in a final DP iteration on the the full G-matrix. This approach is easily generalised to the case where the number of residues is not a whole multiple of the tile size. For example, a protein of 200 residues would be divided into strips of (1:64), (65:136), (137:200), (1:28), (29:100), (101:172), and (173:200).

For large proteins, tiling reduces the computational cost from O(F N 2 ) to O(T f n 2 ) where F is the number of AFPs in the global alignment, f is the average number of fragments in a tile, T is the number of tiles of size n, and where n < N and T << N . Furthermore, tiling allows the possibility of finding good local 3D seeds in the K-alignment which might not appear in the (globally optimal) trace-back path of the full K-matrix. In other words, small conserved motifs may be detected in a tiled alignment whereas they might be completely missed in a global alignment.

Since it is not known in advance whether a given pair of proteins will align well globally or whether a tiled alignment should be applied, Kpax first attempts a global alignment and if the global alignment M-score (see below) is less than 0.5, it applies a tiled alignment in order to find the best seed for a global alignment.

Flexible Alignments using Multiple DP Boxes

Calculating flexible alignments may be considered as a natural extension of the above tiled approach, in which the tiles no longer have a pre-determined size. For Kpax flexible alignments, all residues are initially labelled as belonging to a default segment (Segment 0), and after the optimal rigid superposition has been calculated as described above, all pairs of residues that fall within a given distance threshold (3 Å, by default) are assigned to a new segment (Segment 1). The remaining Segment 0 residues of each structure are then considered as candidates for flexible refinement. More specifically, the first DP problem corresponds to the entire DP matrix, and consecutive runs of four or more residue pairs belonging to Segment 0 are identified and superposed by least-squares fitting. The initial DP matrix is then divided into three or more "boxes", where the first box spans the first and last residue pairs of segment 1, thus leaving one N-terminal and one Cterminal box of the DP matrix, each of which contains candidate residues for flexible refinement. It is also possible for one or more DP boxes to be found within the interior of a superposed segment (i.e. segments do not necessarily contain continuous runs of superposed residues). In any case, each new DP box may be considered as a new SA sub-problem, which may be treated exactly as described above. Thus, by maintaining a list of DP boxes to be processed on a stack, the flexible alignment procedure continues recursively, and it terminates when no more DP boxes remain in the stack.

The evolution of the DP matrix for the first few iterations of this procedure is illustrated schematically in Supplementary Figure 1. Each residue in the moving structure is transformed in 3D space at most two times: once during the rigid body superposition, and possibly once again if it belongs to a DP box that produces a superposition segment. Note, when superposing a moving segment to a rigid template structure, the coordinates of the neighbouring residues relative to the moving segment are taken into account whenever possible. Thus, the final flexibly aligned structure has few obvious gaps or tears in its backbone atom coordinates.

Pivot-First Method of Calculating MSAs

In Kpax, MSAs are calculated using a rigid pivot structure onto which the remaining structures are either rigidly or flexibly aligned and superposed. By default, in both rigid and flexible MSAs, each structure is used in turn as the pivot structure, and the MSA that yields the best overall M-score is retained as the pivot. Thus, for N structures, there will be N -1 alignments between the selected pivot and each of the remaining structures. A MSA may then be induced from the N -1 pair-wise alignments by introducing additional gaps as necessary to re-align the N -2 gapped pivot structures with the pivot of the first alignment. While this induced MSA (IMSA) is not necessarily optimal, it provides a fast way to identify highly conserved columns in the alignment. Any column of the IMSA which is occupied with 60% or more residues is noted as a highly conserved column.

The IMSA is used to provide a look-up table for building a final MSA. Kpax first sorts the N -1 other structures by their G-score similarity to the pivot. It then takes the most similar structure as a seed alignment, and it builds the final MSA by incrementally aligning the remaining N -2 structures in a "pile-up" manner onto this seed. Each in-coming conserved residue is biased to appear close to its position in the initial IMSA by calculating column-wise G-scores for a window of ± 16 residues around "kpax_2016" -2016/10/10 -11:17 -page 4 -#4

i i i i i i i i 4 David W. Ritchie *
the corresponding IMSA position. Thus, the alignment of conserved regions may be further optimised by DP, and in-coming non-conserved residues have the opportunity to be re-aligned against the preceding structures. Consequently, the Kpax MSA algorithm may be considered as a "pivot-first" hybrid between the row-first and column-first approaches.

Although it is conceivable that one could use a guide-tree to align multiple pairs of flexibly aligned structures, this could lead to many false AFPs and many tears in the moving structures. The above pile-up approach helps to ensure that the number of deformations to the flexible structures is small, and that each non-rigid structure is deformed only once. A pile-up approach also allows for easy parallel execution on multiple processors.

What is Wrong with RMSD?

There are two fundamental problems with RMSD as a quality measure. Firstly, it is non-linear in the number of distances concerned. Secondly, RMSDs can hide many sins. Consider one column of a MSA of 10 structures, in which 9 of the Cα atoms superpose perfectly, but one is displaced by 5 Å and therefore probably belongs to an adjacent column (recall that the distance between consecutive Cα backbone atoms is normally 3.8 Å). Suppose we calculate the RMSD with respect the average column coordinate (which for convenience can be assumed to be at the coordinate origin, and where the displacement is along one of the principal axes). This gives ((9 * 0.5 2 + 1 * 4.5 2 )/10) 1/2 = 1.50 Å. Now consider one column of a MSA of 20 structures, with one Cα again displaced by 5 Å. This gives ((19 * 0.25 2 + 1 * 4.75 2 )/20) 1/2 = 1.08 Å. In other words, increasing the number of rows has spread the displacement over more samples and hence apparently reduced the RMSD by about 30%.

But one should not forget that RMSD has units of "distance" / ("number-of-pairs") 1/2 . So it would be more precise to label each RMSD with the number of pairs involved, e.g. RMSD 10 = 1.50 Å, and RMSD 20 = 1.08 Å. This reminds us that it is incorrect to compare or combine RMSD 10 and RMSD 20 directly because they were calculated using different distance scales. Consequently, a fairer way to compare these two cases is to extract the average deviations as RMSD 10 / √ 10 = 0.47 Å and RMSD 20 / √ 20 = 0.24 Å. But this obviously amounts to not using RMSD in the first place. If we wish to score and compare SAs reliably, it seems much more desirable to call a gap for the displaced residue, and to award scores of 9 out of 10 (say) for the first example and 19 out of 20 for second. This is what the M-score aims to achieve.

MSA Quality Assessment Using the M-score

The M-score extends the pair-wise Cα Gaussian overlap model to MSAs of arbitrary dimensions. Firstly, let N be the total number of rows (i.e. structures) in a MSA, C be the total number of columns, T be the total number of residues over all structures, and L be the number of residues in the longest structure. Thus, an "ideal" MSA may be considered to have at least one row of L residues and the remaining (T -L) residues from the other structures would each align perfectly with one of the residues of the longest structure. This implies than many of the MSA columns must contain 2 or more completely overlapping Cα atoms.

A simple way to estimate this mutual overlap is to calculate an average coordinate for each column, and then to calculate the degree of overlap between this average coordinate and each of the column members. Using these ideas, the column overlap score C j with respect to the average 3D coordinate, c j , for the j th column may be defined as:

C j = N i=1 µ ij e -(x ij -c j ) 2 /4σ 2 , ( 4
)
where µ ij is either 1 or 0 according to whether the position at row i and column j contains a residue or a gap, respectively, and x ij is the coordinate of the Cα atom at row i and columns j. Note that one wrongly placed Cα atom will reduce this score by at most one unit. On the other hand, if all of the Cα coordinates in column j coincide with the average coordinate, then C j will be numerically equal to the number of residues in that column. Thus, a global overlap score may be calculated for the MSA as

M = (T -L) + C - C j=1 max(C j , 1), (5)
where the max function ensures that at least one unit is subtracted for each column, and thus deals with columns that contain only a few very poorly superposed Cα atoms, which may occur at the termini of an alignment or in the loop regions around an AFP. If all of the Cα atoms in each column coincide completely with the average column coordinate, then C = L and C j=1 max(C j , 1) = T , and hence M = 0. On the other hand, in the worst possible case of no aligned residue pairs, we have c j = x ij for all i and hence C j = 1 for all j. Therefore, for any MSA, M will have the range

0 ≤ M ≤ (T -L).
However, because it is more intuitive to let large scores correspond to good alignments, the above expression may be scaled and shifted using

M = 1 -M /(T -L) to give M = C j=1 max(C j , 1) -C (T -L) . (6)
With this scoring function, a perfect multiple alignment of any number of structures (e.g. between any protein and one or more identical copies of itself) would give a score of unity. Furthermore, for pair-wise alignments, the factor (T -L) will be equal to the number of residues in the smaller structure, and a score of unity will correspond to the complete overlap of all Cα atoms of the smaller structure with some subset of the Cα atoms of the larger structure. If desired, the M-score may be multipled by (T -L)σ 3 to give an absolute quality measure in Å 3 units.

It may be noted that the M-score contains no gap penalties and has no distance bins or distance cut-offs. Its only potentially adjustable parameter is the Cα radius, which in practice should be fixed at 1.4 Å. Furthermore, by construction, the only way to maximise M is to have as few columns as possible of which some subset must each contain two or more maximally overlapping Cα atoms. This means that the M-score smoothly circumvents the number/RMSD trade-off problem.

In cases where two or more proteins share only a small motif, or where they have long unaligned N or C termini, it can be useful to restrict the M-score to a given block of the MSA. For example, to score the alignment of a motif between columns J and K, it is straight-forward to calculate

M JK = K j=J max(C j , 1) -C JK (T JK -L JK ) , (7)
where C JK and T JK denote the number of columns and total number of residues between columns J and K, respectively, and L JK is the longest chain within that block.

Comparing Pair-wise SA Algorithms

As a first but important comparison, Table 1 shows the flexible alignment results calculated by Kpax, Matt, and FATCAT for for the 10 "difficult" pairs of structures from [START_REF] Fischer | Assessing the performance of fold recognition methods by means of a comprehensive benchmark[END_REF]. Supplementary Table 1 shows the results for the corresponding rigid alignment calculated by Kpax, Matt, FATCAT, 3DCOMB, TMalign, MUSTANG, CE, and DALI. Comparison of the numbers of aligned residues and the RMSDs with the corresponding M-score in these tables shows that the M-score captures rather well the aim to reward a high number of aligned pairs with low RMSD. For example, in the first row of residue identities (I), the alignment RMSD in Å units (R), and the M-score (M). The Matt alignments were made using Matt-1.00 [START_REF] Menke | Matt: local flexibility aids protein multiple structure alignment[END_REF]. The FATCAT alignments were made using jFATCAT [START_REF] Prlić | Pre-calculated protein structure alignments at the RCSB PDB website[END_REF]. † The A chain is used in all PDB structures except for 3hhr (chain B), 1cew (chain I), 1nsb (chain B), and 1bge (chain B).

Flexibly Aligning Multi-Domain Proteins

To my knowledge, the FlexProt algorithm of [START_REF] Shatsky | Flexible protein alignment and hinge detection[END_REF]Shatsky et al. ( , 2004a) ) is the first structure aligner to have been developed to handle structural flexibility in large multi-domain proteins which cannot be rigidly superposed. According to the M-scores in 1 and2, in each case either Kpax or Matt produces a better alignment according to the M-score.

Two Flexible MSA Examples

Figure 1 shows two examples of the rigid and flexible MSA superpositions produced by Kpax using the SABmark superfamily group104 ("sup-group104") and twilight-zone group202 ("twi-group202"). In this figure, the structures are coloured using a simple but informative "rainbow" colour-coding scheme, in which the residues in consecutive columns of a MSA are coloured in consecutive colours of the rainbow. These examples demonstrate the benefit of using each member structure as a candidate rigid pivot structure. For sup-group104 (consisting of the four SCOP structures: d1lgpa_, d1dmza_, d1g6ga_, and d1gxca_), the best rigid MSA pivot is found to be d1dmza_, giving M=0.626 and 68 core columns, while the worst pivot (d1lgpa_) gives M=0.541 and 64 core columns. However, for a flexible MSA with d1lgpa_ as pivot, the above quality measures increase to M=0.779 and 97 core columns. On the other hand, the best pivot changes to d1g6ga_, giving M=0.784 and 98 core columns. This shows that the best pivot for a flexible MSA cannot necessarily be predicted from rigid pair-wise SAs.

While Figure 1 suggests that calculating a MSA for sup-group104 should be relatively easy, it also shows that aligning all 3 members of twi-group202 (SCOP codes d1l2wi_, d1k8kd1, and d1jyoa_) is much more challenging. In this case, the rigid Kpax MSA gives only 14 core columns and M=0.285. In contrast, the Kpax flexible MSA (lower right in Figure 1) gives 51 core columns and M=0.750. This demonstrates the improvement that may be achieved on going from a rigid to a flexible MSA.

The Homstrad, SABmark, and SISY Benchmarks

The version of Homstrad used here (October 2014) consists of 1,032 families of homologous protein families containing a total of 3,458 structures. Each Homstrad family contains from 2 to 25 members, where no two members of any family have more than 90% sequence identity [START_REF] Stebbings | HOMSTRAD: recent developments of the homologous protein structure alignment database[END_REF]. The Homstrad SAs were calculated using STAMP, MNYFIT [START_REF] Sutcliffe | Knowledge based modelling of homologous proteins, part i: three-dimensional frameworks derived from the simultaneous superposition of multiple structures[END_REF], and COMPARER [START_REF] Sali | Definition of general topological equivalence in protein structures: A procedure involving comparison of properties and relationships through simulated annealing and dynamic programming[END_REF]. The SABmark dataset (version 1.65) contains multiple alignments for 3,320 structures belonging of 425 superfamilies with low sequence identity and 1,740 structures belonging to 209 "twighlight zone" families with very low sequence identity. The SABmark dataset has only sequence-based pair-wise reference alignments. It has been shown that these alignments contain mutual inconsistencies [START_REF] Edgar | Quality measures for protein alignment benchmarks[END_REF], hence it impossible to reconstruct MSAs from the reference pairs. Consequently, it is not possible to calculate M-scores or core columns for SABmark.

The SISY datasets were collected by [START_REF] Berbalk | Accuracy analysis of multiple structure alignments[END_REF], who calculated the SISY reference alignments using CE and DALI. SISY-Families consists of 106 families of single chain SCOP domains with less than 40% sequence identity. SISY-pairs consists of 130 protein pairs with low sequence similarity drawn from the SISYPHUS database [START_REF] Mayr | Comparative analysis of protein structure alignments[END_REF]. SISY-RIPC (repetitions, insertions, permutations, and conformational variability) consists of 23 pairs which Berbalk et al. report to be difficult to align. Many of the SISY reference alignments are given not as full alignments but as "trimmed" alignments, consisting of a small list of highly conserved reference pairs or MSA columns. Several of the examples in the SISY-Pairs and SISY-RIPC datasets contain nonsequential alignments or small motifs within large unrelated folds which can be difficult to detect.

As noted by Berbalk et al., some aligners can fail when presented with a PDB file that contains more than one structure, or when a PDB file contains SEQRES records that differ from the ATOM records. Therefore, all SISY structures involving more than one chain per PDB file were excluded, and the remaining PDB files were cleaned by running Kpax with the "-preprocess" option to change HETERO amino acid records to ATOM records, and to remove all other HETERO and SEQRES records. After cleaning, 3DCOMB was found to fail with a memory fault on 2 of the SISY examples and Matt gave a memory fault or failed to terminate on 6 examples. Excluding the above cases left a total of 77 SISY-Families, 111 SISY-Pairs, and 23 SISY-RIPC examples. Applying the same cleaning procedure to the Homstrad dataset detected 7 PDB files with duplicate or corrupted PDB chains, and these were also excluded. Additionally, 5 SABmark superfamily groups and 4 SABmark twighlight zone groups were eliminated because their reference alignments contained no pairs.

MSA Benchmark Quality Comparison

Table 3 presents results for the alignments calculated by Kpax, Matt, and 3DCOMB for the above benchmark datasets. The Homstrad dataset provides superposed structures and a PIR-format MSA alignment file. Hence it is possible to calculate several common quality measures such as the number of fully aligned core columns ("#Core" in Table 3), the number of aligned pairs ("#Pair"), the number of aligned residue identities ("#Identity"), and the average RMSD over all aligned pairs. Additionally, since such global measures do not take into account local alignment quality, it is useful calculate the number of badly aligned pairs, and the number of superposed fragment pairs ("#SFP"), where a SFP is taken to be a run of 3 or more aligned pairs in which the distance between each pair does not exceed 3.5Å, and where any pair-wise distance that exceeds 8.0 Å is counted as a distance violation ("#Viol"). The columns labeled "M 1C " and "M JK " show the M-scores calculated for the full and trimmed MSAs, respectively, where M 1C ranges over all MSA columns and M JK is calculated between the columns containing the first and last pair of aligned residues. The values shown for the number of pairs, identities, distance violations, SFPs, and RMSD values are normalised by the number of pairs of structures involved (i.e. for N structures, the normalisation factor is N (N -1)/2). Supplementary Figures 2-10 show more detailed "boxwhisker" plots [START_REF] Spitzer | BoxPlotR: a web tool for generation of box plots[END_REF] of the results in Table 3. Spread-sheets of numerical results are provided as Supplementary Files.

For this table, the RMSD values for the SABmark reference alignments were calculated using the pair-wise reference alignments to superpose the corresponding structures. The RMSD values for the SISY reference alignments were calculated by iteratively superposing the SISY structures onto an average coordinate for each reference column. Because several of the SISY reference alignments contain non-sequential SAs that cannot be fully aligned by sequential aligners, these datasets are grouped separately in Table 3 (denoted as "seq" and "nonseq", accordingly) because these sub-sets show rather different statistics.

Several revealing observations may be made from Table 3. Firstly, it can seen that some of the reference alignments contain many distance violations, notable the Homstrad and SISY-RIPC-seq datasets. As observed previously for multiple sequence alignments [START_REF] Edgar | Quality measures for protein alignment benchmarks[END_REF], this means that structural reference alignments sometimes contain significant errors, and it is therefore dangerous to treat even well curated reference alignments as a gold standard against with different algorithms should be compared. In particular, for Homstrad, it can be seen that 3DCOMB, Matt, and Kpax all achieve higher M-scores with lower numbers of cores columns and residue identities than the Homstrad reference alignments, while the Homstrad reference alignments give a higher average RMSD and more distance violations than these three methods. For the SISY-RIPC-seq reference set, most of the distance violations come from just 3 pairs, with the pair d2bbma_/4clna_ having 114 distance violations, which is clearly incorrect, and the pairs d1ggga_/d1wdna_ and d1l5ba_/d1l5ea_ having 25 and 19 violations, respectively, which make them of doubtful quality.

The large numbers of core columns calculated by 3DCOMB is impressive, and no doubt reflects 3DCOMB's strategy of aiming to produce many filled columns [START_REF] Wang | Alignment of distantly related protein structures: algorithm, bound and implications to homology modeling[END_REF], but this is often at the expense of higher RMSDs and more distance violations than Kpax and Matt. In contrast, Matt-flex gives amongst the lowest RMSDs and lowest numbers of distance violations on all datasets. However, with the exception of some of the trimmed SISY reference alignments, this is always at the expense of finding fewer aligned pairs and residue identities than the other alignments.

It can be seen from Table 3 that Kpax-flex produces the largest average M-score while also finding the greatest number of residue identities in all datasets except for Homstrad and the small SISY-RIPC-nonseq group. Kpax-flex also finds the largest average trimmed M-score for the Homstrad and SABmark datasets, and it is interesting to note that average trimmed M-score for Kpax-flex is larger than that of the reference alignments for the SISY-RIPC-seq group. Furthermore, Kpax, Matt, and 3DCOMB all find more aligned pairs, residue identities, and AFPs in the SISY datasets than are present in the SISY reference alignments. While the trimmed SISY reference alignments have high M JK values, which is to be expected, the rigid alignments calculated by Kpax and Matt for the full structures Apart from the surprisingly high number of residue identities in the Homstrad reference alignments, the Kpax flexible alignments contain the greatest number of aligned identities in all datasets except for the SISY-Families-nonseq and SISY-RIPC-nonseq groups. Additionally, the final column of Table 3 shows that the flexible Kpax alignments contain the largest average number of SFPs in all datasets. While it might be argued that flexible Kpax is "over-fitting" the 3D structures, it should be noted that the Kpax rigid alignments also contain the most SFPs for all datasets except for the SISY-RIPC-seq and SISY-Pairs-nonseq groups. Because it would seem to be very difficult to match runs of 3 or more low RMSD pairs purely by chance, these results strongly suggest that the Kpax alignments contain more conserved structural regions than the 3DCOMB and Matt alignments, and indeed more than many of the reference alignments. More certainly, the results in Table 3 show that high M-scores are highly correlated with high numbers of SFPs, core residues, aligned pairs, and low RMSDs. Therefore, these results demonstrate the utility of the M-score as a measure of MSA quality. Nonetheless, since current high quality alignments can still contain distance violations, a useful future extension could be to add a term to penalise such occurrences.

It is worth mentioning that Kpax is considerably faster than 3DCOMB and Matt. Aligning the 77 SISY-Families MSAs using Kpax takes around 13 minutes for rigid alignments and 17 minutes for flexible alignments on a workstation with eight E5410 2.3 GHz processors and 12 Gb memory. In comparison, 3DCOMB takes 49 minutes for its rigid alignments and Matt takes 200 minutes for both rigid and flexible alignments.

Conclusions

This article has described a novel approach for calculating flexible protein SAs by recursively dividing a DP scoring matrix into multiple boxes. Each box of the DP matrix corresponds to a rigid SA sub-problem. This allows flexible SAs to be calculated almost as easily as rigid alignments. This approach is easily extended to calculating flexible MSAs by flexibly aligning each mobile structure onto an automatically selected rigid pivot. This article has also described the M-score, a novel Gaussian-based measure of MSA quality. The results presented for three large benchmark datasets show that high M-scores correspond to SAs with high numbers of aligned pairs with low RMSD. These results show that Kpax, Matt, and 3DCOMB generally produce SAs with more aligned pairs, residue identities, and aligned fragments than are present in the Homstrad, SABmark, and SISY reference alignments. It has also been shown that the rigid SAs from Matt and 3DCOMB, as well as several reference datasets, can contain a high proportion of pair-wise distances that exceed 8 Å, indicating the presence of significant alignment errors. While the flexible SAs produced by Matt are free of distance violations, these SAs are often much shorter than those produced by Kpax. In contrast, the full-length Kpax flexible SAs contain more residue identities and structurally aligned fragments than the Matt SAs, while still showing few distance violations and reduced RMSDs. Overall, these results show that high quality SAs may be achieved when structural flexibility is properly taken into account, and they support the utility of the M-score as a measure of MSA quality.

Fig. 1 .

 1 Fig. 1. Kpax rigid (left) and flexible (right) MSAs for SABmark superfamily group104 (top row) and twighlight zone group202 (bottom row), drawn using scripts generated by Kpax. The two groups of superpositions are not shown at the same scale.

Table 1

 1 Table1shows that Kpax produces better flexible alignments than both Matt and FATCAT, while Supplementary Table1shows that Kpax and 3DCOMB produce amongst the best rigid alignments.

	(1fxi/1ubq), Kpax aligns 64 pairs with a RMSD

Table 1 .

 1 Flexible SA comparison between Kpax, Matt and FATCAT for the 10 "difficult" pairs from[START_REF] Fischer | Assessing the performance of fold recognition methods by means of a comprehensive benchmark[END_REF].

	PDBs †	Kpax-flex		Matt-flex		FATCAT-flex
		S A R I	M	A R I	M	S A R I	M
	1fxi/1ubq	4 64 1.6 4 0.734	49 1.0 5 0.607	1 63 3.0 9 0.492
	1ten/3hhr	1 86 1.5 17 0.854	73 0.7 17 0.797	1 87 1.9 16 0.794
	3hla/2rhe	4 81 1.6 11 0.701	58 1.0 9 0.551	3 91 3.2 10 0.497
	2aza/1paz	3 88 1.5 18 0.652	61 0.8 11 0.489	1 89 2.8 15 0.472
	1cew/1mol 3 81 1.0 16 0.811	69 0.6 13 0.715	1 83 2.4 14 0.641
	1cid/2rhe	3 100 1.9 14 0.739	77 1.2 11 0.622	1 103 2.8 13 0.581
	1crl/1ede 11 238 2.2 22 0.606	153 0.9 9 0.470	7 280 3.1 21 0.490
	2sim/1nsb 12 286 1.8 25 0.632	218 1.1 22 0.534	1 305 3.1 30 0.429
	1bge/2gmf 6 105 1.8 18 0.740	59 0.3 11 0.486	1 116 3.0 18 0.488
	1tie/4fgf	4 117 1.8 11 0.800	69 0.6 7 0.543	1 120 3.2 11 0.589
	Mean	5 125 1.7 16 0.727	89 0.8 12 0.581	2 134 2.8 16 0.547

* * Listed are the numbers of segments (S; not available for Matt), aligned pairs (A), aligned

Table 2 ,

 2 Kpax often produces better flexible alignments for the 18 structure pairs of Shatsky et al.

	(average Kpax M-score: 0.794) than both Matt (M=0.732) and FATCAT
	(M=0.689). Although the average numbers of aligned residues and residue
	identities are similar for Kpax and FATCAT (245 and 129 for Kpax
	compared to 254 and 132 for FATCAT, respectively), the lower average
	RMSD for Kpax (1.4 Å compared to 2.3 Å) accounts for its greater average
	M-score. As in Table 1, Matt produces some very low RMSD alignments,
	but this is at the expense of finding fewer aligned pairs and pair-wise
	identities. It may be noted that Matt finds sub-optimal alignments in at
	least three cases (2bbm/1cll, 2bbm/1top, and 1mcp/1tcr) in which the
	M-scores and numbers of aligned residues are significantly lower than
	those for Kpax and FATCAT. In contrast, Kpax finds only one sub-optimal
	alignment (1a21/1hwg) relative to the Matt and FATCAT alignments.
	Although FATCAT finds respectable alignments for all of the examples
	in Tables

Table 2 .

 2 Flexible SA comparison between Kpax, Matt, and FATCAT for 18 pairs from[START_REF] Shatsky | Flexible protein alignment and hinge detection[END_REF]Shatsky et al. ( , 2004a)).

	PDBs †	Kpax-flex		Matt-flex	FATCAT-flex
		S A R I	M	A R I	M S A R I	M
	1wdn/1ggg 2 219 0.7 214 0.970 214 0.4 208 0.963 1 220 1.0 220 0.948
	2bbm/1cll 2 140 1.8 131 0.801	65 0.4 31 0.448 2 144 2.1 141 0.774
	2bbm/1top 7 143 1.6 65 0.837	65 0.4 24 0.436 4 147 2.3 77 0.743
	1ake/2ak3 4 206 1.2 83 0.880 201 0.5 81 0.925 3 206 1.7 84 0.819
	2ak3/1uke 6 184 1.3 49 0.867 168 0.5 52 0.859 1 188 2.5 54 0.654
	1mcp/4fab 2 217 1.2 176 0.921 213 0.6 174 0.952 2 217 1.3 176 0.908
	1mcp/1tcr 5 212 1.6 53 0.841	95 0.7 28 0.419 2 213 2.2 52 0.760
	1lfh/1lfg	2 687 0.8 684 0.960 675 0.2 674 0.973 3 689 1.0 687 0.945
	1tfd/1lfh	3 289 1.2 182 0.915 229 0.5 108 0.769 3 291 1.4 187 0.889
	1qf6/1adj 7 349 1.8 59 0.696 299 0.8 57 0.685 2 353 2.4 61 0.589
	2clr/3fru	8 264 1.3 70 0.888 233 0.5 66 0.853 1 263 3.0 68 0.563
	1fmk/1qcf 9 431 1.0 256 0.933 418 0.4 259 0.945 1 434 2.3 265 0.737
	1fmk/1tki 8 240 1.7 39 0.640 209 0.6 42 0.638 1 242 3.2 42 0.438
	1a21/1hwg 2 89 1.9 18 0.392 120 0.9 17 0.599 2 167 2.9 20 0.490
	1ntr/1dc7 8 122 1.3 122 0.893	90 0.8 90 0.695 1 122 3.1 122 0.565
	1dce/1a17 3 103 2.3 9 0.507 117 0.6 16 0.722 4 150 2.7 14 0.631
	2nac/1psd 4 306 1.3 80 0.736 281 0.6 76 0.737 1 311 3.2 76 0.441
	1chm/1xgs 3 206 1.4 29 0.622 171 0.7 26 0.565 1 214 2.4 28 0.509
	Mean	5 245 1.4 129 0.794 215 0.6 113 0.732 2 254 2.3 132 0.689
	Time/s	7.6		55.6	92.5

* * Column heading are the same as in Table

1

. † The A chain in all PDB structures is used here except for PDB code 1mcp (chain L).

Table 3 .

 3 Kpax, Matt, and 3DCOMB performance on the Homstrad, SABmark, and SISY datasets. Pairs datasets give lower RMSDs than the RMSDs of the SISY reference alignments. Taken together, these observations indicate that a high proportion of the SISY reference alignments contain sub-optimal alignments compared to those calculated by Kpax, Matt, and 3DCOMB. The box-whisker plots in Supplementary Figures 5-10 further support this observation.

	"kpax_2016" -2016/10/10 -11:17 -page 7 -#7

Funding

This work was funded by Inria Nancy -Grand Est and by the Agence Nationale de la Recherche, grant reference ANR-11-MONU-006-02.


TITLE: Kbdock -Searching and Organising the Structural Space of Protein-Protein Interactions AUTHORS: Marie-Dominique Devignes, Malika Smaïl-Tabbone, David Ritchie TEASER: Big data is a recurring problem in structural bioinformatics where even just one experimentally determined protein structure can contain several different interacting protein domains and often involves many tens of thousands of 3D atomic coordinates. If we consider all protein structures that have ever been solved, the immense structural space of proteinprotein interactions needs to be organised systematically in order to make sense of the many functional and evolutionary relationships that exist between different protein families and their interactions. This article describes some new developments in Kbdock, a knowledgebased approach for classifying and annotating protein interactions at the protein domain level.

Protein-protein interactions (PPIs) are fundamental biophysical interactions. Our understanding of many biological processes relies on the 3D-modelling of PPIs, and there is a growing need to be able to classify and analyse the structural repertoire of known PPIs using computational modelling and analysis techniques. Today, there exists over 105,000 experimentally determined 3D structures of proteins and protein complexes in the publicly available Protein Data Bank (PDB; http://www.rcsb.org/pdb/home/home.do). Each entry is composed of the 3D coordinates several thousands (sometimes even millions) of atoms belonging to one or more linear chains of amino-acid residues. By analysing multiple protein structures and their amino-acid sequences, it can be seen that many protein chains are modular, being composed of one or more structural domains. Thus, protein domains may be considered as "knowledge units" because they represent abstract classes which group together proteins with similar amino-acid sub-sequences and similar biological properties. We started to develop Kbdock (Kb for "knowledge-based") to address the problem of PPI modelling and classification at the domain level. In fact, as in many other complex scientific domains, Big Data approaches in the life sciences can only become viable by explicitly considering and making use of prior knowledge. Essentially, Kbdock is a dedicated relational database which combines the Pfam domain classification (http://pfam.xfam.org/) with coordinate data from the PDB to analyse and model domain-domain interactions (DDIs) in 3D space. The detection of a DDI instance relies on the computation of the spatial distance between the surface residues of the domain instances found in each PDB entry. In the latest release of Kbdock, after duplicate or near-duplicate interactions are removed, a total of 5,139 distinct non-redundant DDIs involving two different domains have been identified from nearly 240,000 DDI instances extracted from the PDB. As illustrated in Fig 1A , the Kbdock resource can be queried by Pfam domain to visualise the DDI network involving this domain. Otherwise, Kbdock can return the list of DDI instances corresponding to a given pair of domains contained in two interacting proteins, even when there is little or no similarity with the query proteins. Moreover, calculating 3D super-positions of all DDI instances involving a given domain (Fig 1B) enabled us to perform a spatial clustering of all DDIs involving that domain, thereby identifying a discrete number of so-called "domain family binding sites" (DFBSs) [START_REF] Ghoorah | Spatial clustering of protein binding sites for template based protein docking[END_REF] on the domain of interest. This gives us a new and original kind of knowledge unit, which we find to be essential when studying the structure and specificity of protein binding sites on a large scale [START_REF] Ghoorah | A structure-based classification and analysis of protein domain family binding sites and their interactions[END_REF]. The notion of DFBSs also led us to propose a case-based reasoning approach to the problem of how to retrieve the best available template for modeling protein-protein "docking" interactions. The Kbdock project is run as a collaboration between the Capsid and Orpailleur teams at the Loria/Inria research center in Nancy. It is funded and supported by Inria, the CNRS, and the University of Lorraine, as well as specific ANR ("Agence Nationale pour la Recherche") grants. The Kbdock program is available through its on-line interface (http://kbdock.loria.fr/). It may also be queried programmatically by expert users in order to execute complex or specialised queries. Recent developments to Kbdock make use of a novel protein structure alignment algorithm called "Kpax" (http://kpax.loria.fr) that we have developed [START_REF] Ritchie | Fast protein structure alignment using Gaussian overlap scoring of backbone peptide fragment similarity[END_REF]. This allows queries in Kbdock to span structural neighbours of the retrieved DDIs, thus allowing Kbdock to search over more distant regions of protein structure space and to propose protein docking templates that cannot be found using conventional sequence-based or structure-based comparison techniques.

We are currently working to link KBdock's structural domain binding site classification with the widely used ExPASy Enzyme Classification (http://enzyme.expasy.org) scheme. In order to achieve this, we are developing efficient data-mining approaches to process the millions of sequencefunction associations that are now available in large molecular biology databases such as Swiss-Prot and TrEMBL, which together build the UniProt Knowledgebase (http://www.ebi.ac.uk/uniprot) at the European Bioinformatics Institute.

Useful Links http://kbdock.loria.fr http://kpax.loria.fr http://hex.loria.fr

Contact address:

David Ritchie Tel. +33 (0)3 83 59 30 45 E-mail dave.ritchie@inria.fr


Richard P Brent

Paul Zimmermann
email: paul.zimmermann@inria.fr

Paul Zimmermann 2016 Twelve

TWELVE NEW PRIMITIVE BINARY TRINOMIALS

Keywords: 1991 Mathematics Subject Classification. Primary 11B83, 11Y16; Secondary 11-04, 11N35, 11R09, 11T06, 11Y55, 12-04 GF(2)[x], trinomials, irreducible polynomials, primitive polynomials, Mersenne numbers

published or not. The documents may come

Primitive trinomials of degree up to 32 582 657 were reported in [START_REF]Ten new primitive trinomials[END_REF]. We have completed a search for all new Mersenne prime exponents found by the GIMPS project [START_REF]The Great Internet Mersenne Prime Search, mersenne[END_REF]. Twelve new primitive trinomials were found (see Table 1 Table 1. New primitive trinomials x r +x s +1 of degree a Mersenne exponent r, for s ≤ r/2. For smaller exponents, see references in [START_REF]Ten new primitive trinomials[END_REF] or our web site [START_REF] Brent | Search for primitive trinomials (mod 2)[END_REF].

Our search used the new algorithm [START_REF] Brent | A multi-level blocking distinct degree factorization algorithm[END_REF] relying on fast arithmetic in GF(2)[x], whose details are given in [START_REF] Brent | Faster multiplication in GF(2)[END_REF]. For the squaring of polynomials over GF [START_REF] Brent | Faster multiplication in GF(2)[END_REF][x], we used the new _pdep_u64 Intel intrinsic, which gave a speedup of a factor about 2.5 over the algorithm described in [3, §4]. On a 3.3Ghz Intel Core i5-4590, together with improvements in the gf2x library, we were able to square a degree-74 207 280 polynomial in about 2 milliseconds, and to multiply two such polynomials in about 700 milliseconds. As in [START_REF]Ten new primitive trinomials[END_REF], we produced certificates for non-primitive trinomials, which were checked independently with Magma and NTL (a certificate is simply an encoding of a nontrivial factor of smallest degree). A 3.3Ghz Intel Core i5-4590 takes only 22 minutes to check the certificates of all 37 103 637 reducible trinomials (s ≤ r/2) of degree r = 74 207 281 with our check-ntl program based on NTL [START_REF] Shoup | NTL: A library for doing number theory[END_REF], the largest factor having degree 19 865 299 for s = 9 788 851. Acknowledgements. The authors thank Allan Steel, who independently verified with Magma the twelve new primitive trinomials, and Grégoire Lecerf who verified with Mathemagix the three primitive trinomials of degree 74 207 281. Part of the computations reported in this paper were carried out using a cluster funded by the French ANR CatRel.

  ).

	r	s	Notes
	42 643 801 55981, 3706066, 3896488, 12899278, 20150445 Brent and Zimmermann, 2009
	43 112 609	3569337, 4463337, 17212521, 21078848	Brent and Zimmermann, 2009
	57 885 161	none	Brent and Zimmermann, 2013
	74 207 281	9156813, 9999621, 30684570	Brent and Zimmermann, 2016