Packing topology in crystals of proteins and small molecules: a comparison

We compared the topologies of protein and small molecule crystals, which have many common features – both are molecular crystals with intermolecular interactions much weaker than intramolecular interactions. They also have different features – a considerably large fraction of the volume of protein crystals is occupied by liquid water while no room is available to other molecules in small molecule crystals. We analyzed the overall and local topology and performed multilevel topological analyses (with the software package ToposPro) of carefully selected high quality sets of protein and small molecule crystal structures. Given the suboptimal packing of protein crystals, which is due the special shape and size of proteins, it would be reasonable to expect that the topology of protein crystals is different from the topology of small molecule crystals. Surprisingly, we discovered that these two types of crystalline compounds have strikingly similar topologies. This might suggest that molecular crystal formations share symmetry rules independent of molecular dimension.

The topology of packing in protein crystals may be different from that of other molecular solids, where the packing is, in general, optimal and where there are no internal cavities filled with amorphous material. It is thus interesting to compare the topology of small organic and protein crystals. The most widespread model of molecular crystals is the Kitaigorodskii's model 26 , according to which each molecule in a homomolecular crystal tends to be surrounded by 12 other molecules; other coordination numbers (10 or 14) occur much rarer. However, a comprehensive treatment of all molecular structures stored in the Cambridge Structural Database using a computer-realized Voronoi approach 27,28 showed that coordination number 14 is more abundant than 12. This was explained by the model of deformable molecules, in contrast to the model of rigid ones used by Kitaigorodskii. Since the approach used in references 27,28 can be applied to the molecules of any size and shape it was natural to use it for comparative analysis of packings of protein and small molecules.

Experimental
Data selection. Protein data. All protein crystal structures were taken from the Protein Data Bank 6,7 , according to the following criteria. We discarded structures containing nucleic acids and structures with an average B-factor smaller than 10 Å 2 or larger than 40 Å 2 . We considered only structures with 50-500 amino acids in the asymmetric unit, with resolution better than 3 Å, with working R-factor better than 0.25, and determined at 80-120 K. We discarded structures deposited without their experimental diffraction data. Redundancy was reduced to 30% sequence identity. Structures with missing protein atoms or with protein atoms deposited with zero occupancy were discarded as well as the structures where non-protein and non-water atoms are more than 5% of the total number of atoms.
We note that these are standard criteria for protein crystal structure selection from the Protein Data Bank 29,30 . In particular, resolution and R-factor thresholds ensure that unreliable structures are excluded from the analyzed data sets; moreover, we considered a homogeneous set containing only low temperature crystal structures, which are now the routine, since packing might depend on temperature; furthermore, limitations on the average B-factor and on protein dimension ensure that anomalous structures -much larger or smaller or much more or less flexible than customary proteins -are removed from the data sets; and eventually, redundancy reduction to 30% sequence identity, though very crude, is necessary to avoid distribution biases.
It is also important to observe that the exclusion of structures containing too many (more than 5%) hetero-atoms ensures that crystal contacts due to the presence of co-crystallized small molecules (the so called packing bridges) are minimized. Analogously, the exclusion of incomplete structures, where some of the protein atoms/residues were undetected, ensures that no protein-protein crystal packing contacts are neglected.
We prepared three ensembles of protein crystal structures. In the first, we collected only monomeric proteins that crystallized with only one molecule per asymmetric unit (monomer set; 394 structures). In the second, we grouped dimeric proteins that crystallized with only one dimer per asymmetric unit (dimer set; 207 structures). In the third and last ensemble, we pulled together structures of monomeric proteins with two molecules in the asymmetric unit (double set; 164 structures).
These data sets might seem small if compared to the number of entries deposited in the Protein Data Bank. However, it is necessary to consider that only 31% of the PDB files contain monomeric proteins with one molecule per asymmetric unit; only 7% contain dimeric proteins with one molecule per asymmetric unit; and only 9% contain monomeric proteins with two molecules per asymmetric unit. Moreover, redundancy reduction to 30% maximal sequence identity implies a 6-fold decrease of the number of PDB entries. Eventually, in the large majority of the protein crystal structures deposited in the Protein Data Bank there are residues that have not been seen experimentally 31 . Since they are usually at the protein surface, it is not possible to use these structures to perform reliable analyses of the crystal packing. We note also that the quality of the data extracted from the Protein Data Bank is becoming a crucial issue in structural bioinformatics 29,32 . It is also important to observe that although these are high quality and non-redundant datasets of protein crystal structures, it is impossible to estimate the extent to which they represent the entire protein "universe", a considerable fraction of which is constituted, for example, by membrane proteins that are under-represented in the Protein Data Bank 33 .
The identification codes of all these crystal structures are listed in the Supplementary Information (Table S1).
Small molecules data. Crystallographic data for organic molecular crystals consisting of chemically equivalent molecules were taken from the CSD (release 5.38). We have not analyzed molecular crystals of coordination or organometallic compounds, which have different chemical nature compared to proteins. In total, we have considered all 166,297 structures, which contain one independent molecule and 10,240 structures with two independent chemically equivalent molecules in the asymmetric unit. The entries containing incomplete, erroneous data, or disordered structures as well as those with R f > 10% were excluded; no other restrictions were applied.
Crystal contacts identification. Protein data. The crystal contacts were identified as described previously 13 .
The asymmetric unit was transformed according to all symmetry operations and translated up to three times along the axes a, b, and c in both negative and positive directions. Two molecules -one, the reference molecule, is in the asymmetric unit and the other, the satellite, is symmetry related to the first -were considered to be in contact if at least one atom of one of them is closer than 4.5 Å from an atom of the other one. We refer to this method below as 'Distance' method.
Small molecules data. All interatomic contacts were determined by means of the 'Domains' method 34 implemented into the ToposPro package 35 . This method is based on the Voronoi partition and accounts not only for interatomic distances, but also for atomic surrounding and screening effects. An interatomic contact A-B is assumed to exist if there is a face of an atomic Voronoi polyhedron with a solid angle, Ω > 1.5% of the total solid angle 4π steradian and the line A-B crosses the face (the contact is direct). To distinguish valence and non-valence contacts a set of additional criteria was used as described in reference 34 . The 'Domains' method is designed to process large samples of crystal structures in an automated mode, which was crucial for our study. At the same time, this method was never used for proteins; thus, it was important to check if this method gives similar results for proteins as the 'Distance' method. The molecular coordination number was computed for the molecular Voronoi polyhedron, which is a union of Voronoi polyhedra of all atoms of the molecule ( Fig. 1) 27,28 .
External faces of the molecular Voronoi polyhedron correspond to intermolecular contacts, and two molecules are considered interacting if there is at least one intermolecular contact between them ( Fig. 2) 27 .
This method was also used to compute coordination numbers of some protein molecules for comparison with the simple geometrical approach described above.

Topological analysis.
To analyze molecular packings of both proteins and small molecules we squeezed each molecule into its center of gravity to obtain the so-called underlying net, i.e. a net, which keeps the information on the connection of molecules throughout the crystal, and ignores the information about their geometrical properties and internal structure (Fig. 3) 36 . As a result, the crystal structure is represented as an infinite graph, whose vertices and edges correspond to molecules and intermolecular links, respectively. Such simplification of the structure allows us to put aside differences in size, shape and other geometrical features of molecular ensembles and to focus on the method of their local and overall connection, i.e. on the crystal architecture as a whole. This is a way to find correlations between packings of so different objects as small molecules and proteins.
The packings represented by the same underlying net are topologically equal and belong to the same isoreticular series, which is designated by a topological type symbol. There are several nomenclatures for topological types; here we use two of them: the RCSR three-letter codes 37 and the NDn symbols used in the collection of ToposPro topological databases (TTD) 38 , where N is the series of coordination numbers for all independent nodes, D is one of the letters C, L or T, indicating the dimensionality of the net (C, L, T for chain, layer, and three-periodic framework, respectively) and n is the number of the net with a given ND sequence in the TTD collection. For example, the 16T4 symbol means a three-periodic net with 16-coordinated nodes (i.e. any molecule has 16 neighbors in the packing) and it is the fourth 16-coordinated net in the TTD collection.
Since the strength of intermolecular contacts is quite different, the topology of a molecular packing strongly depends on which contacts are taken into account. Ignoring interactions of a particular level of strength, we generate different underlying nets for the same crystal structure. Each underlying net specifies the packing topology at a given level of interaction and corresponds to a structure representation.
For small molecule crystals, we estimated the strength of an intermolecular interaction by the value of Ω mol that is the part of the sum of solid angles of all contacts, which are formed by a given molecule; this part corresponds to all contacts between a given pair of molecules. ToposPro includes a special procedure, which generates all possible representations by subsequently breaking groups of weak intermolecular contacts with close Ω mol (the difference between neighboring Ω mol in the group does not exceed 1.5% of the total solid angle 4π steradian). As a result the structure is considered at different levels of intermolecular interactions. The underlying net was constructed for each level of intermolecular interaction to obtain the information about the way of molecular connection. Such a multilevel analysis, which we successfully used for elucidation of organic molecular crystals 39 , allows one to separate underlying motifs, which determine the structure architecture. These motifs can be the same in the molecular packings, which differ in details of weaker intermolecular contacts (Fig. 4). In our study we considered only the representations with molecular coordination numbers ≤14; there are 105,549 molecular packings with such connectivity among the 166,297 structures studied. The topologies of molecular packings with a higher connectivity are non-specific and can be usually reduced to a known topology with a smaller coordination number after discarding weak intermolecular contacts.
For protein crystals, strengths of crystal packing contacts were estimated based on their dissociation free energies as described below.
Crystal contacts thermodynamics. Since protein crystallization occurs in aqueous solutions and crystals are largely impregnated by aqueous solutions, the formation and stability of crystal packing contacts cannot be considered out of solvent context. For this reason, we used the program PISA 40 to estimate the stability of individual crystal packing contacts formed by the reference molecule and each of its satellites, which were generated by applying all the symmetry operations to the asymmetric unit and by applying three translations along all the axes in the negative and in the positive directions. Intermolecular contacts were defined as pairs of atoms, one in the reference molecule and the other in a satellite molecule, closer than 4.5 Å.
According to the PISA definitions, the dissociation free energy ΔG diss is defined as: where ΔG int is the binding free energy of the two molecules (the reference molecule plus one of its satellites), T is the absolute temperature, and ΔS is the entropy change upon association of the two molecules. If ΔG diss > 0, the crystal packing contact is thermodynamically unstable in solution and is present only in the crystal as an experimental artifact. ΔG int is defined as:

Results and Discussion
Overall topologies of molecular packings. Since packing contacts were determined with two different approaches in protein crystals ('Distance method' , see Methods) and in small molecule crystals ('Domains method'), it is mandatory to compare the results obtained by following these two different approaches and to verify that they produce identical results. For this purpose we have computed intermolecular contacts for some proteins with the 'Domains method' . In most cases, the results of both methods coincide (Fig. 5a), however, even if they differ, this concerns only weak long contacts, which the 'Domains' method can consider in contrast to the 'Distance' method (Fig. 5b). Since the weak contacts are not considered in the multilevel topological analysis, the topologies obtained by the two methods can be assumed to be comparable. Table 1 shows the 20 most frequent underlying nets for monomeric protein crystals as well as for small molecule crystals. One can see that small molecules tend to dense packings with coordination numbers ranging from 12 to 14, in agreement with previous results 27, 28 . On the contrary, coordination numbers are very often much lower in protein crystals, in agreement with the fact that a considerable fraction of their volume is occupied by liquid water. However, all underlying nets of protein crystals in Table 1 are well known and almost all of them (except nets chb, smt, yfh) are fcu and hcp subnets. Interestingly, the 14-coordinated body-centered cubic packing (bcu-x), which is the most frequent in small molecule crystals is rather infrequent in protein crystals  (Table 1). Thus protein crystals follow the close packing topologies but with some gaps, which are obviously filled by the solvent molecules and decrease coordination numbers of protein molecules.
It is interesting to observe that while few topologies are extremely more frequent than others in small molecule crystals, this is not observed in protein crystals. 30% of the small molecule crystals are associated with the most frequent topology (bcu-x) and 65% of them have one of the four most common topologies (bcu-x, tcg-x, 14T3, gpu-x). On the contrary, only 6% of the protein crystals are associated with the most frequent topology (fcu) and only 19% of them have one of the four most common topologies (fcu, bct, hex, vcs). This suggests that the suboptimal packing of the protein molecules in the crystal state allows a wider number of topologies and none of them can be much more frequent than the others.
Similar results are observed for the crystals of dimeric proteins where there is one dimer per asymmetric unit (Table 1) and for the small molecule crystals where there are two independent molecules in the asymmetric unit ( Table 2). On the contrary, the underlying nets of crystals of monomeric proteins where the asymmetric unit contains two independent molecules set are very diverse and it was impossible to separate any preferable motif ( Table 2).    Table 2. Most frequent underlying net topologies in crystals of protein and small molecules with two molecules in the asymmetric unit.  Table S2) compared to the CN for small molecules. numbers for protein molecules are 7 and 8, while small molecules prefer higher coordination numbers equal to or larger than 14 ( Fig. 6; see also supplementary Table S2). Like for the overall topology analysis, we see that the most frequent coordination number of small molecules (14) is observed in a large fraction of crystals (52%) while the most frequent coordination number for proteins (7 or 8) is observed only in a smaller fraction of crystals (about 20%). Moreover, small molecules can adopt 25 different coordination number values (from 4 to 32) while proteins can have only 14 different coordination number values (ranging from 1 to 16 -coordination numbers 2 and 15 are never observed).
All these observations point out that proteins, which are suboptimally packed in their crystals, can hardly by surrounded by numerous other proteins and that the palette of their coordination numbers is considerably more limited.
Rare overall topologies. Although the overall topologies are quite similar in protein and in small molecule crystals, protein molecules can form packings of a topology that never occurs in packings of small molecules. Among 394 monomeric and 207 dimeric protein structures there are 31 and 10 underlying nets, respectively, which were not found in the packings of small molecules, and all of them have a unique topology (see Table S3 in the Supplementary Material). Among the crystal structures with two independent protein molecules per asymmetric unit, there are 44, which have unique overall topologies never revealed in the small molecule crystals ( Table 2). This means that from 5 to 26% of protein molecules essentially differ from small molecules by their ability to be packed. Multilevel topological analysis. In order to reach a better understanding of the differences between the topologies of protein and small molecule crystals, we performed the multilevel analysis for the molecular packings in 394 monomeric protein structures and 105,549 structures of small molecules. As a result, 1,799 and 1,144,539 structure representations were generated, respectively.
The underlying net topologies are summarized in Table S4 (see Supplementary Material). In general, the underlying motifs that are frequent in protein crystals are also frequent in small molecules crystals. The frequency of the protein topologies and the frequency of the small molecule topologies shown in Table S4 correlate well (Pearson correlation coefficient = 0.625; Spearman correlation coefficient = 0.563) as shown in the histogram in Fig. 7.
It appears that the underlying motifs of packings in both protein and small molecule crystals have simple and quite common topologies, such as sql, hcb, dia, pcu, hex (Fig. 8). The role of close-packed motifs (fcu, hcp, Fig. 8) is expectedly higher for small molecules than for macromolecules. Since the multilevel analysis reveals the basic motifs (skeletons), which are formed by the strongest intermolecular contacts, the results obtained mean that that the skeletons in molecular crystals are mostly independent of the size and shape of the molecule, while the total packing can vary in proteins compared to small molecules.
Topology versus thermodynamics. To find correlations between thermodynamic parameters and molecular connectivity we examined all the 394 protein structures of the monomer set. A total of 3,156 crystal packing contacts were identified. We observed that the number of intermolecular contacts (NC) is strictly related to the ΔG diss . PISA does not compute the number of intermolecular contacts but computes the buried solvent area (BSA; the patch of solvent accessible surface of molecules A and B that become inaccessible to the solvent because of the interaction between A and B). There is a strict relationship between NC and BSA: with Pearson's correlation coefficient of 0.945. This is obviously not surprising and it means that a single intermolecular contact buries nearly 7 Å 2 of solvent accessible surface. This means that about 30% of the surface of an atom of radius equal to 1.4 Å is buriable in a protein-protein interface, while the rest is buried intra-molecularly. Also the relationship between BSA and ΔG diss is strict:

Conclusions
The topologies of carefully selected high-quality protein crystal structures have been determined and compared to those of a large set of small molecule crystal structures. It is known that protein crystals are suboptimally packed, since they contain channels filled with liquid water. Therefore, the coordination numbers in protein structures are considerably lower than in small molecule crystals. Surprisingly, however, we have found their topologies very similar. Although it is impossible, based on the data presented in this communication, to provide a strict and definitive explanation for this surprising similarity, it is clear that some molecular mechanism must exist. It is possible to hypothesize, for example, that during the early stages of nucleation, only relatively few types of assemblies, with a specific stereochemistry, are sufficiently stable to be able to accumulate and to continue to grow into real crystals. At this regard, it might be interesting to remember that three main differences exist between protein and small molecule crystallization [41][42][43][44] . First, the kinetics of protein crystal nucleation and growth are, in general, two or three orders of magnitude slower than for small molecules, because of the protein larger size and lowered diffusivity 45,46 . Second, protein crystals nucleation occurs at very high levels of supersaturation, often two or three orders of magnitude greater than that required to sustain crystal growth. Third, proteins may assume several distinctive solid states that include amorphous precipitates, oils or gels as well as crystals, and most of these are kinetically favored. Further studies are necessary and additional data must be considered to find out the rationale why topology seems to be independent of packing efficiency and crystallization.