Crucial Positively Charged Residues for Ligand Activation of the GPR35 Receptor*

Background: The structure and function of GPR35 are not understood. Results: Using a GPR35 activated state molecular model, we identified crucial amino acid residues required for ligand activation using β-arrestin trafficking, ERK1/2 activation, and calcium imaging. Conclusion: Arginines in TMH3-4-5-6 affected agonist signaling. Significance: Identification of residues for GPR35 agonist signaling is critical for the design of ligands with improved potency. GPR35 is a G protein-coupled receptor expressed in the immune, gastrointestinal, and nervous systems in gastric carcinomas and is implicated in heart failure and pain perception. We investigated residues in GPR35 responsible for ligand activation and the receptor structure in the active state. GPR35 contains numerous positively charged amino acids that face into the binding pocket that cluster in two distinct receptor regions, TMH3-4-5-6 and TMH1-2-7. Computer modeling implicated TMH3-4-5-6 for activation by the GPR35 agonists zaprinast and pamoic acid. Mutation results for the TMH1-2-7 region of GPR35 showed no change in ligand efficacies at the K1.32A, R2.65A, R7.33A, and K7.40A mutants. However, mutation of arginine residues in the TMH3-4-5-6 region (R4.60, R6.58, R3.36, R(164), and R(167) in the EC2 loop) had effects on signaling for one or both agonists tested. R4.60A resulted in a total ablation of agonist-induced activation in both the β-arrestin trafficking and ERK1/2 activation assays. R6.58A increased the potency of zaprinast 30-fold in the pERK assay. The R(167)A mutant decreased the potency of pamoic acid in the β-arrestin trafficking assay. The R(164)A and R(164)L mutants decreased potencies of both agonists. Similar trends for R6.58A and R(167)A were observed in calcium responses. Computer modeling showed that the R6.58A mutant has additional interactions with zaprinast. R3.36A did not express on the cell surface but was trapped in the cytoplasm. The lack of surface expression of R3.36A was rescued by a GPR35 antagonist, CID2745687. These results clearly show that R4.60, R(164), R(167), and R6.58 play crucial roles in the agonist initiated activation of GPR35.

Immunocytochemistry and Microscopy-The expression of wild-type and mutant GPR35 receptors was examined by immunostaining of cells with HA monoclonal antibodies. To achieve surface receptor staining, cells grown on coverslips were incubated over ice for 40 min with a 1:100 dilution of mouse monoclonal anti-HA antibody in blocking buffer (3% BSA in PBS). This was followed by appropriate washes and 40-min incubation with a 1:1500 dilution of Alexa Fluor 568 goat anti-mouse secondary antibody. Cells were then fixed with 4% paraformaldehyde for 20 min at room temperature followed by three washes with PBS. The total receptor immunostaining was done by fixing then permeabilizing cells with 0.2% Triton X-100 before the primary HA antibodies (1:500) application. Following the staining, glass coverslips were mounted on slides and were imaged on a (Nikon E800) fluorescence microscope using a 40ϫ oil objective and 488-nm excitation for GFP. Coimmunostaining of R3.36A and calreticulin (18) was done by co-applying mouse anti-HA antibodies (1:500) and rabbit anticalreticulin antibodies (Cell Signaling; 1:200) followed by Alexa Fluor 568 goat anti-mouse secondary antibody (1:1500) and Alexa Fluor 488 goat anti-rabbit secondary antibody (1:1000).
␤-Arrestin Assay for Determining Receptor Responsiveness-The ␤-arrestin assay was described previously (5). Briefly, stable cell lines were plated onto poly-D-lysine-coated coverslips in 24-well plates. Cells were maintained at 37°C in 5% CO 2 until ready for experiments (80 -85% confluent) and washed once with Hanks' balanced salt solution before drug application, and experiments were performed in Hanks' balanced salt solution. Agoniststimulated redistribution of ␤arr2-GFP was assessed following drug treatment for 40 min at room temperature. Cells were fixed and imaged as described above. Quantification of ␤arr2-GFP aggregates was described previously (5). Concentration-effect curves for agonist-mediated receptor activation were analyzed by nonlinear regression techniques using GraphPad Prism 5.0 software (GraphPad), and data were fitted to sigmoidal dose-response curves to obtain EC 50 or IC 50 values. Statistical analysis was performed using one-way analysis of variance followed by Dunnett's post-test or two-tailed unpaired Student's t test. p values of Ͻ0.05 are considered significant.
In-Cell Western Assay for ERK Activity-The assay was described previously (5). Briefly, cells were grown to confluence in 96-well plates and serum-starved overnight prior to assay. Following the 20-min agonist treatment at room temperature, the medium was removed, and 4% paraformaldehyde in PBS was added to fix cells for 1 h at room temperature. Cells were then permeabilized by 0.1% Triton X-100 in PBS for five washes, 10 min per wash. LI-COR blocking buffer was added, and samples were shaken on a rotator for 1 h. Primary antibodies against phosphor-ERK1/2 (Cell Signaling Technology; 1:100) were applied for 2 h, and then secondary antibodies goat anti-rabbit 800CW (1:800) were applied overnight in a cold room. Sapphire700 (LI-COR, 1:1000) and DRAQ5 (Biostatus Limited, 1:2000) were added together with the secondary antibodies for normalization. The plate was dried and then scanned using a LI-COR Odyssey Infrared Imager set at 169 M resolution, 3 focus offset, and 4.5-6 intensity. Data were analyzed using Excel and Prism 5.0 software.
Calcium Imaging-Intracellular Ca 2ϩ measurements were performed as described previously (19). Briefly, GPR35 and mutant stable U2OS cells were loaded with fura-2/AM and mounted in an open bath chamber (RP-40LP; Warner Instruments, Hamden, CT) on the stage of an inverted microscope Nikon Eclipse TiE (Optical Apparatus Co., Ardmore, PA). Fura-2/AM fluorescence (emission ϭ 510 nm) following alternate excitation at 340 and 380 nm was acquired at a frequency of 0.25 Hz. Images were acquired and analyzed using Nikon NIS-Elements AR 3.1 software (Optical Apparatus Co.). The ratio of the fluorescence signals (340/380 nm) was converted to Ca 2ϩ concentrations.
Quantification of Receptor Surface Expression-We used average fluorescent signal intensity analysis of surface receptor staining on the whole cell surface to measure total expression levels (20). Within the confines of the demarcation of the cell area in each image, we measured the average fluorescent signal to sample receptor expression level by using ImageJ software (1.47v; National Institutes of Health). Sampling the protein expression level of receptors in the easily definable cell area in the images provides an estimate of the total receptor expression that can be used to compare WT and mutant expression levels.
Conformational Memories (CM) Technique for Calculating TMH Conformation-For TMH 1, 2, 4, 5, 6, and 7, a series of possible, sequence-dictated conformations that deviate from the ␤ 2 -AR template were calculated using the CM Monte Carlo/simulated annealing technique (33,34). The CM technique uses multiple Monte Carlo/simulated annealing random walks employing the CHARMM all atom force field (35) in a distancedependent dielectric at 310 K. With this technique, nearly all dihedrals and bond angles are defined as variable, with various restrictions. Bond angle variation ranges were defined as either Ϯ8°for default bond angles or Ϯ15°for side chain angles involving polar hydrogens (i.e., C-O-H: Ser, Thr, Tyr) or the flexible C-S-C bond angle in methionine. The backbone dihedrals ⌽, ⌿, and had standard dihedral variation ranges of Ϯ10°, Ϯ10°, and Ϯ20°, respectively, and the side chain dihedrals had full range of motion. If the sequence contained a known helix deforming residue such as proline, serine (36), threonine (36), or glycine, a larger variation range of up to Ϯ50°on the ⌽ and ⌿ was implemented to allow for increased flexibility. An ideal helix (⌽ ϭ Ϫ62.9°and ⌿ ϭ Ϫ41.6°) with the GPR35 sequence for each TMH was built, and a minimum of 105 conformers were generated for each GPR35 helix independently.
TMH1 in GPR35 has three nonconsecutive glycines located at positions 1.34, 1.39, and 1.46. An ideal helix with the GPR35 sequence was built, and three independent calculations were performed with increased flexibility of the backbone (⌽/⌿ of Ϯ 50°) at the residues mentioned above.
TMH2 in GPR35 has a proline that is shifted by one residue intracellularly as compared with the ␤2-AR sequence: P2.59 in the ␤2AR and P2.58 in GPR35. An ideal helix with the GPR35 sequence was built, and a flexible region (⌽/⌿ Ϯ 50°) was implemented into the P2.58 to L2.54 region to evaluate the flexibility at this proline position.
TMH4 in GPR35 lacks the conserved proline located at the extracellular end of the ␤2-AR TMH4, but it contains a potential helix deforming motif, GS (G4.55-S4.56). Two independent conformational searches were performed from an ideal helix starting structure to determine the impact of this motif. The first explored only the increased flexibility of the backbone at G4.55 (⌽/⌿ of Ϯ 50°), and the second explored both backbone flexibility and 1 of S4.56 restrained to g(Ϫ) (36).
TMH5 in GPR35 has two proline residues; the conserved proline at position 5.50 and an additional proline at 5.43. The flexibly of the proline kink of the conserved proline (G5.45 to P5.50) was restricted to empirically based backbone dihedral ranges, and the proline kink of nonconserved proline (S5.39 to P5.42) was given increased backbone flexibility (⌽/⌿ of Ϯ 50°).
TMH6 in GPR35 lacks the highly conserved CWXP motif but has the conservative substitution, CFXP. The importance of exploring the possible conformational space of TMH6 lies in the documented changes in TMH6 conformation that occur when a class A GPCR is activated (37). In ␤2-AR and rhodopsin, the most dramatic structural change associated with activation is the conformational change in TMH6 about the highly conserved CWXP hinge region (38 -45). Two independent conformational searches were performed with an identical, ideal helix starting structure to be representative of either an R or R* state. The flexibility of the backbone of the conserved proline kink region (V6.46 to P6.50) was restricted to the empirically based ranges of the antagonist bound GPCR crystal structures for the R state and allowed full flexibility (⌽/⌿ of Ϯ 50°) for the R* state. G6.54 was also considered flexible ((⌽/⌿ of Ϯ 50°) in both calculations.
TMH7 in GPR35 lacks the conserved P7.50 in the NPXXY motif and has the nonconservative substitution of DAXXY. To determine whether the GPR35 TMH7 sequence may suggest a flexible region within the same area as the conserved proline kink region, an increased flexibility (⌽/⌿ Ϯ 50°) was introduced in the region of C7.46 to A7.50. The TMH7 conformation chosen for bundle creation was one that allowed D7.43 to hydrogen-bond to Y3.32. This interaction is found in the CXCR4 crystal structure (30).
Construction of the GPR35 Inactive State TMH Bundle-The conformers generated by the CM calculations were used to construct the TMH region of the receptor. The receptor was assembled by aligning the intracellular regions of TMHs 1, 2, 4, and 5-7 with the ␤ 2 -AR template. In the R model, R3.50 and T6.30 form a hydrogen bond that mimics the R3.50/D or E6.30 ionic lock common in class A GPCRs. The bundle was then pulled apart 3 Å and energy-minimized for 1500 iterations utilizing the OPLS2005 all atom force field in Macromodel 9.1 (Schrodinger Inc., Portland, OR) employing a distance-dependent dielectric with extended cutoffs (nonbonded, 8.0 Å; electrostatic, 20.0 Å; and hydrogen bonding, 4.0 Å) and the Polak-Ribiere conjugate minimization scheme to pack the TMHs and relieve side chain clashes. A harmonic constraint was placed on the backbone (, , and ) of 500 kj/mol to maintain helicity.
For the inactive (R) bundle, several highly conserved residue side chain dihedrals were restrained to crystal structure ranges.
Modeling of the Intra/Extracellular Loops and Receptor Termini-Modeler 8.2 (46,47) was used to generate the N and C termini (Met-1 to Pro-16 and Pro-293 to Ala-309, respectively) and the intra-and extracellular loops (IC1: Arg-48 to Glu-54, IC2: Pro-121 to Ser-129, IC3: Ala-199 to Thr-212, EC1: Asp-82 to Thr-86, EC2: Leu-153 to His-168, EC3: Gly-244 to Ala-247). A restrictive distance was implemented between Cys-(162), located on the EC2 loop, and Cys-3.25(89) to emulate the class A highly conserved disulfide bridge between residues at these positions. A restricted distance was also implemented between C(8) in the N terminus and C7.25(248) at the EC end of TMH7 to emulate the disulfide bridge seen in the CXCR4 crystal structure (30). An EC-2 loop conformation was selected that had the second residue past the disulfide bridge (R(164)) pointing into the bundle as this is seen in the ␤ 2 -AR (25)(26)(27) and other x-ray crystal structures. With the TMHs frozen, the final loops were energy-minimized to a gradient of 0.01 kcal/mol in a high dielectric of 80 using the same force field, cutoffs, and conjugate method described earlier.
Modeling of the Human GPR35 R* Receptor-To create an activated GPR35 model (R*), The EC2 and IC3 loops were removed, and a second TMH6 conformer was chosen that was straighter than TMH6 in the R model. This conformer was chosen so that no hydrogen bond was broken. With the loops and termini frozen, the TMH region was energy-minimized for 500 iterations in a distance-dependent dielectric utilizing the same force field, cutoffs, and conjugate method described earlier, with a harmonic constraint was placed on the backbone (, , and ) of 1000 kj/mol to maintain helicity. The EC2 and IC3 loops were then modeled as described above, including the distance restraint in the EC2 loop to mimic the conserved disulfide bridge. The final loops were then energy-minimized, as described above, to a gradient of 0.01 kcal/mol in a high dielectric of 80.
Modeling of the Human R6.58A Mutant GPR35 R* Receptor-To create the R6.58A R* receptor model, residue 6.58 in the WT R* bundle was mutated to Ala, and the bundle was reminimized using the same protocol as described above.
Hartree-Fock ab initio calculations at the 6 -31ϩG* level as encoded by Spartan '08 (Wavefunction, Inc., Irvine, CA) were performed for each compound. For each comformational search, local energy minima were identified by the rotation of subject dihedrals 360°in 60°increments (6-fold search), followed by a Hartree-Fock 6 -31ϩG* energy minimization of each rotamer generated. For this reason, the energy minimum of each tautomeric form was assessed at the 6 -31ϩG* level. Pamoic acid conformational analyses were performed in vacuo and using a SM8 water solvation model as implemented in Spartan '08 (Wavefunction, Inc.).
Docking of Ligands-The global minimum energy conformer of pamoic acid (dianion) and the global minimum energy conformer of the 2,6-dihydro tautomeric state of zaprinast were used for docking studies (see numbering system in Fig. 2). Because of recent experimental data (49), Y3.32(96) was used as the primary interaction site for zaprinast, and two arginine residues located near the binding pocket (R3.36(100) and R4.60(151)) were used as primary interactions for the pamoic acid (dianion). Initially, each compound was manually docked in the binding site of GPR35 R* (or R6.58A R*) model.
Ligand/Receptor Minimization-The ligand-receptor complexes were minimized using an OPLS2005 all atom force in Macromodel 9.9 (Schrodinger Inc.). An extended cutoff (nonbonded, 8.0 Å; electrostatic, 20.0 Å; hydrogen bonding, 4.0 Å) was used in each stage of the calculation. The minimization was broken into two stages: 1) the minimization of the ligand and the TM regions. The first step consisted of 1250 steps of conjugate gradient minimization using a distance-dependent dielectric. Harmonic constraints of 1000 kcal/mol were placed on select dihedrals of the ligand to maintain conformational similarity to the low energy conformer calculated at the Hartree-Fock 6 -31ϩG* level of theory as the receptor. The ligand restraints were released for the last 250 steps to allow the ligand to adapt to the optimized binding pocket. The EC and IC loop regions were frozen during this stage. 2) The extracellular and intracellular loops were minimized for 500 steps in a generalized born/surface area continuum solvation model for water, whereas the TM regions and the ligand atoms were held frozen.
Glide Docking Studies-The automatic docking program, Glide (Schrodinger Inc.) was then used to first score the manual docks and then explore other possible binding conformations or receptor site interactions with flexible docking. Glide was used to generate a grid based on the centroid of the ligand in the binding site (from the manual docks). Any hydrophobic region defined in the grid generation that contacted the ligand was selected as important to the flexible docking procedure. The box for flexible docking was defined large enough to encompass Y3.32(96) as the primary interaction site for zaprinast and R3.36(100) and R4.60(151) as primary interactions for pamoic acid. Standard precision (SP) was selected, and flexible docking was invoked unless scoring a manual dock (50).

GPR35 Receptor Model Development
A homology model of the GPR35 inactive state was developed using the ␤ 2 -AR crystal structure (Protein Data Bank code 2RH1 (25)) as the template. This initial model was refined based on sequence dictated differences between GPR35 and the ␤ 2 -AR crystal structure. The key differences between these two TMH bundles are described below.

TMH1
Despite the presence of three glycines (G1.34, G1.30, and G1.46) in TMH1, conformational memories calculations revealed that TMH1 in GPR35 is similar to TMH1 in ␤ 2 -AR, which leans away from the TMH bundle on the extracellular side.

TMH2
Like the CXCR4 receptor (30), GPR35 contains the XX(C/S/ T)LP motif in TMH2, with proline position at 2.58. This motif creates a tightly wound proline kink where the C/S/T residues in the proline kink region are in g(Ϫ), supporting the tightly wound kink structure. The ␤ 2 -AR has a relaxed broad turn in the 2.59 proline kink region. The different turn ratio between the ␤ 2 -AR and the GPR35 TMH2 proline kink region supports TMH2 fitting into the general seven-TMH GPCR structure, but with a different set of residues facing into the binding site. For example, residue V2.60 faces into the binding pocket of GPR35, whereas residue 2.60 in the ␤ 2 -AR structure (25) is in the TMH2-3 interface and not available from within the binding pocket. Another consequence of the Pro at 2.58 in GPR35 is that the extracellular portion of TMH2 leans toward TMH3. It is important to note that had the GPR35 TMH2 retained the ␤ 2 -AR TMH2 conformation, the GPR35 EC1 loop would not have been able to join the EC ends of TMH2 and TMH3, because the GPR35 EC1 loop is only two amino acids long (SD). In contrast, the EC1 loop of the ␤ 2 -AR is six amino acids in length (KMWTFG).

TMH4
GPR35 lacks the proline at 4.60 found in the ␤ 2 -AR. CM results suggest that the lack of this proline at 4.60 causes a local change at the extracellular end of TMH4, which leans toward TMH5 compared with the ␤2-AR crystal structure (Fig. 3).

TMH5
TMH5 in GPR35 has a proline at 5.50 and at 5.43. Because the ␤ 2 -AR TMH5 has the conserved P5.50 as does GPR35 TMH5, the backbone torsions for that Pro kink were preserved and only the region around P5.43 in GPR35 TMH5 was explored via CM. These calculations suggested that P5.43 leans the extracellular end of TMH5 away from the binding crevice in GPR35. This change is accommodated in GPR35 because the GPR35 EC2 loop is longer than ␤ 2 -AR (GPR35 11 amino acids versus 5 amino acids in ␤ 2 -AR).

TMH6
The CM output gives helices that bend inward toward the binding crevice and also helices that bend toward TMH7. From the helices generated by CM, we picked two helices, one that would maintain the R3.50/T6.30 hydrogen bond TMH6 R (C-␣ to C-␣ distance ϭ 6.15 Å) and one for which this interaction would be broken. TMH6 R had the following bend, wobble, and face shift angles: 32.8°, Ϫ62.1°, and 129.5°. TMH6 R* had the following bend, wobble, and face shift angles: 32.3°, Ϫ71.2°, and 83.3°. The differences in these two conformers are in their wobble and face shift angles. These changes result in the intracellular end of TMH6 in R* pivoting away from the IC end of TMH3.

TMH7
GPR35 TMH7 lacks the highly conserved class A GPCR NPXXY motif, having instead DAICY. The lack of a proline residue in TMH7 resulted in CM output that did not have the usual proline induced bend and deformation seen in crystal Activation of GPR35 FEBRUARY 7, 2014 • VOLUME 289 • NUMBER 6

JOURNAL OF BIOLOGICAL CHEMISTRY 3629
structures of other class A GPCRs. However, GPR35 has a charged residue, K7.40, that faces toward TMH2 in the binding pocket that should draw water into the TMH1-2-7 region. We have shown previously for GPR55 that such charged residues (D7.43 and D7.49 for GPR55) in TMH7 become hydrated, and this water further perturbs the backbone structure on TMH7 by hydrogen bonding to backbone carbonyls. The net result is a region of 3-10 helix from 7.43 to 7.50. Fig. 3 illustrates a comparison of the activated ␤ 2 -AR crystal structure (37) (Fig. 3A) and the GPR35 active state model (Fig.  3B). Here the view is from the extracellular side of each receptor. One hallmark of GPCR activation is the breaking of the ionic lock between R3.50 and E/D6.30 (T6.30 in GPR35), which allows TMH6 to change conformation. The result is an intracellular opening of the receptor, exposing residues that can interact with the C terminus of G␣ (51). The active form of the GPR35 receptor model differs from the inactive form mainly in changes at the intracellular side of the receptor. The intracellular end of TMH6 has changed conformation and moved its IC end, breaking the hydrogen bond between R3.50 and T6.30. The distance between the C␣ carbons of R3.50 and T6.30 is 11.8 Å in the GPR35 R* model, indicating loss of interaction between the IC ends of TMH3 and TMH6. This distance is comparable to the opening formed in the rhodopsin meta II (activated state) crystal structure (R3.50 -E6.30 C␣ distance ϭ 14.7 Å) (52) in which the C-terminal fragment of G-␣ of transducin is inserted. The size of this same opening in the nanobody stabilized ␤ 2 -AR crystal structure (Fig. 3A), is larger (R3.50 and E6.30 C␣ distance ϭ 17.2 Å) (37), because of the physical insertion of the nanobody higher than the C-terminal fragment of G-␣ of transducin, which causes the IC opening to enlarge.

GPR35 R6.58A Model Activated Form (R*)
To create the R6.58A mutant, residue 6.58 in the WT R* bundle was mutated to Ala, and the bundle was reminimized. This mutation did not change the overall conformation of TMH6. The mutant TMH6 had the following bend, wobble, and face shift angles: 32.7°, Ϫ74.4°, and 83.1°. The distance between the C␣ carbons of R3.50 and T6.30 in the R6.58A mutant R* model was 11.9 Å, indicating loss of interaction between the IC ends of TMH3 and 6.

Docking Studies: Pamoic Acid
Pamoic Acid/Wild-type Human GPR35 R*-We used our GPR35 activated state model to probe the nature of the interactions of pamoic acid dianion within the GPR35-binding pocket. Pamoic acid dianion is both charged and highly aromatic. Docking studies suggested that the pamoic acid-binding pocket is in the TMH3-4-5-6 region of GPR35. Fig. 4A illustrates the final docked pose of pamoic acid in the GPR35 R* model (Glide score SP 5.0, Ϫ8.81 kcal/mol). In this docked position, pamoic acid is able to form several highly favorable polar interactions with Arg residues, including several salt bridges and cationinteractions (residues colored magenta). The ligand also forms hydrogen bonds (residues colored orange) and predominantly van der Waals interactions with other residues (residues colored yellow). The geometries of these interactions are discussed below. Pamoic acid establishes a salt bridge, as well as a hydrogen bond with the EC2 loop residue, R(164).  Pamoic Acid/Human R6.58A Mutant GPR35 R*-Despite the fact that pamoic acid is predicted to dock in the arginine rich TMH3-4-5-6 region of GPR35, not all of these arginines appear to be necessary for ligand binding. A case in point is R6.58. Our modeling studies (see above) suggest that this resi-due does not interact with pamoic acid. We created an R6.58A GPR35 R* model to examine possible changes occurring upon mutation. Fig. 4B illustrates the final docked pose of pamoic acid in the GPR35 R* model (Glide score SP 5.0, Ϫ8.95 kcal/ mol). In this docked position, pamoic acid is able to form several highly favorable polar interactions with Arg residues, including several salt bridges and cationinteractions (residues colored magenta). The ligand also forms hydrogen bonds (residues colored orange) and predominantly van der Waals interactions with other residues (residues colored yellow). The interactions that pamoic acid has in the R6.
Docking Studies: Zaprinast Zaprinast-GPR35 R* Complex-The EC 50 of zaprinast is considerably higher than that of pamoic acid. Docking studies suggested that like pamoic acid, the binding pocket for zaprinast is in the TMH3-4-5-6 region of GPR35. Fig. 5A illustrates the final docked pose of zaprinast in the GPR35 R* model (Glide score SP 5.0, Ϫ4.14 kcal/mol). Unlike pamoic acid, zaprinast is not charged and therefore cannot form salt bridges like pamoic acid. In its docked position, zaprinast is able to form hydrogen bonds (residues colored orange) and predominantly van der Waals interactions with other residues (residues colored yel- Residues that form salt bridge and cation-interactions are colored magenta, whereas residues that form hydrogen bonds are colored orange, and residues that contribute mainly van der Waals interactions are colored yellow. Residues that are not part of the ligand-binding pocket are colored gray. low). The CϭO oxygen of the zaprinast fused ring system can form a hydrogen bond with R(164). The N-O distance is 2.64 Å, and the N-H-O angle is 168°. Q3.29 is too far to hydrogenbond with zaprinast, but this residue has significant van der Waals interactions with zaprinast. R4.60 (colored cyan) forms a boundary of the zaprinast-binding pocket but does not have a significant direct interaction with zaprinast. Several other amino acids have predominantly van der Waals interactions with zaprinast, including L3.33, L6.51, and L6.55. Y3.32 has a tilted T aromatic stack with the heterocyclic fused ring of zaprinast. The distance is 4.14 and 5.65 Å for the five-and six-member rings, respectively, and the angle is 49°. R6.58 and R(167) (colored gray) are above the zaprinast-binding pocket and have no direct interaction with zaprinast.
Zaprinast-GPR35 R6.58A R* Complex Fig. 5B illustrates the final docked pose of zaprinast in the R6.58A GPR35 R* model (Glide score SP 5.0, Ϫ5.38 kcal/mol). Our modeling studies suggest that R6.58 in WT GPR35 does not interact with zaprinast but produces crowding in the binding region. Removal of R6.58 via the R6.58A mutation is actually favorable for zaprinast because it acquires additional interactions compared with WT. In its docked position, zaprinast is able to form hydrogen bonds (residues colored orange) as well as predominantly van der Waals interactions with other residues (residues colored yellow). Zaprinast can form a hydrogen bond between Q3.29 and the CϭO oxygen of its fused ring system. The N-O distance is 2.84 Å, and N-H-O angle is 143.5°. Additional contributors to the zaprinast interaction energy come from residues that form predominantly van der Waals interactions with the ligand. These include L6.55, R(164), L3.33, and L6.51. Zaprinast also has two aromatic stacking interactions that add to its receptor interaction, one with F5.47 and one with Y3.32. F5.47 forms a tilted T aromatic stack with the propoxy-benzyl ring of zaprinast. The ring centroid to centroid distance is 6.22 Å, and the angle of the planes of the two rings is 50°. The fused ring system also stacks with F5.47. The ring centroid to centroid distances of the five-and sixmember rings are 5.61 and 5.90 Å, respectively, and the plane to plane angle is 59°. Y3.32 forms a tilted T interaction with the fused ring system of zaprinast. The centroid to centroid distance of the two rings is 4.77 Å, and the plane to plane angle is 72°. R4.60 (colored cyan) forms a boundary of the zaprinastbinding pocket but does not have a significant interaction with zaprinast. R(167) (colored gray) is above the zaprinast-binding pocket and has no direct interaction with zaprinast.
HA immunostaining of the R3.36A and R7.33A (data not shown) mutant U2OS cell lines showed no cell surface expression (Fig. 8A, left panels), but instead showed strong cytosolic expression in permeabilized cells (Fig. 8A, second panels from  left). Immunostaining of HEK293 cells transiently transfected with R3.36A plasmid also showed no surface expression but strong cytosolic expression (Fig. 8A, right two panels). No Residues that form hydrogen bonds are colored orange, and residues that contribute mainly van der Waals interactions are colored yellow. Residues that are not part of the ligand-binding pocket are colored gray, whereas residues that may form a boundary of the binding pocket are colored cyan.
␤-arr2-GFP response was observed in R3.36A mutant cells upon pamoic acid or zaprinast treatment (Fig. 8B). To localize the cytosolic HA-tagged R3.36A mutant, U2OS cells expressing R3.36A alone were permeabilized and immunostained with HA and the ER-resident protein calrecticulin (18). High magnification of confocal images reveals a perinuclear labeling that is characteristic of ER staining (Fig. 8C). Co-localization of the ER protein calreticulin with the R3.36A GPR35 suggests that the R3.36A mutant is retained intracellularly at the ER. This is consistent with the lack of cell surface receptor immunoreactivity for R3.36A (Fig. 8A). Treatment of R3.36A U2OS cells with GPR35 antagonist CID2745687 (5) for 16 h rescued the surface expression of R3.36A, suggesting a chaperone function by this antagonist.
To compare the potencies of these mutants, ␤-arr2-GFP response was quantified as previously described (5). The WT GPR35 has EC 50 values for pamoic acid and zaprinast-induced activation of 52 nM and 1.0 M, respectively ( Table 1). The substitution of Ala4.59 with a glycine did not significantly change the potency for ␤-arr2-GFP response ( Fig. 9 and Table 1). Although modest rightward shifts in the concentration response curves were observed with the K7.40A mutant receptor for pamoic acid and zaprinast (Fig. 9, A and B), they were not significantly different from the WT GPR35 receptor (Table 1). Both R(164)A and the R(164)L mutants showed a robust rightward shift in the concentration response curves for pamoic acid and zaprinast (Fig.  9, C and D, and Table 1). There was no ␤-arr2-GFP response in the R3.36A and R4.60A mutants ( Table 1). The EC 50 value of R6.58A was not obtained because the morphological changes made the quantification of aggregates unreliable.

Agonist-induced ERK1/2 Activation in Mutant GPR35 Cells
In addition to the ␤-arr2-GFP trafficking assay, agonist-induced ERK1/2 phosphorylation was used to examine the functional effects of GPR35 mutations. A downstream effector of GPR35 activation is the G protein-dependent ERK1/2 phosphorylation. It has been shown to be independent of the ␤-arrestin pathway (5). Table 2 summarizes the activities of ERK1/2 by pamoic acid and zaprinast from WT and mutant GPR35 receptors. The EC 50 values of agonists, pamoic acid, and zaprinast on the mutants A4.59G and K7.40A were comparable and not statistically different from the WT receptor. In contrast, a ϳ30-fold increase in the potency of zaprinast was observed in the R6.58A mutant receptor (83 nM) compared with the wildtype GPR35 receptor (2.8 M). A modest increase in potency of pamoic acid in R6.58A (16 nM) was observed as well, compared with WT (62.3 nM).

Calcium Responses in Wild-type and Mutant GPR35 Cells to Agonists
Because GPR35 is a G i/o -coupled receptor, calcium is a downstream effector. Calcium response curves of mutant  R6.58A shifted to the left compared with WT (Fig. 9, E and F).  (Table 1) and ERK1/2 activation (Table 2) results. Quantification of receptor expression levels of WT and mutants (Table 3) showed no significant difference between them, indicating that the expression level is not a contributing    factor for the shifts of response curves of ␤-arr2-GFP, ERK1/2, and calcium.

DISCUSSION
To date, there is no x-ray crystal structure information available for the human GPR35 receptor; thus, in the present study, we have used computational modeling together with site-directed mutagenesis to explore the structural features of human GPR35 involved in agonist-induced receptor activation. The GPR35 sequence contains numerous positively charged amino acids that face into the binding pocket. These residues cluster in two distinct receptor regions: the TMH3-4-5-6 region and the TMH1-2-7 region of GPR35 (Fig. 1).
The most prominent positively charged amino acid in the TMH 1-2-7 region of the GPR35-binding pocket is K7.40, whereas other positively charged residues that face into the binding pocket but are more extracellular include K1.32, R2.65, and R7.33. Mutation of each of these residues to Ala resulted in the retention of WT EC 50 for pamoic acid and zaprinast. Because many of the ligands for GPR35 are anions, a lack of change in potency for any of these TMH1-2-7 residues suggests that the ligand-binding pocket of GPR35 may not lie in the TMH 1-2-7 region of GPR35.
The majority of the positively charged residues that face into the binding pocket are located in the TMH 3-4-5-6 region. These residues include R3.36, R4.60, and R6.58 and the EC-2 loop residues R(164) and R(167). In contrast to results for the TMH 1-2-7-binding pocket region, we found that mutation of many of the positively charged binding site residues in the TMH 3-4-5-6 region affected both pamoic acid and zaprinast signaling. Fig. 10 illustrates the relative positions of pamoic acid and zaprinast in their docked positions as determined by Glide docking studies. In this view from lipid looking between TMH5 and TMH6, it is clear that both ligands bind in the same general area (TMH3-4-5-6 region); however, their positions are shifted relative to each other. Pamoic acid extends further up toward extracellular space than does zaprinast and also moves closer to TMH4.
Many of the interactions identified by docking studies, particularly for pamoic acid, are salt bridge and cationinteractions, which are strong interactions. It is important to emphasize, however, that the presence of water will have a dampening effect on such electrostatic interactions. This would be particularly true for residues in loop regions but is likely true even for transmembrane positively charged residues here, because the high concentration of positively charged residues in the GPR35-binding pocket likely draws waters down into the binding pocket. One therefore cannot expect as profound a change upon mutation and loss of a salt bridge in GPR35 compared with other receptors that may have an analogous salt bridge shielded from water.
Many of the mutations performed here had very striking effects. The EC-2 loop residue R(164) is shown here by mutation and modeling studies to be a key interaction site for pamoic acid, with a significant increase in EC 50 values (R164A 15-fold; R164L, 18-fold) upon mutation to an uncharged residue. R(164) forms a hydrogen bond, rather than a salt bridge, with zaprinast. Mutation of R(164) to an uncharged amino acid would be expected, therefore, to have a less profound effect on the EC 50 of zaprinast when compared with the effect of this mutation on pamoic acid. The 5-fold increase in EC 50 for zaprinast at the R(164)A and R(164)L mutants, therefore, is consistent with this result.
R(167) in the EC2 loop is also shown here by modeling and mutation studies to be a key interaction for pamoic acid, but not for zaprinast. R(167) forms a salt bridge with the part of pamoic acid that is located near the extracellular end of TMH5. This interaction is with the carboxyl group with which R(164) does not interact. The R(167) interaction with pamoic acid is shown to be weaker by mutation studies than the R(164) interaction but still significant with an increase in EC 50 of almost 4-fold for pamoic acid. Zaprinast does not interact with R(167), as shown by modeling, because zaprinast is located further from TMH5. Because the position of zaprinast does not shift toward TMH5 in the R6.58A mutant, no interaction with R(167) is present there as well.
Our modeling studies suggest that at the pamoic acid binding site, R4.60 forms a parallel and a tilted T cationinteraction with both fused rings of pamoic acid. Although R4.60 does not form a very high energy interaction with zaprinast, this residue forms the boundary of the zaprinast-binding pocket and therefore could be anticipated to play an important structural role in the zaprinast-GPR35 R*-binding pocket (Fig. 5). We show here that mutation of Arg4.60 to alanine abolished the agonist-in-  duced receptor activation for both pamoic acid and zaprinast and in two independent assays of ␤-arr2-GFP trafficking and ERK1/2 phosphorylation. This suggests that R4.60(151) is a key residue in the GPR35 R*-binding pocket that may serve to initially align the ligand in the binding pocket either through direct coulombic interactions via a cationinteraction (pamoic acid) or by defining one edge of the binding pocket (zaprinast). R6.58 is located near the extracellular end of TMH6 and points toward the TMH3-4-5-6 region of GPR35. Our modeling studies suggest that this residue does not directly interact with pamoic acid and that the interactions that pamoic acid has in the R6.58A mutant receptor are almost identical, with the exception of the exchange of a direct interaction with Q(157) for R(167). The R6.58A mutation does not cause a large drop in the overall interaction energy of pamoic acid, because interaction with other residues (principally T(166)) improves in the absence of R6.58. Thus, modeling would predict that the R6.58A mutation should have little effect on the efficacy of pamoic acid.
Our modeling studies suggest that R6.58 does not interact with zaprinast either but produces crowding in the binding region. Removal of R6.58 via the R6.58A mutation is actually favorable for zaprinast. The improved receptor interaction results from multiple stronger interactions at L6.55, R(164), L3.33, and L6.51. In each case, the gains are largely due to improved van der Waals interactions. Thus, modeling would predict that zaprinast should have a lower EC 50 in the R6.58A mutant.
Mutation studies showed that R6.58A increases the potency of zaprinast in the pERK assay. On the other hand, it causes a roundup of cells upon agonist treatment in the ␤-arr2-GFP trafficking assay. Because of these morphological changes upon agonist application, a reliable quantification of arrestin trafficking could not be obtained for this mutant. This morphological change upon agonist application was not observed in WT and other mutant cell lines. Interestingly, the calcium response assay showed increased potencies for both zaprinast and pamoic acid in the 6.58A mutant, similar to what was observed in the pERK assay.
Modeling studies suggest that R3.36 is a major contributor to the interaction energy of pamoic acid at GPR35 R*. This residue also interacts with zaprinast, but in a lower energy interaction. An R3.36A mutation can be expected to impact both ligands, but pamoic acid more significantly. Mutation of Arg3.36 to alanine prevented GPR35 receptor expression on the cell surface; it was trapped in the endoplasmic reticulum. This suggests that R3.36A mutation may cause a structural misfolding during the protein maturation process. This is consistent with the ability of the GPR35 antagonist to rescue cell surface expression of R3.36A (Fig. 8D). A similar chaperone function by ligands has been previously documented for V2 vasopressin receptor (18). In the functional assay, there was no ␤arr2-GFP trafficking response for this mutant upon pamoic acid and zaprinast application. Lack of ␤arr2 recruitment response for the mutant R3.36A was reported in another study as well (49). The total loss of ␤arr2 response for R3.36A may be due to the lack of surface expression. Milligan and co-workers (49) have reported that the R3.36A mutation abolished ␤arr2 recruitment in response to several GPR35 agonists. A recent mutant and modeling study by MacKenzie et al. (54) also indicated the contributions of R4.60, R(164), and R6.58 to ligand binding.
The A4.59G mutant is a random mutation created as a control. This residue actually faces lipid in the GPR35 R* model. Mutation studies show here that this mutation did not change the EC 50 of pamoic acid or zaprinast relative to WT.
The Milligan group has reported that Y3.32 is also an important residue for the action of zaprinast (49). Modeling studies reported here suggest that Y3.32 is a part of both the pamoic acid-and zaprinast-binding pockets at GPR35 R*, with its interaction being stronger for zaprinast. In the present study, however, function in response to pamoic acid or zaprinast could not be detected at a human Y3.32L mutant.
Our data strongly support that direct interactions of the ligands with the mutated residues is responsible for the observed functional effects. However, in the absence of radioligand binding experiments, direct interactions are implied and not proven at this point. A recent mutant and modeling study by MacKenzie et al. (54) also indicated the contributions of R4.60, R(164), and R6.58 to ligand binding using functional assays.
Recently it was reported that zaprinast is 38-fold more potent for rat GPR35 than human GPR35 (55). Kynurenic acid also has a higher potency for rat than human GPR35 (15). These findings imply that the ligand-binding site in human versus rat might be distinct. The origins of this difference are currently under investigation.