Phylogenetic Analysis Indicates That Evasin-Like Proteins of Ixodid Ticks Fall Into Three Distinct Classes

Chemokines are structurally related proteins that activate leucocyte migration in response to injury or infection. Tick saliva contains chemokine-binding proteins or evasins which likely neutralize host chemokine function and inflammation. Biochemical characterisation of 50 evasins from Ixodes, Amblyomma and Rhipicephalus shows that they fall into two functional classes, A and B, with exclusive binding to either CC- or CXC- chemokines, respectively. Class A evasins, EVA1 and EVA4 have a four-disulfide-bonded core, whereas the class B evasin EVA3 has a three-disulfide-bonded “knottin” structure. All 29 class B evasins have six cysteine residues conserved with EVA3, arrangement of which defines a Cys6-motif. Nineteen of 21 class A evasins have eight cysteine residues conserved with EVA1/EVA4, the arrangement of which defines a Cys8-motif. Two class A evasins from Ixodes (IRI01, IHO01) have less than eight cysteines. Many evasin-like proteins have been identified in tick salivary transcriptomes, but their phylogenetic relationship with respect to biochemically characterized evasins is not clear. Here, using BLAST searches of tick transcriptomes with biochemically characterized evasins, we identify 292 class A and 157 class B evasins and evasin-like proteins from Prostriate (Ixodes), and Metastriate (Amblyomma, Dermacentor, Hyalomma, Rhipicephalus) ticks. Phylogenetic analysis shows that class A evasins/evasin-like proteins segregate into two classes, A1 and A2. Class A1 members are exclusive to Metastriate ticks and typically have a Cys8-motif and include EVA1 and EVA4. Class A2 members are exclusive to Prostriate ticks, lack the Cys8-motif, and include IHO01 and IRI01. Class B evasins/evasin-like proteins are present in both Prostriate and Metastriate lineages, typically have a Cys6-motif, and include EVA3. Most evasins/evasin-like proteins in Metastriate ticks belong to class A1, whereas in Prostriate species they are predominantly class B. In keeping with this, the majority of biochemically characterized Metastriate evasins bind CC-chemokines, whereas the majority of Prostriate evasins bind CXC-chemokines. While the origin of the structurally dissimilar classes A1 and A2 is yet unresolved, these results suggest that class B evasin-like proteins arose before the divergence of Prostriate and Metastriate lineages and likely functioned to neutralize CXC-chemokines and support blood feeding.


INTRODUCTION
Chemokines are a family of structurally related proteins that have prominent roles in driving leucocyte recruitment -a process that plays a key role not only in the inflammatory response to injury or infection (predominantly by neutrophils and monocyte/ macrophages), but also in immune system homeostasis mediated by trafficking of lymphocytes and dendritic cells (Zlotnik and Yoshie, 2012). The >45 mammalian chemokines are grouped according to the spacing of their N-terminal cysteine residues into CC, CXC, CX3C and XC classes. They bind to a family of 19 G-protein coupled receptors (GPCRs) that are expressed on leucocytes. Activation of GPCR-signalling leads to directed migration or chemotaxis of leucocytes to the site of chemokine expression, for instance, to the site of infection or injury, where they form the first line of defence. The chemokine network is highly redundant with multiple chemokines expressed at the site of injury, many chemokines binding multiple receptor subtypes, and multiple receptor subtypes typically expressed on a single leucocyte subclass (Mantovani, 2018;Bhattacharya and Kawamura, 2020). This redundancy is thought to result in a robust output that is resistant to genetic or environmental variation (Mantovani, 2018).
The site of a tick bite is characterized by early expression of several CC and CXC-class chemokines (Glatz et al., 2017). Following Ixodes ricinus bite in humans, elevated levels of CCL2, CCL3, CCL4, CXCL1 and CXCL8 are observed, and there is early accumulation of macrophages and dendritic cells (Glatz et al., 2017). Following Ixodes scapularis bites in mice, there is local expression of CCL2, CCL7, CXCL2, CXCL5, but not CXCL3, and accumulation of neutrophils, eosinophils, lymphocytes and histiocytes (i.e., tissue monocyte/macrophages) (Krause et al., 2009;Heinze et al., 2012). However, the area immediately adjacent to the tick hypostome is paucicellular -i.e. minimally inflamed, and this is thought to be due to anti-inflammatory effects of tick salivary proteins (Krause et al., 2009). Reinfestation is associated with infiltration of the paucicellular zone with neutrophils, eosinophils, lymphocytes and histiocytes, factors that cause premature tick detachment (Krause et al., 2009). Indeed, basophils appear to have a critical role in acquired resistance to tick reinfestation (Wada et al., 2010).
The nature of the biological arsenal with which ticks may defend against redundant chemokine-induced inflammation has become clearer following the initial report of chemokine-binding proteins in tick saliva (Hajnickáet al., 2001). Ticks produce salivary chemokine-binding proteins that bind and neutralize multiple chemokineswhich would be expected to minimize local inflammation and allow prolonged blood feeding. Three proteins -evasins 1, 3 and 4 (EVA1, EVA3, EVA4) -were cloned from the salivary gland of the dog tick Rhipicephalus sanguineus (Frauenschuh et al., 2007;Deruaz et al., 2008). EVA1 and EVA4 specifically bind CC chemokines, whereas EVA3 specifically binds CXC chemokines. Structural characterization of the CCbinding evasins, EVA1 and EVA4, indicates a four-disulfidebonded core structure (Dias et al., 2009;Denisov et al., 2020). Characterization of the CXC-chemokine-binding evasin EVA3 indicates that it forms a three disulfide-bonded "knottin" motif (Denisov et al., 2019;Lee et al., 2019). While "knottin" structures are reported in multiple organisms, evasins appear to be unique to ticks as sequence homologs have not been reported for any other organisms to date.
In addition to inflammation inhibitors, tick saliva is known to contain inhibitors of blood clotting, complement activation, pain, and itchmany of which have potential therapeutic and vaccine development applications (Simo et al., 2017;Aounallah et al., 2020). The advent of cost-effective RNA sequencing has fuelled a veritable explosion of sequence data from the salivary transcriptomes of a number of tick species (Ribeiro and Mans, 2020). This serves as a valuable resource for identifying novel modulators of human biology that have appeared in ticks through natural selection.
The 50 validated evasins have been distinguished into two classes, class A evasins (CC-chemokine binders) and class B evasins (CXC-chemokine binders), that each share structural and conserved sequence features although their sequence identity may be <30% (Bhusal et al., 2020). Phylogenetic sequence relationships of these biochemically characterized evasins to the many evasin-like proteins annotated in tick transcriptome databases are unclear. Here, using BLAST searches of tick transcriptomes supported by the 50 biochemically characterized evasins (defined here as possessing a chemokine-binding function), followed by sequence-based phylogenetic tree construction, we identify evasin-like proteins (defined here as ones not yet tested for a chemokine-binding function) and show they fall into three distinct classes.

Databases
The transcriptome shotgun assembly database (tsa_nr) and the non-redundant protein database (nr) was downloaded from NCBI (May 2021).

BLAST Search
Databases were searched using command-line BLAST+ blastp (BLAST/2.10.1-Linux_x86_64) on the Oxford BMRC cluster with each of the 21 biochemically characterized class A and 29 biochemically characterized class B evasins (Supplementary  Tables S1, S2, respectively), using default values for blastp and restricting the search to ticks, to retrieve NCBI accession numbers, genus, and species. Distinct sequences with e-value <0.0001 were used for subsequent analysis.

Identification of Mature Protein Sequence
Subsequent steps were performed with customized R scripts. NCBI accession numbers identified by blastp were used to retrieve the corresponding protein sequences using "entrez_fetch" (package rentrez). The first methionine residue in each sequence was used to identify the protein start site. A local installation of the program SignalP 4.1 (downloaded from https://services.healthtech.dtu.dk) was used to identify the signal peptide (Bendtsen et al., 2004). The start of the mature peptide was identified using the parameter posYmax. Sequences lacking a signal peptide and sequences containing stop codons in the predicted mature protein sequence were filtered out. Proteins identified are shown in Supplementary Table S1 (class A evasins/evasin-like proteins) and Supplementary Table S2 (class B evasins/evasin-like proteins).

RESULTS
Sequences of the 50 biochemically characterized evasins (21 in class A and 29 in class B) were used to screen the BLAST databases, which consisted of 42586 argasid, 204166 Metastriate ixodid and 170717 Prostriate ixodid tick sequences (Supplementary Figure  S1A). Using a BLAST e-value threshold of p< 0.0001, we identified 271 class A and 128 class B evasin-like proteins from Prostriate (Ixodes) and Metastriate (Amblyomma, Dermacentor, Hyalomma, Rhipicephalus) ixodid ticks (Supplementary Tables S1, S2). No evasin-like proteins were identified for argasid ticks. A possible reason is that our search criteria are too stringent to identify very distantly related sequences.
The evolutionary relationships between sequences and groups of sequences can be inferred from phylogenetic trees that are rooted in an explicit ancestral node (Wilkinson et al., 2007). To understand the evolutionary history and relationships between the 292 class A and 157 class B evasin and evasin-like proteins of Prostriate and Metastriate ticks, we performed sequence alignment followed by phylogenetic tree construction using the midpoint rooting method (Hess and De Moraes Russo, 2007;Kinene et al., 2016) to locate the ancestral nodes (Figures 1, 2). We were unable to use the alternative outgroup rooting method as no appropriate outgroup sequence, i.e., an evasin from a species that diverged prior to the separation of Prostriate and Metastriate lineages could be identified.
Inspection of the tree created from class A evasins/evasinlike proteins indicates that the two nodes immediately following the ancestral node clearly separate Prostriate and Metastriate lineages (Figure 1 and Supplementary Table S1). We classified the two clades descending from these nodes as class A1 (Metastriate, 205 members) and class A2 (Prostriate, 87 members) respectively. As shown in the heatmap accompanying the tree, class A1 evasins/evasin-like proteins typically have 8 or more cysteine residues that are appropriately spaced to create a Cys8-motif. Class A2 members typically have less than 8 cysteine residues and lack the Cys8-motif. As also seen in the heatmap, several class A1 evasins have less than the minimal number of Cys residues needed to form the four-disulfide bonded tertiary structure, which is thought to play an important role in chemokine binding (Dias et al., 2009;Denisov et al., 2020).
Inspection of the tree created from class B evasins/evasin-like proteins shows that, unlike the situation with class A evasins, there is incomplete separation between proteins from Prostriate and Metastriate lineages (Figure 2 and Supplementary Table S2). The two nodes immediately following the ancestral node define clades B1 and B2. B1 contains proteins from both Prostriate and Metastriate lineages, whereas B2 is exclusively Prostriate. Within B1 there are two subclades, B1.1 being exclusively Prostriate, and the other, B1.2, exclusively Metastriate. As shown in the heatmap, class B evasin-like proteins from B1.1, B1.2 and B2 clades typically have 6 cysteine residues appropriately spaced to create the Cys6-motif, and can be considered as homologous (Figure 2 and Supplementary Table S2). For subsequent analyses we therefore considered class B proteins as a single class. As also seen in the heatmap, several class B evasins have less than the minimal number of Cys residues needed to form the three-disulfide bonded tertiary structure, which plays an important role in chemokine binding (Lee et al., 2019).
Using the output of the above analyses, we next estimated the relative numbers of class A1, A2 and B evasins/evasin-like sequences in Prostriate and Metastriate lineages, and in the species constituting these lineages (Figure 3). Class B evasins/ evasin-like sequences are identified in two Prostriate species (Ixodes ricinus, and Ixodes scapularis) and four Metastriate genera (Rhipicephalus, Hyalomma, Dermacentor and Amblyomma). Class A1 evasins/evasin-like sequences are identified in four Metastriate genera (Rhipicephalus, Hyalomma, Dermacentor and Amblyomma). Class A2 evasins/evasin-like sequences are identified in three Prostriate species (Ixodes holocyclus, I. ricinus, and I. scapularis). The majority of evasins/evasin-like proteins in Prostriate lineages is class B, whereas in Metastriate lineages it is class A1. In keeping with this distinction, the majority of Metastriate evasins characterised to date neutralize CCchemokines, whereas the majority of Prostriate evasins neutralize CXC-chemokines (Figure 4).

DISCUSSION
In this study we have used systematic BLAST searches with uniform selection criteria to identify evasin-like proteins from publicly available transcriptome databases. Comparison with annotations in the NCBI databases (Supplementary Figure S1B) reveal several differences to our analysis. For instance, our analysis indicates there are 195 evasin/evasin-like proteins for I. ricinus and 5 evasin/evasinlike proteins for I. holocyclus (Figure 3). In contrast, the annotations in the NCBI databases suggest that there are 69 for I. ricinus and FIGURE 1 | Neighbour-joining tree of class A evasin and evasin-like proteins. Accession numbers of proteins are indicated at the tree tips, which are colour coded by species. The ancestral node at the centre was defined by midpoint rooting. The outermost circle of the heatmap indicates Prostriate or Metastriate lineage. The next two circles of the heatmap indicate the presence of a C-x(14,17) -C-x(3)-C-x(11,16)-C-x(17,20)-C-x(4)-C-x(4,5)-C-x(8)-C (Cys8-motif) and the numbers of cysteines respectively. The innermost circle of the heatmap indicates the number of CC-chemokines bound by biochemically characterized evasins. Classes A1 and A2 were defined as having diverged from the two nodes immediately following the ancestral node and are indicated in orange and yellow respectively. Distinct subclades of class A1 containing EVA1 and EVA4 proteins are indicated in red and dark orange respectively.  (Figure 3), whereas annotations in the NCBI databases suggest that there are 37 and 11 respectively (Supplementary Figure S1B). The differences with NCBI annotations likely reflect the application of our search criteria and differences in analytical methodology. Specifically, we performed searches using all 50 biochemically characterised evasins and restricted the output to proteins that had a BLAST e-value <0.0001, contained a signal peptide, and had no stop codons within the predicted mature protein sequence. The differences between species may also arise partly from differences in the transcript count per species in the NCBI database (Supplementary Figure S1A). Moreover, composition of the salivary gland transcriptome varies according to the origin of the tick (geographical location, laboratory colony, etc.), developmental stage and gender, whether from fed, feeding, or unfed ticks, the time of feeding relative to attachment, and the host on which the tick is feeding (De Castro et al., 2016;Garcia et al., 2020;Kim et al., 2020). None of these variables are standardised for salivary gland transcriptomes hence, together with differences in transcript quality and quantity, it is perhaps not surprising we found an approximately 100-fold range in number of evasin-like proteins between different species. Also of note are differences between evasin-like proteins identified in salivary gland transcriptomes and in saliva proteomes. For example, putative evasins were not recorded in the saliva proteomes of R. microplus (Tirloni et al., 2014;Garcia et al., 2020). Their apparent absence may be explained by their relatively low expression levels (typically <10%) (Karim et al., 2011;De Castro et al., 2016;Garcia et al., 2020). However, 12 FIGURE 2 | Neighbour-joining tree of class B evasin and evasin-like proteins. Accession numbers of proteins are indicated at the tree tips, which are colour coded by species. The ancestral node at the centre was defined by midpoint rooting. The outer circle of the heatmap indicates Prostriate or Metastriate lineage. The next two circles of the heatmap indicate the presence of a C-x(3)-C-x(6,10)-C-x(3,6)-C-x(1)-C-x(10,11)-C motif (Cys6-motif) and the numbers of cysteines respectively. The innermost circle of the heatmap indicates the number of CXC-chemokines bound by biochemically characterized evasins. Clades B1 and B2 were defined as having diverged from the two nodes immediately following the ancestral node and are indicated in azure and cyan respectively. Subclades B1.1 and B1.2 were defined as having diverged from the two nodes immediately following the node defining clade B1 and are indicated in light and dark azure respectively. putative evasins were annotated in the saliva proteome (and 43 in the salivary gland transcriptome) of A. americanum, providing supporting evidence that putative evasins are secreted proteins . Differences in the numbers of reported putative evasins highlight the need to standardise source material and analytical approaches to unlock the wealth of information in tick salivary gland transcriptomes. Phylogenetic tree construction using evasins/evasin-like proteins showed that class A evasins are separated into two classes A1 and A2 which are exclusive to Metastriate and Prostriate lineages respectively. Classes A1 and A2 are structurally dissimilar, with members of class A1 typically having the capacity of forming a four-disulfide-bonded structure, and members of class A2 typically lacking such capacity. The presence of EVA1 and EVA4 in distinct subclades of A1 is intriguing as both have a common fourdisulfide bonded core structure. However, there are subtle structural differences between the two structures. The EVA4   structure is bulkier and interacts with chemokines via the Nterminal domain, unlike EVA1 which binds through both N and C-terminal domains (Denisov et al., 2020). It is possible that the two subclades reflect these subtle differences. Supporting this idea, E672, a member of the EVA4 subclade, also interacts with chemokines through the N-terminal domain Darlot et al., 2020). In contrast to class A, class B evasins are present in both Prostriate and Metastriate lineages, and members typically have the capacity to form a three-disulfide-bonded "knottin" structure irrespective of lineage. Notably, several class A1 and class B evasins have less than the minimal number of Cys residues needed to form the appropriate Cys-Cys bonded tertiary structure, which is thought to play an important role in chemokine binding (Dias et al., 2009;Lee et al., 2019). It remains to be determined whether such evasinlike proteins are capable of binding chemokines. If neutral evolution and genetic drift are the major reasons for diversity in salivary gland proteins (Mans et al., 2017), many evasin-like proteins may be sub-or non-functional or have new (nonevasin) functionsa hypothesis that needs to be tested. Taken together these results suggest that class B evasins, which are identified in four Metastriate genera and two Prostriate species, arose prior to the divergence of Prostriate and Metastriate lineages, which occurred~200-249 MYA (Mans et al., 2012;Jia et al., 2020), and likely functioned to neutralize CXC-chemokines and support blood feeding. The presence of class A1 in four Metastriate genera suggests that these proteins appeared before Metastriate diversification, which occurred~138 MYA (Jia et al., 2020). The presence of class A2 in three Prostriate species suggests that these proteins appeared before Prostriate diversification, which occurred 217 MYA (Mans et al., 2012). Our analysis indicates that classes A1 and A2 are structurally dissimilar and are exclusive to Metastriate and Prostriate ticks respectively. Possible explanations for this exclusivity include: a) descent with modification from a shared ancestral protein sequence present before the divergence of Metastriate and Prostriate lineages, b) descent from independent ancestral protein sequences with loss of a class [i.e. "gene death" (Nei and Rooney, 2005)] following Metastriate/Prostriate divergence, or c) origin by horizontal gene transfer from another organisma phenomenon that has been shown to occur in arthropod lineages and has functional importance at the tick host interface (Chou et al., 2015;Hayes et al., 2020). These alternative scenarios may be resolved by further analysis of ticks from different lineages.

CONCLUSIONS
Our analysis of publicly available tick transcriptomes indicates that evasins/evasin-like proteins of ixodid ticks fall into three distinct classes. Class A1 is exclusive to Metastriate ticks. Members typically contain a characteristic Cys8-motif that is predicted to form a fourdisulfide-bonded structure. Two subclades containing EVA1 and EVA4 are evident and are likely to have subtle structural differences. Characterised members exclusively bind CC-chemokines. Class A2 is exclusive to Prostriate ticks, members do not contain a Cys8motif, and characterised members bind CC-chemokines. Class B is present in both Prostriate and Metastriate ticks, members typically contain a characteristic Cys6-motif predicted to form a threedisulfide-bonded "knottin" structure, and characterised members exclusively bind CXC-chemokines. Our analysis suggests that class B evasin-like proteins arose before the divergence of Prostriate and Metastriate lineages and likely functioned to neutralize CXCchemokines and support blood feeding. The origin of evasin classes A1 and A2 is yet unclear but may be resolved by transcriptome sequencing from underrepresented tick species. Together with biochemical characterisation of individual proteins such studies will provide further insights into the biological functions of tick evasins, the development of therapeutics that target the chemokine system, and the phylogenetic relationships of these "endless forms most beautiful" (Darwin, 1859).

DATA AVAILABILITY STATEMENT
The original contributions presented in the study are included in the article/Supplementary Material. Further inquiries can be directed to the corresponding author.

AUTHOR CONTRIBUTIONS
SB performed the data analysis. SB and PN wrote the manuscript. All authors contributed to the article and approved the submitted version. Supplementary Table 1 | Class A evasins/evasin-like proteins used to construct the tree in Figure 1. Table shows the identity, mature protein sequences, full protein sequences, genus, species, cysteine count in the mature protein, and the presence of a Cys8 or Cys6 motif.