Motivated by the biophysical constraints of the polypeptide chain (Goncearenco and Berezovsky, 2010, 2011; Berezovsky et al., 2017a; Berezovsky, 2019), the procedure starts from 30-residue long seed sequence fragments of the functional loops represented by a gapless multiple sequence alignment and a position-specific scoring matrix (PSSM) profile. The sequence profile is then iteratively scanned against the non-redundant UniRef database (Hunter et al., 2009) with an expectation–maximization (EM)-like algorithm (Goncearenco and Berezovsky, 2010, 2011) converging to an expanded sequence profile of the descriptor with a functional signature (Berezovsky et al., 2003a,b). The corresponding structures characterizing the functional loop are then obtained by looking for profile matches against the sequences in the Protein Data Bank (Berman et al., 2000; wwPDB consortium, 2019), followed by extraction and encoding of the structural feature in the form of parametrized probability distribution functions. Structural and functional annotations are extracted from the corresponding databases, such as Uniprot, MaCiE, and Conserved Domains Database (CDD) (Bairoch, 2000; Andreini et al., 2009; Fischer et al., 2010; Holliday et al., 2012; Marchler-Bauer et al., 2013; Akiva et al., 2014; Furnham et al., 2014). The structural features include dihedral angles, Van der Waals (VdW) interactions, and hydrogen bonds (H-bonds). Intra-EFL interactions and interactions between the functional loop and the rest of the structure are encoded separately. Thus, a descriptor contains information about the immediate environment of the functional loop in all protein scaffolds and enzymatic functions where it was encountered.

In protein engineering and de novo design, the descriptors of elementary function need to be integrated into a given structural scaffold. Once the elementary function that is required to be incorporated into the protein is selected, the sequence, structure, and interactions that would fit best into the scaffold have to be determined. The objective function scores how well a given structural loop fits in and should be maximized to obtain the best matching implementation of the descriptor with an assumption that the native structure has the best fit. Essentially, the score represents the joint likelihood of all amino acid positions in the grafted loop with respect to distributions parametrized in the descriptor with N-residue positions and M features:  F=∑iN∑jMWiPWijFSij. The weight that is given to a residue position WiP reflects the relative degree of conservation of features in each position: WiP=∑j=0MWijF∑i=0N∑j=0MWijF . Descriptor features j enumerate along a sequence signature (j = a), dihedral angles (j = d), H-bonds (j = h), and vdW interactions (j = v), with the corresponding scores S_ij and weights w_ij. Score S_i,a is the log-odds score for amino acid substitution according to BLOSUM62 (Henikoff and Henikoff, 1993). The weight of the residue feature is given by the relative frequency of the two most frequent residues in the sequence profile Wi,a F=∑k2(argmaxkwk)∑20wk, where w_k is amino acid frequency. Dihedral angles are first clustered as two-dimensional vector quantities using the EM algorithm as implemented in scikit-learn, and their weights and fitting scores are derived from parameters in the trained model with the following equation:

Si,d=1-erf(∑DTCD), where erf refers to the error function, D = X - Y is the displacement of the compared dihedral angle points, normalized by the median of all points in the descriptor, expressed as phi-psi angles, and C is a matrix proportional to the degree of confidence that the point belongs to the distribution. Given the precision matrix Λ=σ⁻² and posterior probability matrix P obtained from the clustered model, C = Λ P.

The isotropic VdW distribution score, S_i,v, is based on counting the number of VdW contacts with non-hydrogen atoms within the 5-Å radius. Hydrogen bonds are additionally split into acceptors and donors: S_i,hA and S_i,hD, respectively. They are scalar quantities and are measured as the number of donor H-bonds present at the residue position. For H-bonds, the effective radius is 3.5 Å. The score is the ratio between the absolute difference in the number of bonds (δ) to the higher number of contacts for either structure (VdW interactions and H-bonds): Si,  [d,v,hA,hD] =abs(δx-δy)max(δx,  δy). The weight of the feature is the SD of the spread of a half-normal distribution fitted onto the data Wi [b,c,d,e]F=σ. Each feature weight is further tuned by an empirically derived scalar factor to account for differences in the spread of absolute values returned by the scoring functions.

Figure 1 contains examples of some of the structural features for three descriptors: the calcium-binding Ca²⁺-binding helix–loop–helix EF-hand motif (Gifford et al., 2007) (with the characteristic signature DxDxD, Figure 1A), the glycine-rich motif of the phosphate-binding in dinucleotide-containing ligands (with the characteristic signature GxGxxG, Figure 1B), and the phosphate-binding P-loop in nucleotide-containing ligands (with the characteristic signature GxxGxG, Figure 1C). For EF-hand (Gifford et al., 2007), dihedral angles on amino acid positions before and after the calcium-binding sites form tight clusters, i.e., the backbone is structurally conserved. Donor-acceptor hydrogen bond pairs spaced three to four residue positions apart are also present in the first and last 10 residues of the EF-section (Figure 1A), representing, together with conserved dihedral angles, two α-helices in the structural motif. Distributions of VdW contacts in the EF-hand motif show a large proportion of the conserved intrinsic contacts flanking the α-helices, whereas the central flexible link between them forms vital external contacts with the rest of the fold. The same patterns of the internal contacts are observed in the second halves of the GxGxxG and GxxGxG motif structures also containing the α-helices. The β-strands in these loops, as expected, interact more strongly with the rest of the structure via VdW interactions (Supplementary Figure 2) and H-bonds (Supplementary Figures 3, 4). In addition to the interactions with ligands, the substrate-binding residues in the descriptor are characterized by multiple interactions with the fold.

The objective function was calibrated and exemplified on three motifs: (i) the EF-hand from Bos taurus calcium-binding protein structure (PDBID 1A29) with the Ca²⁺ ligand; (ii) the phosphate-binding motif (GxGxxG) from Equus caballus oxidoreductase (PDBID 1A71) with dinucleotide-containing NAD ligand; and (iii) the phosphate-binding motif (GxxGxG) from Methanocaldococcus jannaschii ABC transporter (PDBID 1G6H) with nucleotide-containing an ADP ligand. These motifs were replaced with the realization of the corresponding descriptors (see Supplementary Figure 5 and explanations in the figure caption).

Supplementary Figure 6 illustrates how objective function can be used to obtain the best descriptor realization in the corresponding engineering or design tasks. First, a segment of seven consecutive residues (half of the typical functional signature; Berezovsky et al., 2017a) with the highest cumulative fitting score is determined. Then, this segment is extended with residues that contribute the highest scores to the fitting function. There might be insertions/deletions in the descriptor realization, or the reference could be undefined as in the case of missing data, with disordered or structurally unresolved regions in the query structure. In such a case, the algorithm first scores all valid residue positions and finds the appropriate segments to merge if the best match does not come from a singular structure. Next, adjacent segments are extended, possibly from both ends in the case of a gap, to fill in the uncharacterized positions. This returns the best-effort re-engineered structure that partially relies on the input structure and fills in the rest, where data are missing. Supplementary Figure 6 exemplifies the case of descriptor realization that builds the functional loops out of two segments, providing the most optimal score.

To assess the robustness of the derived descriptors and versatility of the fitting function, we selected seven structures with the phosphate-binding functionality, but of different folds, origins, and ligands: NADP-binding Thermotoga maritima lactate dehydrogenase (1A5Z) and Methylobacterium extorquens AM1 methylene-tetrahydromethanopterin dehydrogenase (1LUA), AMP-binding Escherichia coli MoeB-MoaD protein complex (1JWB), flavin–adenine dinucleotide (FAD)-binding Klebsiella pneumoniae udp-galactopyranose mutase (2BI7) and Thermus thermophilus GidA-related protein, and both binding sites for FAD and NADP in E. coli CoA reductase (1PS9) and human dihydrolipoamide dehydrogenase (1ZMC).

The donor functional loop that had to be transplanted is an elementary function of phosphate binding in dinucleotide-containing ligands (GxGxxG). We derived its descriptor from the dataset that excluded the aforementioned structures. Seven proteins representing variations in the structure and sequences of the phosphate-binding functional loops (both in nucleotide- and dinucleotide-containing ligands) were used as the targets or acceptors for the EFL-replacement procedure, using the above descriptor. Figure 2 illustrates the fit between the re-engineered functional loop into the original structure, which binds various ligands, including AMP, FAD, and NAD(P). Additionally, recombinant functional loops can be built from the best matching segments of multiple structures. Taking 1PS9 structure as an example, its NADP binding site is replaced with the loops from structures 1KF6 from E. coli quinol-fumarate reductase; 1VRP from C. sarcosine oxidase; 6GNC from Clostridium acetobutylicum thioredoxin reductase; and 5ER0 from Lactobacillus oxidase. All of these structures contain FAD-binding sites replacing the original NADP-binding site.

To further illustrate the utility of descriptors of elementary functions and the potential of the DEFINED-PROTEINS package in protein design applications, we set up initial conditions for the grafting procedure where we replaced the original functional loop in a given protein with the non-native one representing a different elementary function. We recruited two distinct, but not opposite, elementary functions, for which we have already derived descriptors: binding of the phosphate in dinucleotide-(GxGxxG) and nucleotide-containing P-loop (GxxGxG) ligands (Zheng et al., 2016). Thereby, we cross-grafted nucleotide-binding function into protein folds with native dinucleotide-ligand binding (Figures 3A–C and Supplementary Figures 7A–C), and cross-grafted dinucleotide-binding function into folds with the native mononucleotide binding ability (Figures 3D–F and Supplementary Figures 7D–F), respectively. Although both are elementary functions of the phosphate binding, the latter belong to different dinucleotide- and nucleotide-containing ligands that determine the corresponding diversity of protein functions that these elementary functions evolved into.

The signatures of EFLs and the most frequent interactions with distinct parts of ligands were discussed elsewhere (Berezovsky, 2019), showing both conservatism and diversity depending on the position in the functional loop. From the structural perspective, both EFLs have similar secondary structures of β-turn-α-helix composition and architecture. The latter results in a high number of intrinsic H-bonds and VdW contacts in its second α-helical part, while β-strand elements of both loops form more contacts with the surrounding structure. At the same time, the dinucleotide-binding-GxGxxG-loop is a more compact structure by itself, with more intrinsic contacts between its α and β elements. We replaced the original functional segments with recombinant (Figure 3) and deconvoluted single-loop (Supplementary Figure 7) functional loops sampled from the corresponding descriptors. In all of the cases presented in Figure 3 and Supplementary Figure 7 replacements were done with the highest-scoring matches.

Figures 3A–C and Supplementary Figures 7A–C show results of grafting of the P-loop (GxxGxG) descriptor in places of Gly-rich signature (GxGxxG) of the phosphate-binding in dinucleotide-containing ligands. While obtaining a good match between the original and replacement loops in both unblended single-loop and recombinant cases, the latter was allowed to obtain a better per-residue fit between the original and replacement loops. Recombinant loop replacements consist of several fragments (typically, two to three), covering most of the loop (see examples in Figure 4 and Supplementary Figure 8) with higher scores in the functional positions at the expense of non-function-bearing positions in the secondary structure elements. It agrees with less conserved dihedral angles in the turn segments of the loops that contain functional signatures, which also result in different angles between the secondary structure elements of the loop and difference in the intra-loop contacts in GxxGxG and GxGxxG loops.

Figure 4 and Supplementary Figure 8 show a comparison of secondary structure alignments in several examples of the single-loop and recombinant replacements. The secondary structure is affected by the environment (Minor and Kim, 1996), pointing to the need for adjustments and optimization in length and location after the original descriptor realization is placed instead of the natural ELF. Overall, examples of structural replacements (Figures 3 and Supplementary Figure 7) and alignments (Figure 4 and Supplementary Figure 8) show that realizations of descriptors can be used as a starting point for further optimization of positions and interactions involving functional residues and structures of elementary functional loops to obtain required modifications/design in the context of new protein structures and functions.

It is important to note that, in general, scores obtained with the objective function depend on the sequence/structure characteristics of the elementary function and the type of procedure in which the descriptor is used. The major contributors to the objective function score are sequence conservation and dihedral angles, with weights 0.62/0.52 and 0.3/0.38 for the GxxGxG/GxGxxG signatures, respectively. Weights for VdW interactions and H-bonds are smaller (see Supplementary Table 1 legend for details): the VdW interactions are weak though omnipresent, whereas H-bonds between the loops and rest of the folds are rare, making the weights of both smaller. Nevertheless, contributions of all characteristics to the weight for each individual position are calculated as a sum of their weights normalized by the total feature weights across all residue positions. The score of the reference state for the descriptor, which can be indicative of the latter and can be used as a guideline in design, should be obtained for each descriptor. The straightforward way is to perform the cross-validation experiment, which assesses the robustness of the descriptor and provides the score that can be used as a ground state score. Supplementary Table 1B contains averaged scores obtained in the cross-validation experiment (Figure 2), which can be used as a reference for the engineering and design of corresponding descriptors in other folds. As an example, Supplementary Table 2B shows scores in case of cross-grafting of descriptors of the phosphate-binding signatures in the nucleotide- (GxxGxG) and dinucleotide-containing (GxGxxG) ligands, revealing the deviation from scores in cross-validation that can further increase in case of de novo design of functions based on the descriptors. The utility of averaged (Supplementary Tables 1A, 2A) and per-residue scores is in the information on the relative conservatism (importance) of positions in the descriptor and of the overall match between its realization and the rest of the fold and its capacity to contribute to engineered/designed function. For example, depending on the requirements on the interactions within the ELF and between the loop and the rest of the fold, H-bonds can be introduced in positions with a conservatism level allowing to do so, but the same holds for changing other characteristics.