Identification of L-asparaginases from Streptomyces strains with competitive activity and immunogenic profiles: a bioinformatic approach

Identification and selection of homologous L-Asparaginases

Putative ASNases from Streptomyces were identified through a BLASTp search against the NR database of the NCBI using as seeds the amino acid sequences of EcAII (ID P00805) and Streptomyces coelicolor type II ASNase (ScAII; ID Q9K4F5). The search was restricted to the Streptomyces taxon (Taxonomic ID number: 1883), and an E-value less than 1e−06 was considered significant. Partial proteins and those from unidentified Streptomyces strains were excluded. In a posterior step, the set of protein sequences was filtered at 60% identity as cutoff to avoid redundancy, using the CD-Hit program (http://weizhongli-lab.org/cdhit_suite/cgi-bin/index.cgi) (Huang et al., 2010). Each cluster was analyzed using the HMMER program on the PFAM server (http://pfam.xfam.org/) to determine the protein family to which they belonged (Finn, Clements & Eddy, 2011; Finn et al., 2016).

Phylogenetic analysis

ASNases amino acid sequence alignments were performed using Clustal Omega (Sievers et al., 2011) with default parameters. The quality of the alignments was improved using the model PF06089.11 or PF00710.11 of ASNase, as required. Multiple sequence alignment statistics were computed with AliStat (http://www.csb.yale.edu/userguides/seq/hmmer/docs/node27.html).

Phylogenetic analyses were carried out using the maximum-likelihood method with the program Mega 7. The WAG model was chosen as substitution model, and 1,000 replicates were performed. The best tree was calculated using the majority rule. Additionally, E. coli type I ASNase (EcAI) was included in the phylogenetic analysis of the PF00710.11 cluster. EcAI is closely related to EcAII but it does not have therapeutic potential. For the PF06089.11 cluster, Rhizobium etli type II ASNase (ReAII) was included in the analysis.

Antigenicity prediction

The prediction of the probability of antigenicity of each ASNase was calculated with the server ANTIGENpro (http://scratch.proteomics.ics.uci.edu/) (Magnan et al., 2010). ANTIGENpro is a sequence-based, alignment-free, protein antigenicity predictor with an estimated accuracy of 82%.

HLA class II binding prediction

The amino acid sequence of each candidate ASNase was screened for T-cells epitopes with the MHC II Analysis Resource at the Immune Epitope Data Base (IEDB) server (http://tools.iedb.org/mhcii/). MHC II Analysis Resource parses sequences into 15-mer and assesses the binding potential of each 15-mer to MHC class II molecules of one or more HLA alleles. The IEDB recommended method was used for predictions for a set of 8 HLA alleles that collectively represent >95% world population: HLA-DRB1*01:01, HLA-DRB1*03:01, HLA-DRB1*04:01, HLA-DRB1*07:01, HLA-DRB1*08:01, HLA-DRB1*11:01, HLA-DRB1*13:01 and HLA-DRB1*15:01. The IEDB-recommended method uses the consensus approach, combining NN-align, SMM-align, CombLib, Sturniolo, and NetMHCIIpan (Wang et al., 2010). For each peptide, a percentile rank is generated by comparing the peptide’s score against the scores of five million random 15-mer selected from SWISSPROT database, and the median percentile rank is used to calculate a consensus percentile rank (CPR). Peptides with a CPR < 2 were defined as high-affinity binders and thus selected for epitope density (ED) calculation. Multiple 9-mer cores were identified in overlapped 15-mer peptides. To reduce overestimation of predicted peptides, only the 9-mer cores, predicted by using the Sturniolo method (Sturniolo et al., 1999) and with a CPR < 1, were considered for the analysis. Finally, epitope density (ED) was calculated using the follow equation, modified from (Santos et al., 2013): $$\mathrm{E}\mathrm{D}=\frac{\mathrm{P}\mathrm{r}\mathrm{e}\mathrm{d}\mathrm{i}\mathrm{c}\mathrm{t}\mathrm{e}\mathrm{d}\mathrm{e}\mathrm{p}\mathrm{i}\mathrm{t}\mathrm{o}\mathrm{p}\mathrm{e}\ast (2-\mathrm{A}\mathrm{f}\mathrm{f}\mathrm{i}\mathrm{n}\mathrm{i}\mathrm{t}\mathrm{y}\mathrm{a}\mathrm{v}\mathrm{e}\mathrm{r}\mathrm{a}\mathrm{g}\mathrm{e}(\mathrm{c}\mathrm{p}\mathrm{r}))}{\mathrm{P}\mathrm{r}\mathrm{o}\mathrm{t}\mathrm{e}\mathrm{i}\mathrm{n}\mathrm{l}\mathrm{e}\mathrm{n}\mathrm{g}\mathrm{t}\mathrm{h}\mathrm{s}\mathrm{i}\mathrm{z}\mathrm{e}-\mathrm{E}\mathrm{p}\mathrm{i}\mathrm{t}\mathrm{o}\mathrm{p}\mathrm{e}\mathrm{s}\mathrm{i}\mathrm{z}\mathrm{e}+1}$$

Where Predicted epitope is the number of epitopes with a CPR < 1.

Epitope coverage was calculated as the number of alleles covered by the epitope consensus, according to the following assumption: when a small number of alleles is covered, a lower percentage of the population will develop sensitivity to ASNase.

Protein structure prediction, refinement and quality assessment

The three-dimensional structures of the selected ASNases was modeled by homology using the I-Tasser server (https://zhanglab.ccmb.med.umich.edu/I-TASSER/) (Zhang, 2008). In brief, starting from an amino acid sequence, I-Tasser generates three-dimensional atomic models from multiple threading alignments and iterative structural assembly simulations. A C-score, provided as an estimate of the accuracy of the models generated, typically ranges between −5 and +2, with a higher value indicating higher confidence, and vice versa (Roy, Kucukural & Zhang, 2010).

For each ASNase, the model with the higher C-score was selected and then refined using the ModRefiner server (https://zhanglab.ccmb.med.umich.edu/ModRefiner/). ModRefiner improves the physical quality and structural accuracy of three-dimensional protein structures by a two-step, atomic-level energy minimization (Xu & Zhang, 2011). Finally, the quality of the models was evaluated by PROCHECK (https://servicesn.mbi.ucla.edu/PROCHECK), Qmean (https://swissmodel.expasy.org/qmean/), and Verify3D (http:/servicesn.mbi.ucla.edu/Verify3D).

Molecular docking

The molecular coupling was carried out using Autodock Tools software (Sanner, 1999; Morris et al., 2009). EcAII (PDB ID: 3ECA) was recovered from the PDB protein database (http://www.rcsb.org/) (Swain et al., 1993; Berman et al., 2000). Once refined, selected ASNase structures were prepared using Dock prep at UCSF Chimera and refined using the Gasteiger method (Gasteiger & Marsili, 1978).

The three-dimensional structures of the asparagine and glutamine ligands were obtained from the DrugBank repository (https://www.drugbank.ca/; accession numbers DB00174 and DB00130, respectively) (Wishart et al., 2018). The preparation of the ligands was carried out by the Gasteiger method and, finally, the allocation of the rotation centers was determined (Gasteiger & Marsili, 1978).

For each ASNase, the search box was focused on previously proposed active sites. The box size was defined to cover all residues of the ligand binding site, using a grid size of 0.375 Å.

Blind molecular docking was performed with Autodock 4.2 software, using the Lamarckian genetic algorithm, with 1,000 runs, for a population size equal to 150, with 2.5 × 10⁶ evaluations, a mutation rate equal to 0.02 in 27,000 generations.

In addition, the active site location was predicted by AutoLigand (Harris, Olson & Goodsell, 2008). Briefly, AutoLigand identifies sites of maximum affinity from maps generated by AutoGrid, finding regions with better energy and a lower volume.

L-Asparaginases from Streptomyces cluster into two type families according to its protein architecture

The Blast search against the Streptomyces taxon revealed 296 putative ASNases homologous to EcAII and 703 homologous to ScAII with a significant score. After manual examination of both groups, 136 and 311 complete sequences were kept for EcAII and ScAII groups, respectively. Protein domain analysis using PFAM server showed that 136 sequences are related to the PF00710.11 family of N-terminal ASNases. For sequences homologous to ScAII, PFAM analysis revealed that they belong to the PF06089.11 family of ASNases, a group of enzymes related to ReAII, a thermolabile enzyme induced by L-asparagine and repressed by the carbon source (Moreno-Enriquez et al., 2012; Huerta-Saquero et al., 2013). Representative clusters for PF00710.11 and PF06089.11 families obtained using the CD-Hit suite program were generated at a 60% identity cutoff, with 19 and 7 putative ASNases, respectively (Table 1). ASNases sequences showed similar lengths in both clusters, ranging from 320 to 420 amino acids.

The sequences belonging to the PF00710.11 family have conserved residues located at the ligand binding site necessary for L-asparagine hydrolysis: Thr 12, Tyr 25, Ser 58, Gln 59, Thr 89, Asp 90, and Lys 162 for subunit A; Asn 248 and Glu 283 for subunit C. In this regard, Thr 12–Lys 162–Asp 90 and Thr 12–Tyr 2–Glu 283 are the catalytic triads involved in L-asparagine hydrolysis, where Thr 12 and Thr 89 are involved in the nucleophilic attack of the substrate (Gesto et al., 2013; Sanches, Kraunchenko & Polikarpov, 2016).

Concerning the PF06089.11 family, we identified an N-terminal conserved motif, with sequences NCSGKHxAM, DGCGAPL, SHSGEx(2)H, and PRSx(2)KPxQ probably involved in asparagine hydrolysis. ReAII hydrolyzes L-asparagine at similar levels to Erwinia chrysanthemi, but with lower affinity than L-asparaginases from both E. coli and E. chrysanthemi (Moreno-Enriquez et al., 2012). Furthermore, ReAII is the only ASNase characterized from the PF06089.11 family.

Phylogenetic analysis of ASNases

For the PF00710.11 family, EcAI was added to the multiple sequence alignment in order to know the relationship between this ASNase and the candidate ASNases. EcAI belongs to the same family of proteins as EcAII, but it does not represent a therapeutic option for ALL treatment. It is noteworthy that asparaginases can also be classified according to subcellular localization, (a) periplasmic asparaginases containing secretion signal peptide and, (b) asparaginases with intracellular localization. The former generally have a higher affinity for asparagine. However, according to their architecture, both types of proteins can be found in the PF00710.11 or PF06089.11 families. This is the case of E. coli asparaginases I and II, both belonging to the PF00710.11 family (https://pfam.xfam.org/family/PF00710#tabview=tab1). We found that the ASNase with accession number WP_059134811.1 of Streptomyces alboniger is grouped in the same clade as EcAI, and so it was excluded from subsequent analyses (Fig. 1A).

The phylogenetic reconstruction showed three well-defined clades (Fig. 1A). The first clade includes ASNases from Streptomyces species S. aureocirculatus (WP_078965752.1), S. cattleya (WP_014151616.1), S. thermoautotrophicus (KWW98572.1), S. himastatinicu (EFL23513.1), S. turgidiscabies (ELP65653.1), S. nanshensis (WP_070201703.1), and S. griseus (WP_030748190.1).

The second clade includes ASNases from S. albidoflavus (WP_095730579.1), S. kebangsaanensis (WP_073950513.1), S. fradiae (WP_078649241.1), S. himastatinicus (WP_009718687.1), S. purpureus (WP_078513220.1), and S. paucisporeus (WP_079189481.1). Finally, the third clade contains proteins from S.purpurogeneiscleroticus (WP_053609500.1), S. purpurogeneiscleroticus (WP_053610569.1), S. phaeochromogenes (WP_055617501.1), and S. lavenduligriseus (WP_051815467.1) where EcAII was included, suggesting that proteins clustered in this clade share similar properties to EcAII. In addition, two proteins, WP_053609500.1 and WP_055617501.1, exhibited the largest proportion of antigenic regions, with almost the same probability regions as the EcAII protein.

On the other hand, for the ASNases of PF06089.11, phylogenetic analysis included both the ASNase sequence of R. etli and S. coelicolor (ReAII and ScAII, respectively) (Fig. 1B). The tree defines two clades. In the first one, where the ScAII was included, we also considered ARZ68596.1 from S. albireticuli, SOD64826.1 from S. zhaozhouensis, WP_078645645 from S. varsoviensis, CDR15801.1 from S. iranensis, and WP_020554088 from S. scabrisporus. In the second clade were included the following proteins: WP_044373749 from S. ahygroscopicus and WP_078980718.1 from S. scabrisporus.

Antigenicity predictions

The results for antigenicity showed a likelihood of being antigenic for all ASNases in both sets that was lower than that of EcAII (Fig. 2). Nevertheless, among selected Streptomyces ASNases, the candidates from S. purpurogeneiscleroticus (WP_053609500.1) and S. phaeochromogenes (WP_055617501.1) showed a higher probability of being antigenic, whereas the rest of the ASNases showed very low antigenicity values in comparison with an E. coli ASNase (P00805_EcAII).

T-cell epitope analysis

After antigenicity prediction, the ED, the total number of high-affinity epitopes, the affinity epitopes, and the number of HLA alleles covered by each ASNase were calculated. The results showed that the ASNases with accession numbers WP_053609500.1, WP_053610569.1, EFL23513.1, WP_095730579.1, WP_078513220.1, and WP_052425051.1 have higher EDs than the reference (P00805_EcAII; ED=0.01114; 5 covered alleles) (Fig. 3).

On the other hand, the ASNase with the lowest predicted ED was WP_044373749.1, with an ED of 0.0027 and a coverage of 4 alleles, following by WP_095730579.1 (2 alleles), ELP65653.1 (3 alleles), and Q9K4F5 (3 alleles) (Table 2).

Additionally, the distribution of epitopes was mapped into the sequences of the ASNases (Fig. 3). ASNases of the PF06089.11 family tended to have a lower ED (Table 3) as well as lower allele coverage than those of the PF00710.11 family (Fig. 3).

Next, ASNases with lower allele coverage, lower ED, and lower probability of antigenicity were selected for further analysis. S. coelicolor (Q9K4F5), S. scabrisporus (WP_078980718.1), and S. albireticuli (ARZ68596.1) were selected as promising enzymes.

Protein structure predictions

From selected ASNases, homology-based models were generated (I-Tasser). For the subsequent analysis, the S. scabrisporus asparaginase II model, which had the highest C-value, was chosen (WP_078980718.1 SsAII-2) (Fig. 4). The structural model obtained by I-Tasser (with a C-value of −3.09) was refined with ModRefiner. In addition, the RAMPAGE program (http://mordred.bioc.cam.ac.uk/˜rapper) and Verify3D were used to validate the stereochemical quality of the resulting three-dimensional model. After analyzing the Ramachandran plot, 91.7% and 5.5% of the residues were located in favored and allowed regions, respectively; whereas Verify3D analysis revealed that 80.73% of the residues had an average 3D-1D score <0.2, indicating that the model is compatible with its sequence.

Based on the predicted structure, ASNase WP_0789718.1 (PF06089.11 family) is related in terms of folding to the beta-lactamase family. Beta-lactamases (SCOP data base, entry 56600) consist of a cluster of alpha-helices and an alpha/beta sandwich. This folding is also found in transpeptidases, esterases, penicillin receptors, D-aminopeptidases, and glutaminases (InterPro IPR012338).

Active site prediction

In order to identify the active site residues of the S. scabrisporus ASNase (WP_0789718.1), three approaches were used: genomic comparison, blind molecular coupling simulation, and search for high-affinity binding pockets with AutoLigand (active site). To our knowledge, there is no information regarding the active site of the family PF06089.11 ASNases, so genomic comparison was not possible. Using AutoLigand, two possible high affinity binding sites for L-asparagine were identified (Fig. 5A). The first (site A) had a volume of 121 ų and an energy per volume equal to −0.2149 kcal/mol ų; the second (site B) had a volume of 101 ų and an energy per volume equal to −0. 2136 kcal/mol ų. Site A is located between an alpha-helix in the amino terminal containing the ⁵⁷PRSx(2)KPxQ⁶⁵ motif, and a loop in the central region of the enzyme, containing the ¹⁴¹NCSGKHxAML¹⁵⁰ motif (Table 3). Site B is located in a pocket formed by a set of alpha-helices in the amino terminal of the protein, marked by the presence of the ⁸⁷SHTGQxHFV⁹⁵ motif. On the other hand, by performing AutoDock 4.2 whole-protein molecular coupling simulations, we found that the best ligand-enzyme interaction (L-asparagine-ASNase), with a binding free energy of −4.17 kcal/mol, targeted residues corresponding to the ¹⁴¹NCSGKHxAML¹⁵⁰ motif, which correspond to the site A (Fig. 5B).

Additionally, in order to validate AutoLigand analysis searching active sites in the S. scabrisporus ASNases, a search for binding sites in EcAII was performed. To do this, the monomeric, dimeric, and tetrameric forms of the enzyme (the latter is the catalytically active form) were analyzed using the same conditions used for SsAII-2. It was found that AutoLigand successfully identified the binding site of L-Asn, consisting of Thr 12, Tyr 25, Ser 58, Gln 59, Thr 89, Asp 90, and Lys 162 and also Asn 248 and Glu 283 (Fig. 6), the latter two only for dimeric and tetrameric forms. The sites found (purple squares curves) had energies by volume equal to −0.2119, −0.2242, and −0.2366 kcal/mol ų and a volume of 136, 122, and 102 ų for the monomer, dimer, and tetramer, respectively (Fig. 7). It is relevant that for both the dimeric and the tetramer forms, AutoLigand successfully identified L-Asn binding pockets in EcAII: the pocket formed between the amino-terminal end of subunit A and the carboxy terminal of the subunit C, as well as equivalent pockets for dimer BD. In addition, several other solutions found by AutoLigand (curve with blue or green squares), using up to 90 filling points, converge in the different joint pockets formed by dimers.

Molecular docking

Molecular docking simulations were performed at the putative sites found (Table 4). For EcAII, as the reference ASNase, Thr 12, Tyr 25, Ser 58, Gln 59, Thr 89, Asp 90, Asn 248, and Glu 283 were established as flexible residues; meanwhile, molecular docking for S. scabrisporus ASNase were performed using only the rigid structure of the protein, without defining flexible side chains for L-asparagine binding.

Our results showed a higher affinity for L-asparagine of the S. scabrisporus ASNase site A than site B; however, the affinity was lower than that for EcAII. For S. scabrisporus ASNase site A, the L-asparagine interacts with residues Ser 59, Lys 62, Asn 141, Ser 143, Lys 145, His 146, Gly 237, Lys 255, and Gly 256 (Fig. 8A); for site B, the residues that interact with L-asparagine are Ala 84, Gly 78, Ser 87, Tyr 163, Leu 164, and Asp165 (Fig. 8B). Interestingly, from site A, Lys 62, Asn 141, Ser 143, Lys 145, and His 146 are highly conserved in ASNases of the PF06089.11 family.