Computational Analysis of the 3-D structure of Human GPR87 Protein: Implications for Structure-Based Drug Design

G protein-coupled receptors (GPCRs) are still a large superfamily of membrane bound signaling proteins that hold great pharmaceutical interest as a promising drug target (Costanzi et al., 2012). GPCRs accounts that more than 2% of all genes encoded by the human genome. According to a recent analysis, 40-50% of all current drug targets, with 26 out of the 100 top selling drugs are exposed by targeting different GPCRs family (Hans et al., 2007). GPCRs is also the best target for inflammatory mediators, therefore it’s provides a link between chronic inflammation and cancer disease. In addition, GPCRs intellectual have most important roles in tumor-induced angiogenesis and metastasis, which might involve in GPCR-guided migration of cancer cells to their target organs. Therefore, interfering with GPCRs and their downstream targets its potency provides an opportunity to develop a new technique for cancer diagnosis, prevention and treatment (Robert et al., 2007). GPCRs serve as a key signal transduction conduits by connecting extracellular inputs with diverse cellular responses (Audet et al., 2012). Although molecular biology and bioinformatics approaches builds the identification of orphan GPCRs


Introduction
G protein-coupled receptors (GPCRs) are still a large superfamily of membrane bound signaling proteins that hold great pharmaceutical interest as a promising drug target (Costanzi et al., 2012).GPCRs accounts that more than 2% of all genes encoded by the human genome.According to a recent analysis, 40-50% of all current drug targets, with 26 out of the 100 top selling drugs are exposed by targeting different GPCRs family (Hans et al., 2007).GPCRs is also the best target for inflammatory mediators, therefore it's provides a link between chronic inflammation and cancer disease.In addition, GPCRs intellectual have most important roles in tumor-induced angiogenesis and metastasis, which might involve in GPCR-guided migration of cancer cells to their target organs.Therefore, interfering with GPCRs and their downstream targets its potency provides an opportunity to develop a new technique for cancer diagnosis, prevention and treatment (Robert et al., 2007).GPCRs serve as a key signal transduction conduits by connecting extracellular inputs with diverse cellular responses (Audet et al., 2012).Although molecular biology and bioinformatics approaches builds the identification of orphan GPCRs 1 Department of Pharmacology and Therapeutics, King George's Medical University, (Erstwhile C.S.M. Medical University), Lucknow, Mukta Rani 1 , Anuradha Nischal 1 , Ganesh Chandra Sahoo 2 , Sanjay Khattri 1 * which are responsible for the search of their endogenous ligands which has been a significant challenge.The invention of endogenous ligands has given birth to the reverse pharmacology approach, which uses orphan GPCRs as a target to identify their ligands.This technique is magnificently used for the deorphanization of about 300 GPCRs over two decades (Chung et al., 2007).Due to their characteristic structure and specific ligand binding ability, many GPCRs have been extensively used as a target for therapeutic drug development and designing.Cancer is a deregulated multiplication of cells with the consequence of an abnormal increase of the cell number in particular organs or tissues.At a later stage, these cells are distributed via the hematopoietic and lymphatic systems throughout the body where they can colonize in distant tissues and cause metastasis.
A drug which has been launched in the market has consent through the seven basic steps are disease selection, target selection, lead compound identification, lead optimization, preclinical trials, clinical trial testing and pharmacogenomic optimization (Roderick, 2006).Orphan GPCR like GPR87 (also known as GPR95) is recently discovered orphan receptor and classified in the P2Y 12 subgroup that contains P2Y 12 , P2Y 13 , P2Y 14 , CysLT 1 , and CysLT 2 receptors (Nonaka et al., 2005).GPR87 is located on human chromosome 3q25 region (Wittenberger et al., 2001) in which the given P2Y receptor genes are clustered.GPR87 are found in neoplastic cells of different location with lung carcinomas (squamous), cervical squamous carcinomas, head and neck squamous cell carcinomas (larynx, pharynx, tonsil and tongue), urinary bladder carcinomas (Glatt et al., 2008;Gugger et al., 2008) and also in the placenta (Wittenberger et al., 2001).Human G-protein coupled receptor 87 (hGPR87) plays a crucial role in intracellular signal transduction (Wittenberger et al., 2001) and normal ovarian development as well as ovarian cancer progression.It was shown to be overexpressed in squamous cell carcinoma (SCCs) or adenocarcinoma in the lung and bladder carcinomas (Gugger et al., 2008).The two most promising reason for targeting this protein are first-GPR87 is the only example of the over expressed GPCRs that can be correlated on a mutation-based level to squamous cell carcinoma of the lung.Second-it shows the highest fold level of the over expressed GPCR.GPR87 is a novel target gene of p53 it's over expression in testicular cancer as a probable downstream effect of chemotherapy-induced p53 (Kcrley-Haniilton et al., 2005).Orphan GPR87 plays crucial role in SCCs and it will certainly lead to the identification of novel cancer therapeutic targets.Recent research confirm that the GPR87 is an LPA receptor (Ochiai et al. 2013).
Till date no crystal structure information either from X-ray crystallography or NMR is available for hGPR87 protein are reported.In the present study, different homology models of GPR87 protein with 358 amino acids was constructed by using multiple templates through different homology modeling servers.According to the Mobarec et al. reported that, use of multiple templates for modeling of GPCRs increases its model quality, in terms of higher predictive power in flexible ligand-rigid protein docking experiments (Mobarec et al., 2009).Molecular dynamics (MD) simulations studies with explicit solvent treatments which have been most frequently employed to assess the stability of protein.In silico modification of some substrate ligands has been carried out to search for better inhibitors of GPR87 protein.In the absence of crystal structure, our study provides an early insight into the detailed 3Dstructure of a major drug target or developing new inhibitors for GPR87 which will clearly be a challenging work in the near future.

Sequence alignment
The amino acid sequences of GPR87 from Homo sapiens (GenBank ID: Q9BY21) was retrieved from National Centre for Biotechnology Information (NCBI).Multiple sequence alignment between P2Y receptors, group I: P2Y 1 , P2Y 2 , P2Y 4 , P2Y 6 and P2Y 11 ; group II: P2Y 12 , P2Y 13 and some other pharmacologically related receptors like CysLT 1 and CysLT 2 were retrieved from NCBI.Alignment was performed by using the ClustalW algorithm and cladogram was assessed by TreeView algorithm (Larkin et al., 2007) and maximum likelihood method from www.phylogeny.fr(Dereeper et al., 2008).
Translated protein sequence of GPR87 protein has 358 amino acid was taken for further study.

Molecular modeling
BLASTP (Altschul et al., 1990) search was performed against Protein Data Bank (PDB) to find suitable templates for homology modeling.There was no single template available in PDB that could be used to model GPR87, therefore two templates were selected.Based on the maximum identity with high score and lower e-value crystal structure of squid rhodopsin by intracellularly extended cytoplasmic region with 3.70Å resolution (PDB ID: 2ZIY) (Shimamura et al., 2008) and turkey β1-adrenergic G protein-coupled receptor with 2.70Å resolution (PDB ID: 2VT4) (Warne et al., 2008) were taken as a best templateshave sequence similarity 21% from C to N-terminal ends.By using several homology servers like Geno3D (Combet et al., 2002), Swiss-Model (Arnold et al., 2006), Phyre (Kelley et al., 2009), I-Tasser (Zhang, 2008) and ESyPred3D (Lambert et al., 2002) servers.Dali server (Holm et al., 2008) has been accessed to compare and validate the template for GPR87 and find that the 2VT4 and 2ZIY are best templates and their sequence identity (21%) for constructing the 3Dstructure of GPR87 by using threading method.The threading method uses specialized schemes to relate protein sequences with library of known structures.Even in cases where there is no clear sequence homology but they have been shown to be able to identify the likely protein fold (Shailza et al., 2011).The use of several templates (two) generally increases the model accuracy as it combines information from multiple template structures.The utility of multiple templates in comparative modeling relies on the accuracy of their multiple structure alignment (Pandian et al., 2011).The coordinates obtained from 3D modeled structures were visualized using the program PyMOL (DeLano et al., 2002).All the different modeled structures of GPR87 protein which were generated by different servers are validated through DOPE score (Discovery Studio2.5) which shows the lowest DOPE scores are the best modeled structures which use for further analysis (Table 1).

Model validation
The analysis of conformational correctness and reliability was carried out by using Ramachandran plot through PROCHECK (Laskowski et al., 1993) for validate the modeled structure of GPR87.The initial model that showed problems in conformation of the different loop region, which was further subjected to loop modeling and validated using side chain and loop refinement protocol of DS2.5 (Accelrys Software, 2003) and also by webserver like ERRAT plot (Colovos et al., 1993).ERRAT plot which gives a measurement of structural errors at each residue of modeled structure.This process was repeated iteratively until most of the amino acids were below 95% cutoff value and residue which lies above 95% cutoff value were subjected to loop modeling in MODELLER.Finally, GPR87 model showing the best PROCHECK and ERRAT plot was then subjected to native protein folding energy evaluation using ProSA server (Wiederstein et al., 2007).DOI:http://dx.doi.org/10.7314/APJCP.2013.14.12.7473Computational 3-D Analysis of Human GPR87 Protein of: Implications for Structure-Based Drug Design

Molecular dynamic simulations
Molecular dynamic (MD) simulations were carried out using the CHARMm (Brooks et al., 1993) force field by module of standard dynamic cascade protocol in DSv2.5.The binding of GPR87 model in presence of explicit solvent was analyzed by molecular dynamics simulation.The protein was solvated using explicit spherical boundary with harmonic restraint solvation protocol of DS2.5 was used.The protein atoms were parameterized using the CHARMm force field (Brooks et al., 1993).All bond lengths involving hydrogen atoms were fixed using a SHAKE algorithm (Ryckaert et al., 1977).Simulations were carried out at 300K with 7000 steps (7ns) in steepest descent and adopted basis Newton-Raphson (ABNR) minimization technique, until the root-mean square deviation (RMSD) was less than 0.001 kcal/Mol/Å R.M.S gradient.The whole system was then equilibrated for 500ps at the target temperature of 300K using 'Berendsen Coupling Bath' and finally the production run was carried out for another 500ps for the canonical ensemble (NVT) using Leapfrog Verlet dynamics Integrator.In all the three stages (Heating, Equilibration and Production) of simulation, time step of 1fs was used and only the backbone of GPR87 protein was kept rigid.

Normal mode analyses on GPR87 protein 3-D structure
Normal mode analysis (NMA) was conducted to analyze the intrinsic motions of the GPR87 model.The elNémo online server (http://igs-server.cnrs-mrs.fr/elnemo/index.html)(Neudecker et al., 1988) is a part of the Elastic Network Model, which provides fast simple tool to compute, visualize and analyze low-frequency normal modes of biological macromolecules.The structural model of hGPR87 in dot pdb (*.pdb) format was submitted to NMA analysis; while another 3Dstructure complex was submitted as reference of the desired conformational change.The key parameters which are used in computation that includes: DQMIN=-100, DQMAX=100, DQSTEP=20 and NRBL="auto".A total of 100 normal modes with the lowest frequencies were requested.The normal mode theory is based on the harmonic approximation of the potential energy function, around a minimum energy conformation.This approximation allows an analytical solution of the motion by diagonalizing Hessian matrix.Hessian matrices are used in large-scale optimization problems within Newtontype methods because they are coefficient of the quadratic term of a local Taylor expansion function i.e.

y=f(x+rx)~f(x)+J(x)rx+yrx T H(x)rx
Where J is the Jacobian matrix, which is a vector (gradient) for scalar-valued functions.The full Hessian matrix can be difficult to compute in practice; in such situations, quasi-Newton algorithms have been developed that use approximations to the Hessian.The most wellknown quasi-Newton algorithm is the BFGS algorithm (Binmore et al., 2001).

Active site prediction
Pocket-Finder or Q-site finder is a molecule-binding site prediction server based on a Ligsite algorithm (Hendlich et al., 1997).It works by scanning a probe radius 1.6Å along all gridlines of grid resolution 0.9Å surrounding the protein.The probe also scans cubic diagonals.Grid points are defined to be part of a site when the probe is within range of protein atoms followed by free space followed by protein atoms.Grid points are only retained if they are defined to be part of a site at least five times (Hendlich et al., 1997).The different active sites in hGPR87 protein also predicted through CASTp server (Dundas et al., 2006) and compare their results with the pocket finder result.The different binding sites and functional residues were identified and stored for further ligand protein interaction analysis.

Docking with substrate ligands
Protein-ligand interaction i.e. docking study was done by LigandFit protocol of DS2.5 (Venkatachalam et al., 2003).The substrate ligands were used i.e. 4-{ 2-carbonitrile (P32), 2-Hydroxymethyl-6-octylsulfanyl -Tetrahydro-Pyran -3,4,5 triol (SOG), palmitic acids (PLM) and retinal (RET) were obtained from the template crystal structure (PDB ID: 2VT4, 2ZIY) (Shimamura et al., 2008;Warne et al., 2008) to dock with GPR87 refined model.The LigandFit docking algorithm combines a shape comparison filter with a Monte Carlo conformational search to generate docked poses consistent with the binding site shape.Among these poses, the most suitable docking mode with a high dock score and a high fitness score from consensus functions was finally selected for analysis.Interaction study with some other protein i.e.EMR2 in H sapiens with chondroitin 4-sulphate is reported (Rani et al., 2011), TLR4 in H sapiens (Das et al., 2012) and GAPDH protein in Leishmania were already been reported (Sahoo et al., 2013).GPCR-membrane models would provide a significant tools for studying receptor behavior and its modulation by small drug-like ligands, a relevant issue for drug discovery (Sadiq et al., 2013).

Sequence alignment and phylogenetic analysis
The translated amino acid sequences of G-protein coupled receptor 87 from Homo sapiens consists of 358 amino acids (Q9BY21) taken from NCBI.Multiple sequence alignment between two groups of P2Y receptors, group I: P2Y 1 , P2Y 2 , P2Y 4 , P2Y 6 and P2Y 11 and group II: P2Y 12 , P2Y 13 and some other pharmacologically related receptors like CysLT 1 and CysLT 2 were obtained from NCBI.Multiple sequence alignment was performed to calculate a distance matrix by using ClustalW and cladogram was constructed by using TreeView algorithm is shown in Figure 1A.By maximum likelihood method clustering of GPR87 of Homo sapiens with Mus musculus showed 89% sequence similarity and with P2Y receptors, CysLT 1 and CysLT 2 receptors producing similarity in between 17-38% Figure 1A.In sequence alignment there are total eight conserved residues were examined in two groups of P2Y receptors and CysLT receptors are C114, R139, P146, W166, Y226, P267, H269 and P311 predicted by ClustalW server.
For the finding of the best template by using BLASTP, three significant hits i.e. 2VT4|A, 2ZIY|A and 3KJ6|A are found, which shows 16-21% sequence identity with GPR87 are shown in Figure 1B.It shows that total 41 amino acids are highly conserved with these three templates sequences.In GPR87 protein, it comprises two different domains that consist of 1-226 residues in the first domain, similarly 227-358 in second domain.The disulfide bond lies in the first domain of GPR87 predicted by Dipro/SCRATCH protein predictor program which also shown in Figure 1B.In disulfide bonds there are 78 amino acids involved in between 114-192 amino acids in length.

Model building and validation
The 3Dstructure of orphan GPR87 protein was built with threading method on the basis of two different templates i.e.PDBID 2VT4 and PDBID 2ZIY by using different homology modeling serves.The comparative analysis of various 3Dstructures which were built by using Geno3D, Swiss-Model, Phyre, I-Tasser and ESyPred3D are summarized in Table 1.The model generated by I-TASSER (model E) is the only model which lays the entire 358 amino acids.For model validation, Ramachandran plot (the Φ/Ψ distribution of the backbone conformation angles for each residue of the 3Dstructure) statistics of model 'C' revealed that the highest 92.1% residues in the core region or 6.8% in allowing region, 0.4% lies in generously allowed region and 0.8% residue in disallowed region is approximate equal to the validation of template structure i.e. 2VT4, but the model structure lay only 286 residues instead of 358.Based on the scores of the Ramachandran plot (validation), the exact length of residues and RMSD value of 3Dstructure in comparison with template main chain model 'E' is the best model are plotted in Figure 2. In model 'E', 2.4% residues lie in disallowed regions which includes Phe3, Leu29, Asn33, Thr74, Lys110, Asp149, Gln183, Asp187, Trp203, Ile241, Glu296 and Lys290 were refined through loop refinement and side chain refinement protocol of DSv2.5 and their RMSD of the main chain with their template is 1.049Å.All the models were validated by DOPE scores (Discrete Optimized Protein Energy) of DSv2.5 and model 'E' has the highest DOPE score is -46755.238281which also proof that it is the best model structure of hGPR87.The DOPE score is the scoring function of protein which can be viewed as a conformational energy and measures the relative stability of a conformation with respect to other conformations in the same protein (Shen et al., 2006).
The secondary structure of GPR87 proteins (α-helices and β-sheets) and their comparisons with two templates (2VT4|A, 2ZIY|A) which contains 7 transmembrane helices (7TM) and 8 small α-helices, 7TM and 10 small α-helices in 2VT4 and 7TM and 12 small α-helices and 2 β-sheets have been observed are represented in Figure 3A.The superimposition of the hGPR87 model (orange) over with its highest closest homolog (2VT4) (cyan) resulted in a root mean square deviation (RMSD) of 1.286Å represented in Figure 3B.Therefore, we influenced that we had reasonable and reliable confirmation of modeled 3Dstructure of GPR87 protein for further analysis.
Different models of GPR87 were further validated by comparison with structural features with their template by using Verify-3D; ProSA and ERRAT servers were used for determining the stereo-chemical parameters of the model structure.ProSA is widely used to check 3D models of protein structures for potential errors.ProSA energy profile of modeled GPR87 in comparison with 2VT4 is presented in Table 1.Verify-3D shows that 53.20% of the residues have an averaged 3D -1D score greater than 0.09 and ERRAT plot gives an overall quality factor of 89.714 for the modeled structure of GPR87 proteins (Table 1).The quality assessment of the model via ProSA revealed that models 'B' and 'E' built by Swissmodel and I-Tasser are matched with the NMR region of the plot with Z-score -1.96 and -1.97 respectively which indicates good quality of the model.Furthermore, we considered I-Tasser modeled structure for further study which consists of 7 transmembrane helices (7TM) and 8 small α-helices prepared by using DS2.5 is shown in Figure 4.In GPR87 there are sequence conflicts in amino acids Ser154 which usually of unknown origin and they may originate from sequencing errors or may be yet uncharacterized natural polymorphism (SNP) reported in UniProtKB (Ota et al., 2004).In Figure 5, schematic representation of 3Dstructure of GPR87 protein is shown using the one- Note that the C-and N-terminus ends are shortened and also the loop between helix 2 and 3 has been shortened.The black letter amino acids indicates conserved Asp-Arg-Tyr (DRY) motif is located at the beginning of the second intracellular loop and at the end of helix 3 where the G-protein binds and thus activation signals are transduced letter residue abbreviation in snake plot.Furthermore, the DRY-motif (an Asp-Arg-Tyr sequence) at the end of helix 3 is highlighted (138-140aa) it's a significant position where the G-protein binds and thus activation signals are transduced.The DRY-motif of GPR87 is located at the beginning of the second intracellular loop which is useful to find novel lead compounds which can be designed on the basis of ligand protein interaction (docking).The DRY-motif patterns are also present in all three different templates i.e. 2VT4|A, 2ZIY|A and 3KJ6|A are shown in the sequence alignment (Figure 2).Similarly, the same hypothesis which have already reported that the C-terminus tail which provide the main driving force for stable interaction between the receptor and β-arrestin and the DRY-motif in these receptors appears to play important roles in β-arrestin binding as well as in maintaining the receptor conformations (Kyeong et al., 2008).

Transmembrane helices (TMs) prediction
There are twelve different transmembranes prediction servers i.e.DAS, TMHMM, HMMTOP, TMpred, TopPred, SOSUI, SPLIT, Predictprotein (PHD), MEMSAT, Wave TM, HMM-TM and SMART servers are used to predict the position and number of transmembrane regions in hGPR87 protein are summarized in Figure 4 and 5.The transmembrane domain of orphan GPCR87 is composed of 6-7 transmembrane helices (TMs) are listed in Table 2.These seven TMs were then modeled by using the known structure of β1-adrenergic GPCR (PDBID 2VT4) as a template.The extracellular loops on the transmembrane domain of GPR87 are relatively long and therefore its 3Dstructure is modeled on the basis of appropriate templates.The comparative analysis of transmembrane helices prediction programs showed that the lowest range and higher range of transmembrane helices in first TM is 42-71 residues, 74-99 in second TM, 111-149 in third TM, 155-180 in fourth TM, 206-232 in fifth TM, 252-278 in sixth TM and 292-322 in seventh TM.In case of TopPred and SOSUI servers both predict only six transmembrane helices at position second and sixth there is missing one transemebrane residues of GPR87 protein respectively.

Protein function assignment of GPR87 protein
Drug design can be applied to arrest the posttranslational modification of GPR87 which will inactivate functional characteristics of that particular protein.It shows six different N-glycosylation sites revealed by using MOTIF SCAN but on PROSCAN server only one site at 129-132 (NITV) is present are shown in Table 3.In GPR87 protein of H sapiens, one N-myristoylation site is present at 224-229 amino acids (GCYIAI).The GPCR87 protein belongs to Zn metallo-β-lactamase family having an E-value of 0.00015.It is a member of rhodopsin-like GPCR superfamily which varies in between 44-322 residues and rhodopsin-like GPCR superfamily (7tm_1) which lays involvement of 59-314 amino acids (Table 3).The closest structural domain lies in lactamase_B domain i.e. metallo-β-lactamase superfamily which shows 10-158 residues by MOTIFSCAN and 5-148 through the Intrepro tool of SWISSMODEL.The family metallo-β-lactamases which degrade all β-lactam drugs including carbapenems and requires zinc as a cofactor (Varsha, 2008).These enzymes are not inhibited by routine β-lactamase inhibitors, so these compounds are being evaluated to inhibit these types of enzymes.

Molecular dynamic simulations
Based on intrinsic dynamic, structural stability and improved relaxation of the modeled GPR87 protein, energy minimization by steepest descent methods of standard dynamic cascade protocol of DSv2.5, which is an important step for the convergence of free MD simulations, were calculated.The effect of solvent on the stability of binding of the GPR87 was studied by MD simulation of a fully hydrated model is shown in Figure 6.The radius of the sphere is 40Å and salt concentration is 0.145 is fixed.The total energy and temperature were found to remain steady with little fluctuation during the production stage of 500ps.During the production stage of 500ps following the initial heating and equilibration phases, total energy and simulation temperature were found to remain steady with slight fluctuation.The snapshots of steepest descent dynamics trajectory at 0 to 8000 of the production runs are shown in Figure 7.The result indicates that the variation in total energy is -15782.28127KJ/Mol to -19888.33249KJ/Mol calculated by the steepest descent method at 6500ps dynamics trajectory was retrieved and further minimized to full convergence.Van der Waals energy did not show much variation both by steepest descent method (-2344.01480KJ/Mol to -2312.75325KJ/Mol).The

Normal mode analyses on GPR87 3Dstructure
Normal Mode Analyses (NMA) program was conducted to analyze the intrinsic motions (vibrations and thermal properties) of the GPR87 protein modeled structure.The structural model of GPR87 in apo-form was submitted for NMA analysis; while the GPR87/2VT4 model complex was submitted as a reference of the desired conformational change.Essential features of the top ten low-frequency normal modes, including its frequency, amplitude, collectivity of atom movements and the overlap with observed conformational changes are illustrated in Table 4.The two lowest-frequency normal modes (mode 7 and 8 in Table 4) are illustrated in Figure 8 one feature showing bending of the N-terminal domain towards the transmembrane domain, the other feature showing rotation of the N-terminal domain on the top of the transmembrane domain.In order to confirm the outcomes of NMA analysis, we also conducted the principal component analysis (PCA) on GPR87 in apo-form to detect its major conformational motions during MD simulation.Each snapshot of the GPR87 structure was fitted with a reference configuration to remove translational and rotational motions.

Protein-ligand interaction study
After knowing the active sites in GPR87 protein, we analyze the specific substrate ligands of the template were docked with both the 3D model and crystal structures optimally to estimate the site of reaction on the ligand molecule and their binding affinity with the target protein.
There are sixteen different binding sites were detected in GPR87 protein and also in 2VT4 but 2ZIY have fourteen sites detected by using find sites from receptor cavities protocol of DSv2.5.
The highest ligandfitness score has been detected as 61.414 with P32 and low as 56.952 with palmitic acids (PLM) having binding affinity for GPR87 protein are listed in Table 5.The ten different conformations (poses) have been observed during docking study in DS (Accelrys).The groove contains several positively charged side chains that interact directly with the negatively charged groups on P32, SOS, PLM and retinal are shown in Figure 10 A-D.The GPR87 belongs to the members of the opsin family, which have been considered to be typical members of the rhodopsin superfamily which share several motifs mainly 7TMs, GCPRs of the rhodopsin superfamily.Mostly the opsin that binds with chromophore, such as 11-cis-retinal but in case of retinal the binding affinity is 60.675 and there is an absence of H-bonding.The most significant function of opsins other than the photoisomerases is split into two steps: light absorption and G-protein activation.The interactions between all the four substrate molecules and GPR87 protein, the frequently involve amino acids in the formation of H-bond are Arg115, Ser118, Tyr122, Pro198, Gly200, Val201, Lys202, Tyr268 and Lys296 and their atom involves OG, HH, O, HN, HN, O, HZ1, HZ3, OH and HZ3 revealed ion-pair interactions with substrate ligands atoms are O18, O20, O2, O1, N5, N17 and O16 respectively.For the same the crystal structure interaction shows the 69.723 dock score with P32 with PDBID 2ZIY and 45.726 with PDBID: 2VT4 and interacts frequently with Lys305, Pro312, Trp332, Val333 and Arg284, Arg345 and Arg350.Our docking results suggest that, GPR87 can bind close to the active site with similar binding energy with template crystal structures so, we concluded that these substrate ligands and their various analogues can be the best suitable lead compounds for the development of a novel class of selective drugs for cancer in order to

Discussion
In conclusion, this is the first comprehensive study highlighting the structural features of a potential drug target i.e. orphan GPR87 protein.The structural aspects of the active site and the different predicted motifs help to understand the functional details of the protein.On the basis of modeled structure, predicted binding sites and functions, high throughput screening of various compounds may be carried out to find appropriate leads against squamous cell carcinoma, which can aid in the design and development of novel drug molecules against this drug target GPR87.In the absence of crystal structure, our study provides an early insight into the structure of the major drug target GPR87 and facilitating the inhibitor design.As the 3-D structure of GPR87 protein is known from this study, novel lead compounds can be designed on the basis of ligand protein interaction (docking) of available anti-cancer drugs and designing various analogues of the presently available drugs or defining a novel candidate on the basis of interaction with different binding sites of GPR87 protein in H sapiens.It would be helpful in the future if computational methods could provide guidelines as to what chemical manipulations of ligands will be useful to control protein function.In conclusion, substrate-binding prediction by docking studies and molecular dynamics may contribute to the identification of metabolic routes to drug development and thus aid in the designing of new anti-cancer drugs.

Figure
Figure 1.A) Phylogenetic tree of hGPR87 of H sapiens with P2Y receptors.Group I: P2Y1, P2Y2, P2Y4, P2Y6 and P2Y11; and group II: P2Y12, P2Y13 and Some Other Pharmacologically Related Receptors Like CysLT1 and CysLT2 were Obtained from http://www.phylogeny.fr/with Maximum Likelihood Method.The tree is formed with 0.99 maximum sequence difference values; B) Multiple Sequence Alignment of 3 Different Templates their PDB IDs: 3KJ6|A 2R4R|A and 2VT4|A Taken from BLASTP Which Shows 16-21% Sequence Similarity with GPR87 Protein and Find 41 Conserved Residues.The disulphide bridges are comprised in 114-192 residues represented in blue color

Figure 2 .Figure 4 .Figure 5 .
Figure 2. Ramachandran Plot of Modeled 3Dstructure (model E) of GPR87 Protein in H sapiens Calculated by the PROCHECK Program

a
Only the 10 normal modes with lowest frequencies are displayed here; b The level of collectivity indicates the percentage of residues that are involved in a certain normal mode; c The level of overlap measures the similarity between a desired conformational change and that of a certain normal mode

Figure 6 .
Figure 6.Overview of GPR87 Receptor Solvated in the Hydrated Lipids Bilayer.The backbone of the receptor is represented in red helix, water is in red/white and the internal water molecules are displayed as spheres

Figure 7 .
Figure 7. Comparisons and calculate energy versus time plot by using Steepest descent of energy minimization protocol of DSv2.5 software.The X-axis is Time (ps) and Y-axis Total energy in KJ/Mol

Figure 8 .Figure 9 .
Figure 8.The Two Different Normal Modes (i.e., modes 7 and 8) of GPR87 Protein with the Lowest Frequencies.Left: bending of the N-terminal domain towards the transmembrane domain.Right: rotation of the N-terminal domain on the top of the transmembrane domain

Figure 10 .
Figure 10.Interaction of GPR87 Protein with Substrate Ligands of Template Crystal Structure (PDB ID: 2VT4 and 2ZIY).In the above figure, A) represents P32; B) is SOG; C) is palmitic acids; and D) is retinal.It shows ten different conformations (poses) have been observed during docking study.The green color is the substrate ligands, protein is in atomic color, green dotted lines represented as H-bonding and GRASP model of GPR87 protein shown in electrostatic potential energy red is positive and blue is negative and white are neutral.All figures generated by using DSv2.5

Table 1 . The Comparative Analysis of Various 3-D Models of GPR87 Protein Which were Built by Different Servers Like Geno3D, SwissModel, Phyre, I-Tasser and ESyPred3D
*Predicted number of TMs; **Numbers in bracket are the lengths of TMs