Berberine potently inhibits protein tyrosine phosphatase 1B: Investigation by docking simulation and experimental validation

Abstract

Berberine was investigated as an inhibitor of human protein tyrosine phosphatase 1B (h-PTP 1B) in an attempt to explain its anti-hyperglycemic activitiy. The investigation included simulated docking experiments to fit berberine within the binding pocket of h-PTP 1B. Berberine was found to readily fit within the binding pocket of h-PTP 1B in a low energy orientation characterized with optimal electrostatic attractive interactions bridging the isoquinolinium positively charged nitrogen atom of berberine and the negatively charged acidic residue of ASP 48 of h-PTP 1B. Experimentally, berberine was found to potently competitively inhibit recombinant h-PTP 1B in vitro (Ki value = 91.3 nM). Our findings strongly suggest that h-PTP 1B inhibition is at least one of the reasons for the reported anti-hyperglycemic activities of berberine.

Keywords::

• Berberine
• tyrosine phosphatase
• h-PTP 1B
• diabetes
• docking simulations
• inhibition

Introduction

Diabetes Mellitus can result in frequent and serious complications including macrovascular and microvascular complications. Insulin independent diabetes mellitus (Type 2 diabetes) accounts for more than 90% of diabetic cases [1]. Resistance to the biological actions of insulin in tissues like muscle, liver and adipocytes is a major feature of the pathophysiology in human obesity and in type 2 diabetes [2].

Insulin resistance in diabetes can result from a defect in the insulin signaling pathway at an early step in the signal transduction cascade post to insulin binding to its receptor. The interaction of insulin with its receptor leads to tyrosine phosphorylation of specific intracellular proteins and also within the insulin receptor kinase (autophosphorylation). The reversible tyrosine phosphorylation of the insulin receptor and its cellular substrate proteins play a central role in the mechanism of insulin action [3].

Berberine (Figure 1), which is an isoquinoline alkaloid of wide distribution in nature and widely used in traditional eastern homeotherapy, particularly in treating gastrointestinal infections 4-6, has been extensively reviewed in the literature as an interesting natural compound of wide potential medicinal applications 7-12.

Figure 1 Structure of berberine.

Berberine has been reported to possess potent anti-diabetic activity 13-17. The hypoglycemic effects of berberine were accidentally discovered when it was administered to a diabetic patient with diarrhoea [18]. Since then berberine has often been used as an anti-hyperglycemic agent by many Chinese physicians.

From studies on mice and hepatocyte cell lines, berberine has been shown to antagonize the hyperglycemic action of glucose and the gluconeogenic action of alanine [16]. It seems to improve insulin resistance by raising insulin sensitivity [15,19] without affecting insulin secretion [13,14]. Interestingly, some recent reports suggested that berberine might inhibit the intestinal absorption of glucose [20]. Nevertheless, most accumulated scientific evidence suggests that the hypoglycemic action of berberine is related to post-insulin receptor mechanisms 13-1619.

Protein tyrosine phosphatase 1B (h-PTP 1B), a cytosolic non-receptor PTPase, has been implicated as a negative regulator of insulin signal transduction. h-PTP 1B directly catalyzes the dephosphorylation of cellular substrates of the insulin receptor kinase resulting in a down regulation of insulin action 21-23.

In a recent h-PTP 1B knockout study, it has been reported that mice lacking functional h-PTP 1B showed an enhancement of tyrosine kinase activity and increased insulin sensitivity [24]. The levels of h-PTP 1B expression in muscle and adipose tissues in humans were strongly correlated to insulin resistance states [25]. In another animal-based approach, treatment with anti-sense oligonucleotide specific for h-PTP 1B resulted in normalization of blood glucose and insulin levels in animal models of type 2 diabetes [26]. Accordingly, h-PTP 1B inhibitors could increase insulin receptor tyrosine phosphorylation, mimic cellular and in vivo actions of insulin, and lower plasma glucose in diabetic animal models 27-30.

The anti-hyperglycemic action observed for berberine, combined with recent interest in human h-PTP 1B as an antidiabetic target, prompted us to evaluate possible binding interaction(s) between berberine and h-PTP 1B as a new suggested mechanism for the hypoglycemic action of this interesting natural compound. We decided to evaluate the binding interactions by employing computer-aided molecular docking and scoring. Eventually, the findings of the docking simulation study were experimentally validated by in vitro bioassay against human recombinant h-PTP 1B. Additional kinetic studies were conducted to evaluate the type of h-PTP 1B inhibition by berberine, i.e. competitive, non-competitive or uncompetitive.

Simulated molecular docking is basically a conformational sampling procedure in which various docked poses/conformations are explored to identify the correct one. This process can be a very challenging problem given the degree of conformational flexibility at the ligand-macromolecular level 31-36. Docking consists of two parts, namely, (i) prediction of the conformation, orientation and position (pose) of the bioactive compound into the binding pocket, and (ii) estimation of the tightness of target-ligand interactions (scoring) to guide conformational sampling [37]. The final docked conformations are selected according to their scores. We decided to conduct the docking study utilizing the program Ligandfit® [38] which was recently reported to illustrate good overall performance, particularly in virtual high-throughput screening experiments (vHTS) [36,39].

Materials and methods

Materials

All of the chemicals used in these experiments were of reagent grade and obtained from commercial suppliers: phosphopeptide Asp-Ala-Asp-Glu-phosphoTyr-Leu-Ile-Pro-Gln-Gln-Gly (from BIOMOL, USA), assay buffer components (HEPES, NP-4, DTT and EDTA, BIOMOL, USA), recombinant h-PTP 1B (BIOMOL, USA), glycerol (SIGMA, USA), BSA (BIOMOL, USA), Berberine (SIGMA, USA), DMSO (SIGM, USA), standard h-PTP 1B inhibitor (RK-682, BIOMOL, USA).

Molecular modeling

Hardware and software

Docking and scoring modeling studies were performed using the CERIUS2® suite of programs (version 4.8, Accelrys Inc., San Diego, California, www.accelrys.com) installed on a Silicon Graphics Octane2 desktop workstation equipped with a 600 MHz MIPS R14000 processor (1.0 GB RAM) running the Irix 6.5 operating system.

Preparation of h-PTP 1B crystal structure

The 3D coordinates of h-PTP 1B were retrieved from the Protein Data Bank (PDB code: 1g7f)[40]. The selected structure is of the best 3D resolution (1.80 Å) compared to other available h-PTP 1B structures. Hydrogen atoms were added to the protein utilizing CERIUS2® templates for protein residues. Gasteiger charges were assigned to the protein atoms as implemented within LigandFit® [38]. The protein structure was utilized in subsequent docking experiments without energy minimization. Explicit water molecules were kept in the structure.

Preparation of berberine structure

The chemical structure of berberine was sketched in Chemdraw Ultra (6.0), saved in MDL molfile format and imported into CERIUS2®. Berberine atoms were assigned partial charges using the default Gasteiger method (polygraph 1.0) implemented in CERIUS2® [41]. Finally, the structure was energy-minimized employing the UNIVERSAL force field® (version 1.02, default settings) implemented within CERIUS2®.

Docking simulations

LigandFit® considers the flexibility of the ligand and treats the receptor as rigid. There are two steps implemented in the LigandFit® process:

(1) Defining the location(s) of potential binding site(s) by shape-based search for cavities in the protein. The algorithm for cavity detection calculates a rectangular grid enclosing the protein, cavity regions, and explicit water molecules around the complex. The protein is mapped on the grid. All grid points occupied by the protein (or crystallographically explicit water molecules) are not available in the site search. The unoccupied grid points inside the protein are potential binding sites. However, if a certain ligand is co-crystallized with the targeted protein then it is possible to generate the binding site from the docked ligand by collecting all grid points that lie within the radius of any atom of the ligand to form the binding site [38].

In the current docking experiments the binding site was generated from the co-crystallized ligand with the targeted protein. The grid resolution was set to 0.5 Å, the radius of hydrogen atoms in the co-crystallized ligand and the protein was set 2.0 Å, while the radius of heavy atoms (carbon, nitrogen, oxygen and sulfur) in the co-crystallized ligand and the protein was set to 2.5 Å. Figures 2 and 3 show the relative size of the binding pocket and the co-crystallized ligand (Pnu177496) that was employed for defining the binding pocket [40].

Figure 2 Structure of Pnu177496.

Figure 3 (A) Perspective cartoon view of h-PTP 1B co-crystallized with Pnu177496 (in red). (B) Cartoon view of the h-PTP 1B showing the binding pocket as dotted red area.

(2) Docking the ligands in the binding site. In LigandFit® docking is composed of few major substeps: [38] (i) conformational search of the flexible ligand employing Monte Carlo randomized process. (ii) pose/conformation selection based on shape similarity with the binding site. (iii) candidate conformers/poses exhibiting low shape discrepancy are further enrolled in calculation of the dock and interaction energies. The dock energy is composed of two terms, namely, the internal energy of the ligand and the interaction energy with the receptor, summarized by van der Waals and electrostatic energy terms. To improve the time-consuming computation of the interaction energy, an approximation by grid-based interpolation is employed in Ligandfit®. A grid encloses the site, and at each point of the grid the potentials are computed for the active site. The potentials at the ligand atom locations are subsequently interpolated. (iv) each docked conformation/pose is further fitted into the binding pocket through a number of rigid-body minimization iterations, i.e., minimization of the interaction energy via molecular rotations and translations of the docked ligand and (v) docked conformers/poses that have docking energies below certain user-defined threshold are subsequently clustered according to their RMS similarities. Representative conformers/poses are then selected, further energy-minimized within the binding site and saved for subsequent scoring.

In the current docking experiment the following docking configuration was employed

• Monte Carlo search parameters: number of trials = 10000; search step for torsions with polar hydrogens = 30.0 degree.

• The RMS threshold for ligand-to-binding site shape match was set to 2.0 employing a maximum of 5.0 binding site partitions.

• Interaction energy parameters: the interaction energies were assessed employing CFF force field, dielectric constant = 1, non-bonded cutoff 10.0 Å. An energy grid extending 3.0 Å from the binding site was implemented. The interaction energy was estimated by tri-linear interpolation value using soft potential energy approximations [38].

• Rigid body ligand minimization parameters: 100 iterations of rigid body minimization (molecular translational and rotational movements) were applied to every orientation of the docked ligand.

• The docked conformations/poses of calculated interaction energies ≤ 30.0 Kcal/mol were clustered using the complete linkage algorithm in CERIUS2® with RMS similarity threshold of 2.0 Å. The best member within the cluster was selected and was further energy-minimized within the binding site for a maximum of 100 rigid-body iterations and 500 flexible conformation iterations. Eventually, a maximum of 20 optimal conformers/poses were saved for each molecule for subsequent scoring. However, this docking procedure yielded 11 plausible docked conformers/poses for berberine.

Scoring of docked conformers/poses

The best-docked conformers/poses (11 structures) were scored using the following scoring functions: LigScore1, LigScore2, LUDI, PLP1, PLP2, PMF and JAIN. However, only the best molecular orientation ranked by a consensus score based on the 7 scoring functions was evaluated. The consensus function assigned a value of 1 for any molecular pose ranked within the highest 40% by the particular scoring function, otherwise it assigned zero for the pose, i.e., if it was within the lowest 60%. Subsequently, the consensus function summed up the scores for each molecular pose/conformer and ranked the molecular orientations accordingly. The followings describe briefly the different scoring functions employed in the current study:

(1) LigScore1 and LigScore2: [38,42]. These are slightly different scoring functions implemented in CERIUS2®. According to the description in the CERIUS2 user's manual, these scores represent the sum of three terms, Where vdW is a softened Lennard-Jones 6–9 potential, C+pol is a count of the buried polar surface area between the complex involving attractive protein-ligand interactions, Totpol2 is the square of the buried polar surface area between the complex involving both attractive and repulsive protein-ligand interactions.

In our study, LigScore1 and LigScore2 scores were calculated employing CFF force field (version 2002) and using grid based energies with a grid extension of 5.0 Å across the binding site. The CERIUS2 user's manual mentions that these scoring functions were calibrated by fitting to known protein-ligand binding affinities. Thus, they fall into the empirical scoring function category.

(2) PLP1 and PLP2 [43,44]. These are two closely related versions of the same scoring function (PLP) as they use slightly different algorithms and parameters sets. Generally, PLP scoring function is a sum of pair-wise linear potentials between ligand and protein heavy atoms with parameters dependent on interaction type. It can be expressed conceptually as:

Ligand and protein heavy atoms are classified as hydrogen bond donors, acceptors, donors/acceptors, or nonpolar. Each pair of interacting atoms is then assigned one of the three interaction types: hydrogen bonding between donors and acceptors (EH-bond), repulsive donor-donor and acceptor-acceptor contacts (Erepulsion), and generic dispersion of other contacts (Econtact). Both the hydrogen bonding and repulsive terms are modulated by a scaling factor that imparts a crude distance and angular dependence. Small (fluorine and metal ion), medium (carbon, oxygen and nitrogen), and large (sulfur, phosphorus, chlorine, and bromine) atoms are assigned atomic radii of 1.4, 1.8 and 2.2 Å, respectively. These parameters are derived from interatomic distances observed from a large number of high-quality crystal structures.

(3) PMF: this knowledge-based potential of mean force scoring function is based on the work of Muegge et al. 45-47 who analyzed 697 protein-ligand complex structures from the Protein Data Bank and derived a set of distance-dependent interaction potentials (dij) for various atom pairs. Both enthalpic and entropic effects are assumed to be implicitly included in this potential. The protein-ligand interaction energy is then defined as a sum of potentials over all heavy atom pairs (Aij) between the complex (ignoring hydrogen atoms): In this study, the PMF scores were calculated employing cutoff distances for carbon-carbon interactions and other interactions of 12.0 Å.

(4) LUDI: this empirical scoring function was developed by Böhm and is one of the pioneering empirical scoring functions [48,49]. It dissects protein-ligand binding free energy as

The first two terms account for the hydrogen bonds and ionic interactions formed between the complex. The contribution of each hydrogen bond is scaled by a distance- and angle-dependent function (R and α, respectively) in order to penalize the deviations from an ideal geometry. The third term accounts for the hydrophobic effect, which calculates the buried hydrophobic molecular surface (A_hydrophobic). The fourth term counts all the rotable single bonds (N_rotor) in the ligand, which is supposed to be related to the torsional entropy loss of the ligand upon protein-ligand complexation. The last term is a regression constant. This scoring function was calibrated by fitting known dissociation constants of 87 protein-ligand complexes.

(5) JAIN: it is an empirical scoring function, which is a sum of five interaction terms: (i) lipophilic interactions, (ii) polar attractive interactions, (iii) polar repulsive interactions, (iv) solvation of the protein and ligand and (v) an entropy term for the ligand.

Only close protein-ligand atoms are considered for the pair-wise interaction terms. The lipophilic and polar interaction terms are captured by a weighted sum of Gaussian-like and sigmoid functions. This functional form is very short-ranged and leads to a pronounced maximum at close surface contacts. It also results in a penalty for large overlaps between protein and ligand atoms [50].

Measurement of PTPase inhibition

The detection of free phosphate released into the medium by the phosphatase enzyme is based on the classic Malachite green ammonium-molybdate method. The assay as described previously [51] and adapted for the plate reader, was used for the nanomolar detection of liberated phosphate by recombinant h-PTP 1B.

The assay used the phosphopeptide Asp-Ala-Asp-Glu-phosphoTyr-Leu-Ile-Pro-Gln-Gln-Gly as substrate. This phosphorylated peptide corresponds to the 988-998 catalytic domain of epidermal growth factor receptor (EGFR), and it is one of the most efficient peptide substrates known for h-PTP 1B [52]. The phosphopeptide was diluted with the assay buffer, pH 7.2, containing 50.0 mM HEPES, 1.0 mM DTT, 1.0 mM EDTA, and 0.05% NP-40 to obtain 150 μM substrate concentration.

The recombinant h-PTP1B was diluted with the buffer solution pH 7.2, glycerol (50 mg/ml), and of BSA (1 mg/ml) to yield a final enzymatic solution with approximate activity 30 pmol/min.

Berberine was dissolved in DMSO and diluted with the assay buffer. The DMSO concentration was less than 1% in all experiments and controls.

The diluted enzyme (5 μl) was pre-incubated with 35 μl of the assay buffer with or without berberine for 15 min at 25°C. Then the reaction was initiated by the addition of 5 μl of peptide substrate solution to get a final concentration of 150 μM peptide substrate.

The mixture was equilibrated to 25°C and incubated for 30 min. The reaction was terminated by the addition of 100 μl malachite green ammonium molybdate-tween 20 reaction terminating reagent. Negative controls were prepared by adding the substrate after the addition of reaction terminating reagent. The Controls were prepared using the pre-incubated recombinant enzyme without the addition of the berberine.

A standard h-PTP 1B inhibitor, RK-682, was included as enzyme activity control for screening inhibition. RK-682 is a potent tyrosine phosphatase inhibitor isolated from Streptomyces sp. 88-682 [53]. The color was allowed to develop at room temperature for 30 min, and the sample absorbance was determined at 630 nm using a plate reader (Bio-Tek instruments ELx 800, USA). Samples and blanks were prepared in duplicates.

Inhibition of recombinant h-PTP1B by berberine was calculated as a percent activity of the uninhibited phosphatase control.

The percent activity was plotted against the logarithmic transformation of the corresponding berberine concentrations for determining the IC₅₀ value. The standard inhibitor RK-682 was used at a final concentration of 100 μM and the h-PTP 1B shower about 20% activity.

Mode of Berberine Inhibitory Action

The Michaelis-Menten parameters were determined from a direct fit of the double reciprocal plot of velocity versus substrate concentration data to the Michaelis-Menten equation (i.e., Lineweaver-Burk plot). The reaction mixture consisted of different concentrations of h-PTP 1B substrate (2–80 μM) in the presence of different berberine concentrations (0, 60, and 110 nM). The inhibition constant (Ki) was determined for h-PTP 1B catalyzed hydrolysis by fitting the data to the following equation: where, V_o and V_max are the initial and maximum enzymatic velocities, respectively, [S] and [I] are the substrate and inhibitor concentrations, respectively, Km and Ki are the Michaelis-Menten and inhibition constants, respectively.

Results and discussion

Despite accumulating literature regarding berberine, the exact mechanism by which it produces its antihyperglycemic activity is still controversial [14,20,54,55]. Accordingly, the fact that h-PTP 1B inhibitors are potent hypoglycemic agents prompted us to evaluate any possible inhibitory effect for berberine on h-PTP 1B.

Our efforts commenced by evaluating the possibility of binding via computer-aided molecular modeling techniques. Accordingly, we docked berberine into the binding pocket of h-PTP 1B (PDB code: 1g7f). The binding site was defined from the crystallographic structure of a bound high-affinity ligand (Figures. 2 and 3, see Docking Simulations under Methods). Docking simulation suggested eleven discrete binding modes for berberine within h-PTP 1B. Table I shows the scores produced by the above mentioned 7 scoring functions for the different molecular poses of berberine within h-PTP 1B.

Table I.  the score values for each molecular pose suggested by the LigandFit® docking engine.

	Scoring Functions
Pose	Ligscore 1	Ligscore 2	-PLP1	-PLP2	JAIN	-PMF	LUDI	Consensus Score
1	3.87	4.46	64.74	64.1	− 0.33	73.19	473	7
2	3.57	3.69	47.99	49.83	− 0.98	83.00	416	6
3	2.62	4.14	55.38	53.1	− 0.93	69.38	538	5
4	3.20	2.59	32.85	39.06	− 0.81	61.71	400	4
5	3.97	3.30	39.44	42.67	− 0.97	42.77	392	3
6	2.93	3.58	37.87	35.37	− 0.8	40.31	336	2
7	1.93	3.41	38.58	38.14	− 0.01	46.37	285	1
8	2.44	3.47	35.8	39.13	− 0.83	54.99	353	0
9	2.37	3.43	19.06	17.46	− 1.73	44.48	259	0
10	1.22	3.29	34.95	35.51	− 1.23	38.94	266	0
11	1.74	3.08	24.61	26.06	− 2.06	11.71	143	0

The molecular interactions of the highest ranking binding mode can be summarized in Figures 4 and 5. Clearly from the figures, the cationic nitrogen of berberine interacts with the negative charged carboxylate moiety of ASP 48. On the other hand, the dimethoxylated aromatic ring in berberine is situated at around 4.0 Å from the guanidine moiety of ARG 24 suggesting significant charge transfer attraction between the electron-rich aromatic ring and the positive charged guanidine group. Finally, the fact that the dioxymethylene substituted aromatic ring of berberine is nearly 4.0 Å away from the phenyl ring of TYR 46, suggests the existence of van der Waals' aromatic stacking attraction between them. Overall, the three attractive interactions cooperate in stabilizing the proposed complex.

Figure 4 The highest ranking binding mode of berberine as suggested by the consensus scoring function. (A) Berberine structure docked into the water accessible surface (Connolly's surface) within the binding pocket of h-PTP 1B. The surface is represented as yellow dotted cavity (B) Perspective cartoon view of h-PTP 1B and the docked berberine structure (in red).

Figure 5 Detailed view of the docked berberine structure and the corresponding interacting amino-acid moieties within the binding site of h-PTP 1B.

However, to validate our docking-scoring procedure, we employed the same conditions to dock a well known h-PTP 1B inhibitor (Pnu177496 [40], Figure 2) into the binding pocket of this enzyme. Our docking simulation resulted in a very close model to the crystallographic structure, which supports our conclusions regarding berberine/h-PTP 1B binding (Figure 6). Interestingly, ASP 48 seems to play a central role in the binding of Pnu177496 within the h-PTP 1B binding site, as it forms two hydrogen bonds with the central amido nitrogens of Pnu177496.

Figure 6 The highest ranking binding mode of inhibitor Pnu177496 as suggested by our docking-scoring conditions. (A) Detailed view of the docked structure and the corresponding interacting amino-acids. (B) Comparison between the docked conformer/pose of inhibitor Pnu177496 (green) as produced by the docking simulation and the crystallographic structure of this inhibitor within h-PTP 1B (red, PDB code: 1g7f)

Further validation of the docking procedure comes from that fact that a similar docking-scoring approach has been recently successfully implemented to construct self-consistent and predictive protein-aligned comparative molecular field analysis models (CoMFA) [56].

The proposed inhibitory action of berberine was experimentally validated against recombinant h-PTP 1B using a phosphotyrosyl decapeptide (Asp-Ala-Asp-Glu-phosphoTyr-Leu-Ile-Pro-Gln-Gln-Gly) as substrate, since it corresponds to the 988-998 catalytic domain of the epidermal growth factor receptor. The enzymatic reaction progression was monitored through the release of inorganic phosphate. The in vitro activity was expressed as the concentration of berberine that inhibited enzyme activity by 50% (IC₅₀). Figure 7 clearly illustrates the potent inhibitory action of berberine on h-PTP 1B leading to an IC₅₀ of 156.9 nM.

Figure 7 The effect of variable berberine concentrations on the relative activity of h-PTP 1B. Data are expressed as means of three replicates ± standard deviation of measurements.

To investigate the type of inhibition of berberine, we evaluated h-PTP 1B kinetics at different berberine concentrations. Subsequently, the enzymatic velocities were calculated and plotted against substrate concentrations within a Lineweaver-Burk plot (Figure 8) [57]. As evident from the figure, berberine increased the Km values of h-PTP 1B without altering the Vmax value, which clearly suggests a competitive mode of inhibition such that berberine blocks the catalytic pocket of the enzyme thus supporting our docking results regarding the site of berberine binding within h-PTP 1B, i.e., within the enzyme's catalytic site. The Ki value of berberine was estimated to be 91.3 nM.

Figure 8 A Lineweaver-Burk plot of berberine inhibition of h-PTP 1B. The inhibition kinetics were determined at various fixed concentrations of berberine: 0 nM (•), 60 nM (▪), 110 nM (▴). Each point represents the average of three replicates.

Conclusion

We have unequivocally proved through experimental testing and theoretical docking simulations that berberine is a potent competitive inhibitor of h-PTP 1B, which explains, at least partially, the reported use of berberine as an effective anti-hyperglycemic agent. Furthermore, this finding does not exclude other possible hypoglycemic mechanisms of berberine. However, further studies are probably required to evaluate the significance of this mechanism in vivo.

Acknowledgements

This project was sponsored by the Deanship of Scientific Research at the University of Jordan (grant number 19/2003–2004). The authors wish to thank the Deanship of Scientific Research and Hamdi-Mango Center for Scientific Research at the University of Jordan for providing funds towards purchasing O2 and Octane2 Sgi workstations and CERIUS2® software package.