Stringent in-silico identification of putative G-protein-coupled receptors (GPCRs) of the entomopathogenic nematode Heterorhabditis bacteriophora

In order to predict putative GPCRs in H. bacteriophora and C. elegans, we followed a combination of different bioinformatic approaches, including both motif-based and alignment-free methods, as well as structural similarity–based approaches. The sequence alignment methods were comprised of BLASTX, BLASTP (Altschul et al., 1990), and Pfam (Bateman et al., 2000). The alignment-free method included GPCRPred (Bhasin and Raghava, 2004), whereas TMHMM (Krogh et al., 2001), GPCRHMM (Wistrand et al., 2006), and GPCRPipe (Theodoropoulou et al., 2013) come under machine learning or statistical algorithm-based methods. The entire pipeline for GPCR identification was composed of three key steps, as illustrated in Fig. 2.

Figure 2:

The initial step was to filter the large number of protein sequences based on their length. The proteomic dataset (Bai et al., 2013; Mclean et al., 2018) of H. bacteriophora and C. elegans (C. elegans Sequencing Consortium, 1998) containing 21,699, 15,218, and 28,447 predicted protein sequences, respectively, were filtered based on the length of sequences to discard proteins with less than 200 amino acid residues. These sequences were then checked for the number of transmembrane helices. Four different transmembrane detectors were used for this purpose: (1) TMHMM 2.0 (https://services.healthtech.dtu.dk/service.php?TMHMM-2.0), a membrane protein topology prediction method based on hidden Markov model (HMM) (Krogh et al., 2001); (2) Phobius (https://phobius.sbc.su.se/), an HMM-based signal peptide predictor that detects transmembrane protein topology also (Käll et al., 2004; Käll et al., 2007); (3) HMMTOP2 (http://www.enzim.hu/hmmtop/), an automatic server for predicting transmembrane helices and topology of proteins based on HMM (Tusnády and Simon, 1998; Tusnády and Simon, 2001); and (4) TOPCONS (https://topcons.cbr.su.se/), a membrane protein topology prediction tool based on a complex algorithm (Tsirigos et al., 2015). TMHMM2 was run using the “one line per protein” option. Phobius was run in the Normal prediction mode with the short output format mode selected. HMMTOP2 was run in the advanced mode with the parameters: FASTA format, Single Sequence type, Reliable prediction type, text output, and the results in one line. TOPCONS was operated with default parameters.

Seven transmembrane sequences, detected by any of the three transmembrane detectors, were then subjected to four different GPCR prediction tools to identify putative GPCRs from those sequences. These were GPCRHMM (https://gpcrhmm.sbc.su.se/index.html), an HMM-based GPCR prediction tool (Wistrand et al., 2006); GPCRPipe (http://aias.biol.uoa.gr/GPCRpipe/search.php), which detects GPCRs using pHMM library and HMM for GPCRs (Theodoropoulou et al., 2013); GPCR Prediction Ensemble (GPCR-Pen) (https://gpcr.utep.edu/), which determines GPCR based on four different algorithms (Begum et al., 2020); and GPCRPred (https://webs.iiitd.edu.in/raghava/gpcrpred/), a support vector machine (SVM) based approach (Bhasin and Raghava, 2004). In case of GPCRPipe, it detects GPCRs in two consecutive steps. In the first step, it uses Hidden Markov Model and then a library of 39 Pfam profile HMMs (specific to different families of GPCR). Analyses using GPCRHMM were performed with the local scoring option turned on. We used the GPCRPipe “AND” method as it allows prediction of GPCR only when it is confirmed by two methods. This selection resulted in a reduced number of GPCRs and a limited number of false positives. GPCRPred was run in its webserver without any manipulation in default parameters. While using GPCR-Pen, we restricted the prediction algorithms to Pfam and BLAST to get results only based on sequence similarity, and we fetched those sequences which were confirmed by both of these methods.

The outcomes of all the bioinformatic tools were then compiled and visualized using the Venny 2.1 software (https://bioinfogp.cnb.csic.es/tools/venny/). Venn diagrams were created to compare and to reflect the number of seven transmembrane sequences and predicted GPCRs made by different combinations of tools.

Apart from predicting GPCR, GPCRPred was also used to know the families and subfamilies of GPCRs based on the dipeptide composition of proteins. It categorizes GPCRs into class A (rhodopsin-like), B (secretin and adhesion), C (metabotropic glutamate), D (fungal pheromone receptors), E (cAMP receptors), or F (frizzled), as per International Union of Basic and Clinical Pharmacology (IUPHAR) nomenclature (http://www.guidetopharmacology.org/nomenclature.jsp; Attwood and Findlay, 1994; Kolakowski, 1994; Sadowski and Parish, 2003).

All GPCR candidates predicted by the pipeline were analyzed by another alignment-free method, GPCR-CA (http://218.65.61.89:8080/bioinfo/GPCR-CA), whose algorithm is based on cellular automation (CA) (Xiao et al., 2009). GPCR-CA is a bilayer predictor: the first layer identifies a query protein as GPCR or non-GPCR; if it is a GPCR, then a second layer classifies it into an A-F classification system, based on sequence homology and functional similarity (Attwood and Findlay, 1994).

In order to validate the pipeline, the putative GPCRs were then checked for the presence of an extracellular N terminus and an intracellular C terminus using HMMTOP2 and TMHMM2.

We screened all these GPCRs through Pfam (http://pfam.xfam.org/) (Punta et al., 2012), a database of protein families and an HMM-based (Eddy, 1998) protein classifier to get a complete and accurate classification of protein families and domains. Apart from this, InterPro (https://www.ebi.ac.uk/interpro/) was also used to analyze GPCR sequences functionally and to classify them into families, and ultimately to predict the presence of signature domains and important sites. We revalidated all these analyses through conserved domain analysis by NCBI CDD search (https://www.ncbi.nlm.nih.gov/Structure/cdd/wrpsb.cgi).

Putative H. bacteriophora GPCRs were classified into GPCR families according to the families to which their orthologous proteins in C. elegans and Ancylostoma caninum were assigned. For this purpose, orthology was detected between these resultant sequences and other organisms with the widely used reciprocal BLAST approach in WormBase Parasite (https://parasite.wormbase.org/index.html). This is done by using each GPCR sequence detected from H. bacteriophora as a query to the Caenorhabditis elegans proteome and selecting the highest scoring C. elegans protein sequence. This C. elegans sequence was then used as the query in a BLAST search against the proteome of H. bacteriophora. If the C. elegans sequence reciprocally identified the original query sequence as the highest scoring hit with the BLAST, these two sequences were then considered as an orthologous pair, which implies that they may possess similar functions and biological roles and fall into the same GPCR subfamily. H. bacteriophora GPCRs were assigned to families and subfamilies based on their orthologous GPCRs in other model organisms, such as C. elegans and A. caninum. In all BLAST comparisons, we used the BLOSUM62 scoring matrix and a cut-off threshold of 1e-1. Additionally, motif prediction of identified GPCRs was performed in MEME Suite 5.4.1 (https://meme-suite.org/meme/). Motifs were functionally annotated in HHpred (http://toolkit.tuebingen.mpg.de/hhpred) using default parameters. Functional motifs were identified in the GPCR sequences using the MOTIF search tool (https://www.genome.jp/tools/motif/).

PRED-COUPLE2 (http://athina.biol.uoa.gr/bioinformatics/PRED-COUPLE2/) (Sgourakis et al., 2005) was used to predict the coupling specificity of GPCRs to the four families of G-proteins, with a stringent cutoff of 0.3 to discriminate between positive and negative predictions. Therefore, results below this limit were not considered as positive predictions and were discarded.

Out of the total 21,699 predicted proteins of H. bacteriophora, 21 sequences were identified as GPCRs using the pipeline presented in Fig. 2. At the initial stage of identification, 7,613 sequences with more than 200 amino acid residues were short-listed. Four different transmembrane predictors identified varying numbers of sequences containing transmembrane domains. TMHMM2, a widely used transmembrane predictor, identified 97 sequences containing 7-TM helices out of 7,613 sequences, of which 14 7-TM sequences were uniquely identified by TMHMM2 (Fig. 3). Another tool, Phobius, detected 139 7-TM sequences, of which 42 were unique (Fig. 3). In addition, 153 sequences were predicted by HMMTOP2, with 59 unique sequences. TOPCONS could identify only 123 proteins containing 7-TM helices, but only 25% of those were uniquely identified. A total of 69 sequences were predicted to have 7-TM helices by at least three out of four transmembrane predictors used for this study (Fig. 3A, supplementary Table 1).

Figure 3:

The common 69 heptathetical transmembrane sequences were passed through four GPCR prediction tools to identify the GPCRs among these sequences. GPCRHMM, GPCRPipe, GPCR-Pen, and GPCRPred identified 22, 21, 15, and 62 sequences as GPCRs, respectively (Fig. 3B). Interestingly, GPCRHMM, GPCRPipe, and GPCR-Pen did not identify any unique GPCR sequence, but GPCRPred identified 36 unique GPCRs not predicted by any other tool. In summary, 21 sequences were considered as putative GPCRs as they were identified by at least three of the four tools used for GPCR prediction, and further confirmed by detailed analysis of the sequences (Fig. 3B, Tables 1–3, supplementary Tables 2 – 4).

Using the same pipeline, 27 GPCR sequences were identified from the proteomic dataset of H. bacteriophora, published by McLean and coworkers in 2018. Out of these sequences, 13 were found to be similar, with the GPCRs fetched from the earlier version of annotation, published in 2013 (Table 4).

A total 1,252 GPCRs (supplementary Table 5) were identified out of the 28,447 predicted proteins of C. elegans using the pipeline proposed in this study (Fig. 2). Initially 21,372 sequences with >200 amino acid residues were screened for 7-TM domains. Out of four different transmembrane predictors, Phobius detected the highest number (1782) of 7-TM sequences, with 116 unique ones. TMHMM2 identified 1,566 sequences containing 7-TM helices, of which 80 sequences were unique. In addition, 1,710 sequences were predicted by HMMTOP2, with the highest number (178) of unique sequences. TOPCONS was able to identify only 791 protein coding genes with 7-TM helices, but <1% of those were unique. A total of 1,404 sequences were predicted to have 7-TM helices by at least three out of four tools (Fig. 3C). The 1,404 resultant sequences were screened using four GPCR predictors. GPCRHMM, GPCRPipe, GPCR-Pen, and GPCRPred identified 1,255, 1,393, 258 and 1,388 sequences as GPCRs, respectively (Fig. 3D). Interestingly, GPCRHMM and GPCR-Pen could not identify any unique GPCR sequence. Finally, 1,252 sequences were considered putative GPCRs as they were identified by at least three of the four GPCR detectors (Fig. 3D, supplementary Table 5).

All the sequences identified as GPCRs above were validated by the GPCR-CA tool, as it uses Cellular Automaton (CA) images to reveal the features hidden in complex protein sequences. It designated all these GPCRs as Class A Rhodopsin-like GPCR, except Hba_17528, which was identified as Class D Fungal pheromone GPCR (Table 1). Further, the predicted GPCRs were classified based on function and protein family, and four chemosensory GPCRs (Hba_07805, Hba_18427, Hba_18743 and Hba_17528), ten 7-transmembrane receptors under the rhodopsin family, one rhodopsin-like GPCR with transmembrane domain, one frizzled type, and one secretin type GPCR were identified (Table 1, supplementary Table 1). In addition, the reciprocal BLAST was used to find the orthologous sequences of the predicted GPCRs from the proteomic dataset of closely related model organisms. This approach also validated the pipeline and confirmed that all the identified protein sequences were GPCRs. Sequence similarities found a frizzled type (Hba_19080), a secretin type (Hba_20566) and the other 19 as rhodopsin types of GPCR (Table 1). The GPCRPred tool classified all the 21 shortlisted sequences as Class A rhodopsin-like GPCR and further classified them into three different subfamilies: peptide (15 sequences), biogenic amine (3 sequences), and lysospingolipid (1 sequence) (Table 1). A search against InterPro database confirmed the families of all the 21 proteins as G protein-coupled receptors except Hba_17528. All the GPCRs were suggested to be involved in G protein-coupled receptor signalling pathway except Hba_19080, which was annotated to have a role in cell surface receptor signalling. Additionally, Hba_18878 was suggested to be involved in the neuropeptide signalling pathway, and Hba_18906 in pheromone responsiveness. Most of them are rhodopsin-like GPCR, except Hba_19080 (frizzled/secreted, frizzled-related protein) and Hba_20566 (secretin-like). All these GPCRs were found to be integral components of the cell membrane by gene ontology (GO) analysis, except Hba_18906 (Table 2).

Additionally, NCBI conserved domain analysis revealed that four identified proteins were members of the class A seven-transmembrane GPCRs and belonged to FMRFamide (Phe-Met-Arg-Phe)-like receptors and related proteins. Eight of the proteins were rhodopsin receptor-like class A family of the seven-transmembrane GPCR superfamily, which constitutes about 90% of all GPCRs. They include light-sensitive rhodopsin and receptors for biogenic amines, lipids, nucleotides, odorants, peptide hormones, and various other ligands. Six of the proteins were broadly classified under the 7tm_GPCRs superfamily. Among these, Hba_18743 was a serpentine type 7TM GPCR chemoreceptor under the Srsx family, the only family among the various superfamilies of chemoreceptors. Another serpentine type of chemoreceptor GPCR (Hba_17528) was found under the srx family, which is a part of the Srg superfamily of chemoreceptors. Interestingly, Hba_20566 was a pigment-dispersing factor receptor (PDFR), a member of the B1 subfamily of class B seven-transmembrane GPCRs, also referred to as the secretin-like receptor family. Hba_18203 was found to be an FMRFamide receptor and a member of the class A family of seven-transmembrane G protein-coupled receptors. Hba_09978 was a cholecystokinin receptor and came under the class A family of seven-transmembrane GPCRs. This group represents four GPCRs that are members of the RFamide receptor family, including cholecystokinin receptors (CCK-AR and CCK-BR), orexin receptors (OXR), neuropeptide FF receptors (NPFFR), and pyroglutamylated RFamide peptide receptors (QRFPR). Hba_08130 was an amine receptor of the class A family of GPCRs, which include adrenoceptors, 5-HT (serotonin) receptors, muscarinic cholinergic receptors, dopamine receptors, histamine receptors, and trace amine receptors (supplementary Table 2).

Lastly, in the 21 GPCR sequences, 18 different motifs were identified, which include 7-transmembrane receptor (rhodopsin family), serpentine type 7TM GPCR chemoreceptor Srw, serpentine type 7TM GPCR chemoreceptor Srsx, frizzled/smoothened family membrane region, 7-transmembrane receptor (secretin family), serpentine type 7TM GPCR chemoreceptor Srx, serpentine type 7TM GPCR chemoreceptor Srt, frizzled/smoothened family membrane region etc (Table 3, Fig. 4).

Figure 4:

PRED-COUPLE2 analysis showed possible interaction of 21 GPCRs to the different families (Gi/o-, Gq/11-, Gs- and G_12/13) of G-proteins (Table 1). Among these, Gi/o is the most abundant kind of G-protein, which can bind to all the fetched GPCRs. There are few GPCRs which have coupling specificity only to Gi/o. Additionally, there are promiscuous GPCRs, except Hba_14891, Hba_19161, Hba_20566, and Hba_17528, which can couple to members of more than one G-protein subfamily. Hba_18878, Hba_13948, and Hba_19080 are the only three GPCRs which can interact with all the four types of G-proteins. No coupling specificity is observed with any of the G-proteins; only in the case of Hba_20096. Though Hba_18878 and Hba_19080 have very low sequence similarity, they have coupling ability to members of the same subfamily of G-proteins. While Hba_10668 and Hba_14446 belong to the same amine subfamilies of rhodopsin like GPCR, they often couple to members of distinct G-protein subfamilies. Several sequences under the peptide subfamily have considerably high sequence similarity but they differ in coupling specificity with different families of G-proteins, except G_i/o (Table 1).