Interaction of letrozole and its degradation products with aromatase: chemometric assessment of kinetics and structure-based binding validation

Abstract Letrozole is one of the most prescribed drugs for the treatment of breast cancer in post-menopausal women, and it is endowed with selective peripheral aromatase inhibitory activity. The efficacy of this drug is also a consequence of its long-lasting activity, likely due to its metabolic stability. The reactivity of cyano groups in the letrozole structure could, however, lead to chemical derivatives still endowed with residual biological activity. Herein, the chemical degradation process of the drug was studied by coupling multivariate curve resolution and spectrophotometric methodologies in order to assess a detailed kinetic profile. Three main derivatives were identified after drug exposure to different degradation conditions, consisting of acid-base and oxidative environments and stressing light. Molecular docking confirmed the capability of these compounds to accommodate into the active site of the enzyme, suggesting that the sustained inhibitory activity of letrozole may be at least in part attributed to the degradation compounds.


Introduction
The emergence and progression of hormone-responsive breast cancer largely depend on the endocrine oestrogenic activity. The treatment of this form of cancer is currently pursued following two main approaches. The first one, more effective during premenopausal events when the gonads are the primary producers of oestrogen, is based on the use of antagonists targeting oestrogen receptors. In this case, tamoxifen is the most prescribed drug 1 . A second approach is preferred in post-menopause, when the production of oestrogen is limited to extra-glandular tissues and consists in the administration of aromatase inhibitors. This enzyme, belonging to the CYP450 superfamily, is expressed in several tissues including gonads, placenta, brain, bone, and adipose tissue. Aromatase acts as a monooxygenase, catalysing the conversion of androgen to oestrogen; thus, its inhibition results in a near complete oestrogen deprivation. Several compounds have been proposed as aromatase inhibitors, leading to a first, second, and third generation of specific agents, the latter being the most prescribed nowadays. Fadrozole and vorozole, which show a nonsteroidal structure, were once commonly prescribed but are no longer utilised today, whereas the use of the irreversible steroid inhibitor exemestane is limited to patients pre-treated with tamoxifen in adjuvant therapy 2,3 .
Letrozole (LTZ) and anastrazole, two compounds showing structural similarities, are today the most used reversible non-steroidal inhibitors 4 . The scaffold of both LTZ and anastrazole is based on a 1,2,4-triazole ring substituted in position 1 with different benzylic residues bearing two cyano groups, whose presence seems necessary for the anti-enzymatic activity. These peculiar features suggest a considerable chemical reactivity, which has been studied under different conditions 5,6 . In particular, LTZ is slightly more potent and longer-acting than anastrazole, and has been shown to easily undergo hydrolysis of its cyano groups in both acidic and alkaline medium, resulting in the formation of the corresponding dicarboxylic acid or dicarboxyamide, respectively ( Figure 1). On the other hand, under oxidising condition, the N-oxide in position 2 of the triazole was the main product generated. These three main degradation products have not yet been fully studied with regard to their intrinsic inhibitory activity, which could lead to a prolongation of the drug action in vivo.
Previous studies have shown that LTZ is highly sensitive to extreme pH values and oxidative conditions, whereas it is stable under light irradiation [5][6][7] . Many experimental studies have been focussed on LTZ and its metabolites or degradation products, which have been investigated by various analytical techniques. The chromatographic approach was successfully applied using different detector systems such as UV/Vis, fluorescence, and mass spectroscopy 5,[8][9][10] . Spectroscopic methodologies were also used for assaying the drug in pharmaceutical formulations [11][12][13] .
In this work, we describe in detail the degradation profile of LTZ under different experimental conditions. It is well known that nitrile catabolism can follow two distinct pathways: (1) conversion to carboxylic acid catalysed by nitrilase, and (2) amide formation mediated by nitrile hydratase/amidase. These reactions can be easily reproduced using different hydrolytic chemical conditions 14 . A further factor leading to LTZ degradation is represented by non-hydrolytic oxidative conditions, when an easy formation of the corresponding 1-(bis(4-cyanophenyl)methyl)-1H-1,2,4-triazole-2-oxide (LNO) is observed. The formation of other minor degradation products has been reported 6 , but their presence has not been considered here. This study aimed to define the kinetics of the degradation products of LTZ, monitoring them by UV/Vis spectrophotometry and processing the spectral data by Multivariate Curve Resolution-Alternating Least Squares (MCR-ALS) methodology. This chemometric procedure was used because it is particularly effective in following the chemical transformation processes, allowing to resolve the spectra and concentration profiles of the components involved. An independent HPLC-DAD method was defined to validate the results obtained by the multivariate resolution of the UV kinetic studies. The UV spectra from the HPLC-DAD detector were compared with the spectra predicted by the multivariate procedure, demonstrating a significant overlap of the spectral curves related to the degradation products. Furthermore, docking experiments were performed to verify whether the degradation compounds were able to accommodate in the aromatase active site with binding mode conformation and affinity comparable to the parent drug. To this aim, LTZ and its three main degradation compounds were docked in the crystal structure of aromatase and the results were compared with those obtained for androstenedione (ASD), the endogenous ligand of the enzyme.

Chemicals
All samples were prepared using distilled water and analytical grade reagents. LTZ standard was purchased from Merck (Milan, Italy), with certified purity >98%, and was used as received. Hydrochloric acid, sodium hydroxide, hydrogen peroxide, HPLC grade water and acetonitrile were also supplied by Merck (Milan, Italy).

Instruments
A Crison pH-meter GPL 22 (Barcelona, Spain) was used to measure the pH values of the samples during degradation. Photodegradation experiments were carried out in a chamber Suntest CPSþ (Atlas, Milan, Italy), equipped with a xenon lamp, using light irradiation according to the ID65 standard of the ICH rules (International Council for Harmonisation of Technical Requirements for Pharmaceuticals for Human Use) 15 . Spectrophotometric measurements were recorded using an Agilent 8453 Diode Array spectrophotometer (Agilent Technologies, CA, USA). HPLC analysis was performed by using an Agilent 1100 series chromatograph (Agilent Technologies, CA, USA), equipped with a binary pump delivery system and a diode array UV-Vis detector.
Spectrophotometer ChemStation software (Agilent Technologies, CA, USA) was used for acquiring experiment data and converting raw UV spectra files (.sd) to human-readable files (.csv), suitable to be imported directly into MATLABV R (The MathWorks, Inc., MA, USA). All chemometric analyses were performed under MATLAB V R computing environment, where the MCR-ALS procedure (GUI version 2.0) has been implemented 16 .

Experimental procedures
A standard stock solution of LTZ (1 mg/mL) in ethanol was properly diluted to obtain seven samples (10 mg/mL) at pH values close to 1, 3, 5, 7, 9, 10, 12 by addition of HCl or NaOH. The pH adjustment was verified not to cause any significant dilution. These samples, transferred in perfectly stoppered quartz cells, were monitored by UV/Vis spectrophotometric analysis in the wavelength range 220-400 nm, starting just after sample preparation (t ¼ 0 min) and every 10 min, up to 24 h.
A second sample set was prepared by combining different pH values and oxidative conditions. Appropriate aliquots of 30% (v/v) H 2 O 2 , HCl and NaOH were added to the LTZ stock solution to obtain three LTZ samples with concentration 10 mg/mL containing 0.5% (v/v) of H 2 O 2 at pH values close to 2, 7 and 12. These samples were monitored by UV/Vis analysis (200-400 nm) every 10 min, up to 24 h.
Photolability of LTZ was tested by exposing the LTZ sample (10 mg/mL) to a xenon lamp under the following conditions: light radiation between 300 and 800 nm, power 550 W/m 2 (33 kJ/(min x m 2 ), constant temperature of 25 C 15 . LTZ behaviour to light exposure was monitored by recording UV spectra at time t ¼ 0 min and every 10 min along the photodegradation experiments, up to 12 h.
HPLC analysis was performed on a reverse phase Gemini LC column (250 Â 4.60 mm, 5 lm C18, Phenomenex, Torrance, CA) using phosphate buffer (pH 5.8) and acetonitrile in the ratio 80:20 (v/v) as mobile phase, pre-filtered through a 0.45 mm filter. The injection volume was 20 mL and the flow rate of the mobile phase was 1.0 ml/min at room temperature of 20 C. HPLC was performed on LTZ ethanol solutions (1 mg/mL) exposed to hydrolytic degradation with HCl and NaOH and oxidation with H 2 O 2 immediately after preparation (t ¼ 0 min) and after 24 h.

Multivariate curve resolutionalternating least squares (MCR-ALS)
Multivariate Curve Resolution (MCR) defines a family of methodologies aimed at resolving the chemical contributions to the outcome of an experiment, described through a data matrix. They have been applied to study different types of multivariate or multicomponent chemical systems 17 . MCR method decomposes the experimental data matrix (D) into a reduced set of contributions of chemical species (in our study, LTZ and its degradation products), using a bilinear model obtained rewriting the Lambert-Beer's law in multiway mode 17,18 : where D (n,m) is the experimental data matrix, C (n,k) is the matrix that includes the concentration evolution of the k components, S T (k,m) represents the pure spectra of k species, and E (n,m) is the unexplained variance in the model.
The number of components involved in the matrix D (chemical rank) can be estimated by Principal Component Analysis (PCA) or Singular Value Decomposition (SVD) algorithms 19 . The chemical rank assumes that the species contributing to the measured spectra give singular values larger than the other signal contributions, such as experimental or instrumental noise. A rank deficiency may occur in some cases, leading PCA to identify a chemical rank lower than the actual value, but multi-experiments and the subsequent matrix augmentation obtained by processing all experiments can be used to avoid this drawback 16 . MCR allows simultaneous processing of data from different origins, via rowwise, column-wise, or both row-and column-wise data matrix augmentation.
Once the number of components is known, the ALS iterative algorithm uses a preliminary estimation of S T or C matrices, and a series of constraints to optimise the MCR model. The application of constraints such as non-negativity, unimodality, and concentration closure allows one to optimise the results according to a chemical meaning.
The quality and reliability of the multivariate resolution can be assessed using the explained variance (%R 2 ) and the lack of fit (%) 20 . These figures of merit allow to evaluate the dissimilarity between the experimental matrix D and the data modelled by MCR-ALS elaboration. The equations are listed below: where d ij is an element of the matrix D and e ij is the related residual value obtained from the difference between the matrix D and the calculated CS T matrix. In the present work, two different MCR-ALS techniques were applied to the datasets. The standard soft-MCR-ALS (S-MCR-ALS) was applied for the preliminary analysis of all data matrices and the hard and soft MCR-ALS (HS-MCR-ALS) was used to analyse in depth the degradation kinetics.

Molecular docking
Molecular docking was performed starting from a crystallographic structure of aromatase. The Protein Data Bank (PDB) contains 38 structures of aromatase and, among them, 27 belong to Homo sapiens. The entry selected was the one recently deposited with PDB ID: 5JKV 21 , and it was chosen over the others because it has a more complete amino acid sequence and contains the endogenous ligand ASD in the binding site. It also contains in a secondary allosteric binding sites a molecule of pentaethylene glycol (PEG), which acts as a weak inhibitor of aromatase. Both ASD and PEG were removed from the protein, whereas the haem cofactor was left in its position. The molecular structures of the ligands (LTZ, DBA, DBO, LNO; see again Figure 1) were built by using the modelling software Avogadro 22 .
Docking calculations were performed by using AutoDock Vina 1.1.2 23 . Preliminary conversion of the structures from the PDB format was carried out by using the graphical interface AutoDock Tools 1.5.6 24 . During the conversion, polar hydrogens were added to the crystallographic enzyme structures, whereas apolar hydrogens of the ligands were merged to the carbon atom they are attached to. Full flexibility was guaranteed for the ligands, resulting in 3, 5, 7 and 3 rotatable dihedral angles for LTZ, DBA, DBO and LNO, respectively. A single simulation run was carried out in each case at very high exhaustiveness, 16 times larger than the default value 25,26 . The binding modes of the ligands were analysed through visual inspection, and intermolecular interactions were evaluated by using the automated protein-ligand interaction profiler PLIP 27 .
Initially, S-MCR-ALS analysis was performed for each experiment by defining two augmented matrices from D pHi and D oxj subsets and the single matrix D k . Augmented were built by merging data subsets in column-wise mode, considering the data obtained in different experiments monitored by the same analytical procedure. The rank analysis and soft modelling of the UV data resolved four components in the acid-base degradation process, and three components in the oxidative degradation. A fair stability of the drug in solution was observed when it was exposed to light for up to 360 min, as evident in Figure 2, because the spectral signal did not undergo significant changes.
The rank analysis of matrix D pHaug estimated 4 species distributed within all subsets. C pHaug and S pHaug T were optimised by the S-MCR-ALS algorithm with satisfying data fitting in term of %LOF ¼ 1.45% and R 2 ¼ 99.9%. The band boundaries of feasible S-MCR solutions were evaluated by using the MCR-BANDS routine contained in MCR-GUI 18 . The band boundaries of feasible S-MCR-ALS solutions evaluated for LTZ acid-base degradation experiment are shown in Figure 3. The rotational ambiguities were evident in matrix C and S T for all resolved species. The matrices C pHaug and S pHaug T described the concentration evolution during the degradation process (time direction) and the spectral features correlated to the species involved in the experiments arranged in D pHaug . In the acid-base degradation process, two species were resolved for LTZ at time 0, depending on the different pH. LTZ is characterised by a pK a ¼ 2.17 and two forms, the ionised LTZ A and the non-ionised LTZ B , were detected with different spectral shape in S pH T (Figure 3). C pHaug resolved the component corresponding to the first degradation product at acidic pH values, which formed due to the partial hydrolysis of the cyano moiety to amide DBA (C pH1 ). The hydrolysis process was pH-dependent and occurred only at extreme acidity values. The drug was moderately stable close to neutral pH conditions, as no significant spectral variations were observed at pH 5, 7 and 9. In contrast, a complete hydrolysis from cyano to carboxylic groups was obtained at very high pH values.
In the experiment at pH 10, only the signal corresponding to DBA was resolved, whereas at pH 12 both the amide and carboxylic hydrolysis products DBA and DBO were detected (C pH12 ). The acid-base degradation process of LTZ was found to follow two sequential first-order reactions at different pH values, and the proposed reaction process is LTZ > DBA > DBO (see again Figure 3).
The matrices C oxaug and S oxaug T described a different behaviour of LTZ when exposed to an oxidative environment. MCR-ALS analysis of the matrix D oxaug (% LOF ¼ 0.87% and R 2 ¼ 99.9%) resolved three components, which encompassed both forms of LTZ (LTZ A and LTZ B ) and the oxidation product was characterised by the N-oxide at position 2 in the triazole of LNO, obtained from a first-order reaction LTZ > LNO (Figure 4). The different evolution of the chemical species in the various experiments highlighted the dependence of the oxidation process on the pH.
In some cases, MCR processing can be limited by the failure to guarantee a unique result, as the technique suffers from some intensity and/or rotational ambiguity (Figures 3 and 4). In the case of kinetic studies, the scale of C or S T is known, and the total mass balance within the system is constant. In a closed reaction kinetic system, intensity ambiguity can be solved by a normalisation procedure or applying the closure constraint. In case of rotation ambiguities, the adoption of more constraints to MCR modelling can be useful in reducing them. Aside from the natural constraints, such as no-negative, unimodality and closure constraints, the most effective tools for reducing rotation ambiguities are the use of local rank and/or selective constraints, the simultaneous analysis of multi-datasets and the use of hard modelling routines 18 .
In a subsequent step, matrices D pHaug and D oxaug were processed by HS-MCR-ALS, which was used to analyse in depth the kinetic degradations. HS-MCR-ALS uses a non-linear curve fitting routines to iteratively fit the concentration profiles resolved by MCR-ALS with the optimal kinetic rate, leading to an improved resolution of the concentration profiles of the species formed during the acid-base and oxidation experiments. It also reduces the rotational ambiguities associated with the S-MCR-ALS results, allowing to define chemical models explaining the data variance in the matrices and to evaluate the rate constants 28 . HS-MCR-ALS results are shown in Table 1 and Figures 3 and 4, Figure S1 (Supplementary Material).
The comparison of the estimated rate constants allowed us the evaluation of the influence of pH on LTZ hydrolysis and oxidation. The drug was more affected by basic rather than acidic conditions and it was stable in an environment close to neutrality. In a basic environment, the drug degraded quickly, reaching high k values at pH 12 (k pH12a ¼ 5.63 Â 1 0 À5 s À1 ) with a complete hydrolysis. Rapid degradation in oxidative environment was also observed, with a rate constant k ox7 of about 1 0 À6 s À1 . Under mixed conditions, when the oxidative environment was associated with either acidic or basic pH, the stability of LTZ was reduced when the oxidative environment was coupled to a pH of 12 (k ox12 ¼ 3.17 Â 1 0 À5 s À1 ), whereas it increased at acidic pH (k ox1 ¼ 4.0 Â 1 0 À8 s À1 ).
The MCR-ALS results obtained in the analysis of the timedependent UV data were confirmed by chromatography investigation. HPLC-DAD analysis of the samples before and after the degradation experiments confirmed the formation of the same degradation products elaborated by the MCR-ALS procedures.  HPLC analysis of acid-base degradation samples provided two additional peaks of DBA (RT 6.22 min) and DBO (RT 5.81 min) in addition to the LTZ signal (RT 9.12 min), as well as a signal for the LNO product was detected during oxidative degradation (RT 7.34 min). The peaks processed by DAD had the same spectral profiles as the products resolved in the previous spectrophotometric study (see Figure S2 and S3 in Supplementary Material).

Docking simulations
Aromatase is a monomeric peripheral membrane protein with 503 amino acid residues, organised in 12 a-helices joined by coil regions and surrounding a haem group, as shown in Figure 5. The active site (see inset of Figure 5) lies next to a channel through which water and oxygen molecules involved in the catalysis, as well as substrates, can enter from the outer face of the membrane. At the channel entry, the couple of residues Arg192/Glu478 acts as a gate-system, whereas the inner portion of the channel ending with the binding site is characterised by the presence of both polar and ionic residues (Asp309, Thr310 and Ser478), and hydrophobic side chains (Phe221, Val313 and His480). The substrate binds at the androgen-specific haem-distal pocket, interacting with the key residues Arg115, Met374, Thr310 and Asp309. In the case of interaction with ASD, the side chains of the first two amino acids form polar contacts with the 17-carbonyl oxygen of the ligand, whereas the last two are involved in the binding with the 3-carbonyl oxygen. In particular, the protonated side chain of Asp309 seems to be critical for substrate binding orientation and catalysis 29 .
The side chain of Arg115, Ile133, Phe134, Phe221, Trp224, Ala306, Thr310, Val370, Val373, Met374 and Leu477 together with the haem porphyrin portion form a haem-distal catalytic cavity in which the ligands, including ASD, easily accommodate, interacting by means of van der Waals forces. In the aromatase structure there is also a haem-proximal cavity, located on the opposite side of the haem group with respect to the main binding site, which could host a second ligand. In the molecular complex we used as a model (PDB ID: 5JKV), this site is occupied by PEG, a pentameric fragment of polyethylenglicole, which is a weak inhibitor of aromatase able to participate in the modulation of the enzymatic activity of the protein ( Figure 5). A third site acts as an access channel and connects the haem-distal binding site with the endoplasmic reticulum via the outer surface of the protein, allowing the entry and exit of steroids. Multiple additional hotspots potentially suitable to accommodate ligands have been identified within the enzyme protein, but their significance towards the protein catalytic activity remains uncertain 21 .
A quantitative structure-activity relationship (QSAR) study 30 predicted that the reversible aromatase inhibitors LTZ should bind the active site by forming interactions through the two cyano groups and the triazole moiety. In our docking experiments, LTZ was observed to anchor to aromatase in the main protein binding site, with a binding energy of À7.9 kcal/mol. This value corresponds to a dissociation constant in the low micromolar range, which indicates a good affinity. Furthermore, the ligand occupies the same location of ASD, as shown in Figure 6. The N 4 atom of the 1,2,4-triazole coordinates the Fe 2þ of the haem prosthetic group, whereas the two cyano groups of the LTZ act as a hydrogen bond acceptor by interacting with Met374 and Ser478 residues, similarly as the 3-and 17-carbonyl oxygens of the physiological ligand ASD. A further hydrogen bond forms between N 2 atom of the triazole and Thr310. The diphenylmethane scaffold of LTZ mimics the backbone of steroids, preserving for this compound a geometry similar to that of endogenous ligands inside the binding pocket.
Furthermore, we investigated whether the degradation compounds, deriving mainly from the transformation of the cyano groups of LTZ, could interact with aromatase in the same active site as the parent compound, and could consequently participate to the overall biological activity. Thus, LTZ degradation derivatives DBA, DBO and LNO were docked into the main catalytic site (Figure 7). In particular, as detailed in Table 2, DBO forms the same interactions as LTZ, including hydrogen bonds with the protein residues Thr310, Met374 and Ser478, whereas DBA lacks the interaction with Met374. On the other hand, a fourth hydrogen bond was observed between Arg115 and a cyano group of LNO. Hydrophobic interactions with several key residues of the site were also observed in all cases. The binding energy values observed range from À8.2 and À8.3 kcal/mol, which are slightly more favourable compared to that of the parent compound LTZ.
In summary, the simulation findings concur to predict the association of LTZ and its three main degradation products (DBA, DBO, and LNO) within the binding site of aromatase. In spite of some differences in their structural details, including charge delocalisation and ability to form hydrogen bonds, all these compounds showed clear similarities in the docking affinity and anchoring modes in the protein binding pocket. These features were evidenced in blind docking experiments, and without the need to introduce a bias to help accommodating the ligands in the protein binding site. Furthermore, the binding of these compounds is consistent with pharmacophore-based QSAR features dictated by the active site of aromatase, which may be considered an independent confirmation of our docking results. Finally, overlap in the binding modes between our compounds and the physiological ligand ASD, in spite of differences in their structure (LTZ and the degradation products assume a triangular pyramid shape, whereas ASD is a steroid with a more linear arrangement), further support the reliability of the results obtained.

Conclusions
Degradation studies indicated that LTZ, although stable to light irradiation, generates three main derivatives under stressing chemical conditions, such as changes in pH or oxidising environment. The analytical data collected by UV/Vis spectrophotometry during the kinetic experiments were processed by chemometric methodologies. The soft multivariate curve resolution of the data matrices confirmed the role of pH and oxidants in favouring the formation of the diamide (DBA), diacid (DBO) and N 2 -oxide (LNO) degradation products. The estimation of the rate constants of the  degradation processes by HS-MCR showed that the parent drug is more sensitive to basic pH rather than acid or neutral pH, leading to a faster degradation and a complete transformation into the hydrolysis product DBO at pH 12 or LNO in an oxidative environment, respectively. Molecular docking experiments demonstrated that LTZ and its degradation compounds have the ability to adapt to the aromatase active site by interacting with the same key amino acid residues involved in the binding of physiological ligand ASD. The results obtained show that the metabolism of LTZ could lead to compounds still effective, suggesting that this could prolong the drug activity.