Hydroxylated derivatives of NPC1161: theoretical insights into their potential toxicity and the feasibility and regioselectivity of their formation.

For antimalarial 8-aminoquinoline (8-AQ) drugs, the ionization potential (energy required to remove an electron) of their putative metabolites has been proposed to be correlated in part to their hemotoxicity potential. NPC1161 is a developmental candidate as an 8-AQ antimalarial drug. In this work, the ionization potentials (IPs) of the S-NPC1161 (NPC1161a) hydroxylated derivatives, which are possible metabolites derived from action of endogenous cytochrome P450 (CYP450) enzymes, were calculated at the B3LYP-SCRF(PCM)/6-311++G**//B3LYP/6-31G** level in water. The derivative hydroxylated at N1' (8-amino) was found to have the smallest IP of ∼ 430 kJ/mol, predicting that it would be the most hemotoxic. The calculated IPs of the derivatives hydroxylated at the C2 and C7 positions were ∼ 475 and ∼ 478 kJ/mol, respectively, whereas the calculated IPs of those hydroxylated at all other possible positions were between 480 and 490 kJ/mol. The homolytic bond dissociation energies (HBDEs) of all C-H/N-H bonds in NPC1161a were also calculated. The smaller HBDEs of the C-H/N-H bonds on the 8-amino side chain suggest that these positions are more easily hydroxylated compared to other sites. Molecular orbital analysis implies that the N1' position should be the most reactive center when NPC1161 approaches the heme in CYP450.

Introduction

Malaria is a mosquito-borne infectious disease that is widespread in the tropical and subtropical regions, causing extensive mortality and morbidity. Global malarial deaths increased to a peak of 1 817 000 in 2004, decreasing to 1 238 000 in 2010, mostly affecting those under the age of 5 years.[1] Primaquine (PQ, Chart 1), an 8-aminoquinoline (8-AQ) first synthesized in 1946,[2] is the only drug approved by the U.S. Food and Drug Administration (FDA) for the radical cure of relapsing malaria.[3] It is also being used as a prophylactic against all major forms of human malaria[4] and was found to have significant sporontocidal and gametocytocidal activity.[5] The greatest life-threatening problem for PQ and other 8-AQ antimalarials is the aggravated methemoglobinemia[6] (in which hemoglobin is oxidized to form methemoglobin, an Fe3+ protein that is unable to carry oxygen) and hemolysis[7] that occurs in patients deficient in glucose-6-phosphate dehydrogenase (G6PD). The hemolysis of 8-AQ analogues can be induced by their hydroxylated metabolites, which are formed by the action of cytochrome P450 isoforms.[8−13] Among a number of metabolites of PQ shown in Chart 1,[14−21] it has been found that 5-hydroxyprimaquine (5-OH-PQ), 5,6-dihydroxyprimaquine (5,6-2OH-PQ) and 6-methoxy-8-(N-hydroxy)aminoquinoline [6-MeO-8-(N-OH)-AQ] are more toxic than the parent PQ molecule.[20,21] In contrast, carboxyprimaquine (cPQ) and 6-methoxy-8-aminoquinoline (6-MeO-8-AQ) are less toxic because they generate less methemoglobin than PQ.[13,21,22]

Chart 1

8-[(4-Amino-1-methylbutyl)amino]-5-[3,4-dichlorophenoxy]-4-methylquinoline (NPC1161, Chart 1), an 8-AQ analogue, is a developmental candidate with better antimalarial efficacy and lower toxicity than PQ.[23,24] It is also significantly active against Leishmania[9,24−26] and Pneumocystis carinii infections.[9,24] In addition, its (−) enantiomer (NPC1161b) showed more markedly reduced general toxicity in mice and reduced hematotoxicity in the dog model of methemoglobinemia than its (+) enantiomer (NPC1161a).[9,24] However, the properties of the hydroxylated NPC1161 metabolites, particularly their potential to cause methemoglobinemia, is not clear. In this work, we evaluated the potential toxicity as well as the feasibility and regioselectivity of formation of hydroxylated metabolites of NPC1161 using computationally intensive approaches. The properties we studied in this work are identical for the two enantiomers, so the calculations we report were performed on just one enantiomer, NPC1161a.

Computational Methods and Software Packages

The MMFF94s force field in Schrödinger’s MacroModel[25] was employed for the conformational search of neutral NPC1161a with an energy window of 130 kJ/mol. The 25 conformers with the lowest energies were then subjected to geometry optimization and single-point energy calculation at the B3LYP[27−29]/6-31G** and B3LYP/6-311++G**//B3LYP/6-31G** levels in the gas phase and at the B3LYP-SCRF(PCM)/6-311++G**//B3LYP/6-31G** level in water. Twenty conformers were found after all of the B3LYP/6-31G** geometry optimizations, and their conformational distributions were calculated using the total energies obtained at the above levels, including zero-point vibrational energies (ZPVEs) calculated at the B3LYP/6-31G** level. The IPs of the predominant conformers of NPC1161a were computed at the above levels. The IPs of metabolites singly hydroxylated at all possible positions were also calculated using the above protocols, and the most predominant conformer of neutral NPC1161a was selected as a model. Potential energy surfaces (PESs) were scanned at the AM1 (semiempirical molecular model) or B3LYP/6-31G** levels to locate the energetic minima of the metabolites of NPC1161a hydroxylated at the C2, C7, and N1′ positions. Homolytic bond dissociation energies (HBDEs) of the selected NPC1161a conformer were calculated at the B3LYP/6-31G** and B3LYP/6-311++G**//B3LYP/6-31G** levels in the gas phase. All computations at the semiempirical and quantum mechanical levels were performed using the Gaussian 09 software package.[30] If not mentioned elsewhere, the results from the solvation calculations at the B3LYP-SCRF(PCM)/6-311++G**//B3LYP/6-31G** level are presented below and discussed, whereas those at all other levels are provided in the Supporting Information.

Results and Discussion

Conformational Analysis of NPC1161a

In our previous studies,[31−34] we proposed that the ability to lose an electron is correlated in part to the hemotoxicity of antimalarial 8-AQ drugs. In this work, we calculated the ionization potentials (IPs) of antimalarial candidate NPC1161a to evaluate its possible hemotoxicity. To locate the most favorable conformers to be used to calculate the IPs for this very flexible compound, a Monte Carlo random conformational search at the MMFF94s level of theory was performed, using Schrödinger’s MacroModel software package,[35] yielding 453 conformers within an energy window of 130 kJ/mol (Figure S1, Supporting Information). The 25 conformers with the lowest energies within an energy cutoff of 8.4 kJ/mol (∼2 kcal/mol) were submitted to full hybrid density functional theory (DFT) geometry optimization, leading to 20 conformers located at the B3LYP/6-31G** level in the gas phase (Figure 1). Harmonic vibrational frequencies were calculated at the same level to confirm that they were minima on the potential energy surface, and single-point energies were computed at the B3LYP/6-311++G**//B3LYP/6-31G** level in the gas phase and at the B3LYP-SCRF(PCM)/6-311++G**//B3LYP/6-31G** level in water. The rotation of the phenolic group about C5–O (C1″) is the major geometric difference between conformers 1–10 and 11–20, in which the dihedral angles of C6–C5–O–C1″ are about −77° and 77°, respectively. This rotation is the only significant geometric difference between the conformers 01 and 11, 02 and 20, 03 and 12, 05 and 14, and 10 and 19. Among the above conformers, weak CH···N and/or NH···N hydrogen bonds were found in conformers 03, 05, 10, 12, 14, and 19. The distances of N6′···H(C7) and N6′···H(C2′) in conformers 03 and 12 are 2.67 and 2.54 Å, respectively; those of N6′···H(N1′), N6′···H(C3′), and N1′···H(C5′) in 05 and 14 are 2.23, 2.65, and 2.62 Å, respectively; and those of N6′···H(C2′) in 10 and 19 are each 2.42 Å. Conformational analysis (Table 1) indicated that these conformers with hydrogen bonds are negligible in water, as evidenced by the fact that they accounted for only 4% at the B3LYP-SCRF(PCM)/6-311++G**//B3LYP/6-31G** level in water. The conformers with the most extended side chain, specifically 02 and 20, were found to be the major conformers, accounting in total for 45% at the B3LYP-SCRF(PCM)/6-311++G**//B3LYP/6-31G** level in water.

Figure 1

Optimized geometries of conformers of neutral NPC1161a at the B3LYP/6-31G** level in the gas phase.

Table 1

Conformational Analysis of Neutral and Ionized NPC1161a at the B3LYP-SCRF(PCM)/ 6-311++G**//B3LYP/6-31G** Level in Water

	neutral NPC1161a		ionized NPC1161a
conformer no.	ΔEa	P%b	ΔEa	P%b
01	2.94	7.0	4.09	5.4
02	0.00	23.0	0.24	25.5
03	7.80	1.0
04	5.11	2.9	6.11	2.4
05	10.01	0.4
06	9.64	0.5
07	19.00	0.0
08	1.47	12.7	1.83	13.5
09	7.59	1.1
10	8.95	0.6
11	3.22	6.3	3.93	5.8
12	8.42	0.8
13	13.50	0.1
14	9.49	0.5
15	7.35	1.2	6.76	1.8
16	1.14	14.5	1.79	13.7
17	3.79	5.0	5.02	3.7
18	9.73	0.5
19	9.44	0.5
20	0.15	21.6	0.00	28.2

Relative energy (kJ/mol), including zero-point energy at the B3LYP/6-31G** level in the gas phase.

Conformational distribution of individual conformers.

Ionization Potentials of NPC1161a

The geometries of all conformers found in water with a Boltzmann contribution greater than 2% were employed as the starting points to locate the ionized conformers of NPC1161a at the B3LYP/6-31G** level in the gas phase. Single-point energies and the conformational distribution were computed using the same protocols described above. The ionized conformers 02 and 20 were again found to be the major ones, with contributions of 26% and 28%, respectively, in water at the B3LYP-SCRF(PCM)/6-311++G**//B3LYP/6-31G** level (Table 1). The ionization potentials were calculated at the B3LYP-SCRF(PCM)/6-311++G**//B3LYP/6-31G** level in water and are listed in Table 2. The greatest difference in IP values was found to be only 3.14 kJ/mol at the B3LYP-SCRF(PCM, H2O)/6-311++G**//B3LYP/6-31G** levels, indicating that the IPs in this case are relatively insensitive to conformation. Thus, if not mentioned elsewhere, hereafter we employed conformer 02 (NPC1161a02) as a conformer model for the further computations because it was found to be the most predominant one in water.

Table 2

Theoretically Calculated Relative Ionization Potentials (kJ/mol) of Predominant Conformers of NPC1161a at the B3LYP-SCRF(PCM)/6-311++G**//B3LYP/6-31G** Level in Water

conformer no.	ΔIPa
01	2.74
02	1.70
04	1.99
08	0.94
11	3.14
15	0.00
16	1.55
17	1.74
20	1.23

Zero-point energy at the B3LYP/6-31G** level in the gas phase included.

Ionization Potentials of Hydroxylated NPC1161a02

Because hydroxylated metabolites of 8-AQ have been postulated to be responsible for the 8-AQ hemotoxicity,[20,21] ionization potentials of hydroxylated metabolites of the conformer NPC1161a02 were calculated at B3LYP-SRCF(PCM)/6-311++G**//B3LYP/6-31G** level in water. Potential energy surfaces were scanned at the AM1 level in the gas phase, and the predominant conformers of metabolites hydroxylated at the C2, C7, and N1′ positions (Figures S2–S4, Supporting Information) were used for the geometry optimization. Figure 2 depicts the optimized geometries of metabolites of neutral conformer NPC1161a02 hydroxylated at all possible positions at the B3LYP/6-31G** level in the gas phase. Their relative energies (Table 3) and ionization potentials (Table 4) were computed at the B3LYP-SCRF(PCM)/6-311++G**//B3LYP/6-31G** level in water. The metabolite hydroxylated at the C2 position was found to be the most stable, at least 20 kJ/mol more stable than any other, whereas those hydroxylated at 2′Me, 4Me, N6′, and N1′ were found to be the least stable, for both the neutral and ionized states. The calculated IPs indicated that the introduction of a hydroxyl group into NPC1161a02 would have little effect on its ability to lose an electron, except if the hydroxylation were at the N1′ position. Compared to its parent compound, the hydroxylation of NPC1161a02 at the N1′ position would dramatically lower its energy barrier to donating an electron by 44.37 kJ/mol, implying that this hydroxylated metabolite should contribute the most to the potential hemotoxicity of NPC1161a. Interestingly, the hydrogen atom in the hydroxyl group attached at N1′ in the neutral hydroxylated metabolite (OH@N1′ in Figure 2) automatically shifted to the N1 position while being converted into a radical cation, as shown in the potential energy surface scanned at the B3LYP/6-31G** level in the gas phase (Figure 3). The N1′–O bond length decreased from 1.34 Å in the neutral form to 1.29 Å in the radical, and the (N1′)O–H bond dissociated (the distance increased from 1.09 Å in the neutral form to 1.61 Å in the radical), whereas the hydrogen shifted to form the (N1′O)H–N1 bond (the distance was decreased from 1.40 Å in the neutral form to 1.05 Å in the radical). The hemotoxic contribution from the metabolites hydroxylated at the C2 and C7 positions was predicted not to be significant. The hydroxylated metabolites at all other possible positions were found to have calculated IPs indistinguishable from those of their parent.

Figure 2

Optimized geometries of the metabolites of neutral conformer NPC1161a02 hydroxylated at all possible positions at the B3LYP/6-31G** level in the gas phase.

Table 3

Relative Energiesa of Neutral and Ionized Metabolites of Conformer NPC1161a02 Hydroxylated at All Possible Positions at the B3LYP-SCRF(PCM)/6-311++G**//B3LYP/6-31G** Level in Water

	ΔE (kJ/mol)			ΔE (kJ/mol)
OH@b	neutral	ionized	OH@b	neutral	ionized
2	0	0	2′	33.81	49.25
3	37.54	42.46	2′Me	54.27	63.14
4Me	52.08	68.72	3′R	29.4	37.42
2″	41.45	46.88	3′S	39	44.73
5″	33.3	39.97	4′R	32.15	38.01
6″	35.84	42.9	4′S	32.01	36.55
6OMe	31.25	39.93	5′R	20.49	28.82
7	55.31	59.16	5′S	23.53	31.85
N1′	207.83	188.53	N6′	176.69	183.38

Including zero-point energy at the B3LYP/6-31G** level in the gas phase.

Position at which the hydroxylation occurs. Numbering as shown for NPC1161a in Chart 1.

Table 4

Relative Ionization Potentialsa of Metabolites of Conformer NPC1161a02 Hydroxylated at All Possible Positions at the B3LYP-SCRF(PCM)/6-311++G**// B3LYP/6-31G** Level in Water

OH@b	ΔIP (kJ/mol)	OH@b	ΔIP (kJ/mol)
parent	0	2′	7.59
2	–6.36	2′Me	0.55
3	–1.63	3′R	0.76
4Me	8.78	3′S	–0.84
2″	–1.71	4′R	1.16
5″	–0.01	4′S	–1.54
6″	–0.28	5′R	1.6
6OMe	1.93	5′S	1.83
7	–3.27	N6′	0.06
N1′	–44.37

Including zero-point energy at the B3LYP/6-31G** level in the gas phase.

Position at which hydroxylation occurs. Numbering as shown for NPC1161a in Chart 1.

Figure 3

Potential energy surface of the hydrogen shift in the metabolite of NPC1161a02 hydroxylated at the N1′ position while being converted to a radical cation (scanned in the gas phase at the B3LYP/6-31G** level; key distances are shown in Ångstroms).

Feasibility and Regioselectivity of the Formation of Hydroxylated Metabolites of NPC1161a

Cytochrome P450 isoforms are the major enzymes involved in drug metabolism, accounting for about 75% of the total number of different metabolic reactions,[36] and have been shown to be responsible for causing oxidative stress and hemotoxicity of 8-AQs,[8] specifically by catalyzing the hydroxylation of 8-AQs. Both experimental observations[37] and theoretical calculations[38] indicated that the barrier to H-abstraction, the first and rate-controlling step in the hydroxylation carried out in CYP450, exhibits a linear correlation with the homolytic C–H bond dissociation energy (HBDE), implying the possibility and reliability of predicting the feasibility and regioselectivity of hydroxylation of substrates in the protein environment of P450. Therefore, in this work, we also calculated the HBDEs of C–H/N–H bonds at all possible positions in NPC1161a to evaluate the feasibility and regioselectivity of forming hydroxylated metabolites in NPC1161a, at the B3LYP/6-311++G**//B3LYP/6-31G** level in the gas phase.

The conformer NPC1161a02 with an extended side chain was selected as a model for the HBDE computations. The calculated HBDEs are listed in Table 5, from which one can see that the methyl group at C4, the methoxyl group at C6, and the positions on the side chain are the easiest to be hydroxylated, specifically in the order of (C4)Me, C2′, C5′, N1′, (C6)OMe, N6′, and C4′ at the B3LYP/6-311++G**//B3LYP/6-31G** level. The positions in the quinoline and phenolic rings are predicted as relatively difficult to be hydroxylated. The molecular orbital analysis of neutral conformer NPC1161a02 (Figure 4) indicated that the electronic density in the highest occupied molecular orbital (HOMO) was localized more at N1′ (and at C5) whereas the electronic density in the lowest unoccupied molecular orbital (LUMO) was localized more at N1 (and at C4), implying that the N1′ atom is more likely than all other atomic centers in the molecule to pair with an electron from the oxygen attached to the iron in the active site of CYP450. The biochemical interaction/reaction between NPC1161 enantiomers and CYP450 isoforms should also be taken into account to interpret the enantiomeric differences in hemotoxicity contributions from the hydroxylation of NPC1161.

Table 5

Calculated Relative C–H/N–H Homolytic Bond Dissociation Energies (HBDEsa) in Conformer NPC1161a02 at All Possible Positions at the B3LYP/6-311++G**//B3LYP/6-31G**++ Level in the Gas Phase

H@b	ΔHBDE (kJ/mol)	H@b	ΔHBDE (kJ/mol)
4Me	0	2′Me	61.89
2′	8.41	2	80.43
5′	12.8	7	101.51
N1′	27.51	3	105.31
N6′	46.65	5″	108.86
6OMe	40.45	2″	109.36
3′	46.82	6″	116.12

Including zero-point energy at the B3LYP/6-31G** level in the gas phase.

Position at which the hydrogen is attached.

Figure 4

Calculated frontier molecular orbitals of neutral conformer NPC1161a02 at the B3LYP/6-31G** level in the gas phase.

Conclusions

Conformational analysis of the very flexible antimalarial candidate NPC1161a showed that an extended conformation of the side chain was found in 45% and 54% of the Boltzmann-weighted structures for the neutral and ionized conformers, respectively, with conformers 02 and 20 predominating, at the B3LYP-SCRF(PCM)/6-311++G**//B3LYP/6-31G** level in water. The calculated IPs of NPC1161a indicated their conformational insensitivity. The metabolite hydroxylated at the N1′ position was suggested as contributing the most to the potential hemotoxicity of NPC1161a, as evidenced by its lowest IP of ∼430 kJ/mol at the B3LYP-SCRF(PCM)/6-311++G**//B3LYP/6-31G** level in water. The calculated homolytic bond dissociation energies demonstrated the preference for hydroxylation on the side chain. The results in this work should apply equally to both enantiomers NPC1161a and NPC1161b, whereas differences in their activities and toxicities could arise from interactions with proteins such as transporters or other metabolizing enzymes. Further theoretical investigations of biochemical interactions and reactions between NPC1161 enantiomers and CYP450 isoforms are currently in progress to elucidate the source of toxicity differences between enantiomers.