Structural Impact Assessment of Cytochrome P450 2A13 Polymorphisms Using Molecular Dynamics Simulations

Abstract

One of the members of CYP, a monooxygenase, CYP2A13 is involved in the metabolism of nicotine, coumarin, and tobacco-specific nitrosamine. Genetic polymorphisms have been identified in CYP2A13, with reported loss or reduction in enzymatic activity in CYP2A13 allelic variants. This study aimed to unravel the mechanism underlying the diminished enzymatic activity of CYP2A13 variants by investigating their three-dimensional structures through molecular dynamics (MD) simulations. For each variant, MD simulations of 1000 ns were performed, and the obtained results were compared with those of the wild type. The findings indicated alterations in the interaction with heme in CYP2A13.4, .6, .8, and .9. In the case of CYP2A13.5, observable effects on the helix structure related to the interaction with the redox partner were identified. These conformational changes were sufficient to cause a decrease in enzyme activity in the variants. Our findings provide valuable insights into the molecular mechanisms associated with the diminished activity in the CYP2A13 polymorphisms.

INTRODUCTION

CYP2A13 is a highly expressed enzyme in the lungs and mucosa and is involved in the metabolism of various chemicals such as nicotine, coumarin, and tobacco-specific nitrosamines.^1–4) Additionally, CYP2A13 contributes to genotoxicity by metabolically activating aflatoxin and 4-(methylnitrosamino)-1-(3-pyridyl)-1-butanone (NNK), which are typical carcinogens found in tobacco. The enzyme also holds significance in tumorigenesis in human respiratory organs. Immunohistochemical analysis revealed the expression of CYP2A13 in pancreatic alpha islet cells, suggesting its involvement as a molecular factor contributing to the significantly higher risk of pancreatic cancer in current smokers.^5,6) Furthermore, CYP2A13 has been reported to be expressed in various other tissues, including the liver, testis, breasts, uterus, ovary, and coronary endothelial cells, with a notable presence in the bladder.^7,8)

CYP2A13 comprises 494 amino acid residues and shares a 93.5% similarity with CYP2A6.⁹⁾ They encompass helices A–L and β1–4, featuring six substrate-recognition sites (SRS1–6).^10,11) SRS-1 consists of B′ and the flanking area (residues 102–112), SRS-2 includes the C-terminus of helix F (residues 197–204), SRS-3 consists of the N-terminus of helix G (residues 232–239), SRS-4 consists of the N-terminus of helix I (residues 288–296), SRS-5 comprises the β2 area (residues 358–368), and SRS-6 contains the central region of β5 (residues 468–476). Due to the similarity between CYP2A6 and CYP2A13, it is presumed that the substrates of CYP2A6 can be metabolized by CYP2A13. Despite their similar three-dimensional structures, the active site of CYP2A13 is larger than that of CYP2A6, even under identical substrate binding conditions.^9,10) The active site CYP2A13 includes Phe107, Phe111, Ala117, Phe118, Phe209, Leu296, Asn297, Phe300, Ala301, Glu304, Thr305, Met365, Leu366, Leu370, and Phe480. Furthermore, the key residues of the active site are Ala117, Ser208, Phe300, Ala301, Met365, and Leu366.⁹⁾ In contrast, those of CYP2A6 are Val117, Ile208, Ile300, Gly301, Val365, and Ile366. CYP2A13 hydroxylates testosterone,¹²⁾ and its catalytic efficiencies for nicotine, cotinine, and NNK are one or two orders of magnitude higher than those of CYP2A6.¹³⁾ The redox partner binding to CYP is required for electron transfer and is predicted to involve the helices C, K, and L of CYP.¹⁴⁾

CYP polymorphisms are well known to affect the metabolism of drugs and pollutants. Genetic polymorphisms of CYP2A13 have been identified in Japanese and French populations.^15–17) The kinetic parameters of the CYP2A13 polymorphisms were characterized using nicotine and coumarin as substrates¹⁸⁾ (Table 1). CYP2A13.7 has no enzymatic activity due to an immature stop codon (Arg101-Stop). All polymorphic mutations resulted in reduced catalytic activities, and in the case of CYP2A13.4, even a single residue mutation affected the activity to such an extent that it could not be determined. There were no significant differences in K_m values, particularly for nicotine, but relatively large effects on catalytic activities were observed.¹⁸⁾ These polymorphisms may contribute to individual variations in cancer development. However, because the experimental structure of the CYP2A13 variant has not been determined, the detailed mechanism of reduced metabolic activity remains unclear.

Table 1. Kinetic Parameters of CYP2A13 Variants for Nicotine C–Oxidation and Coumarin 7-Hydroxylation¹⁸⁾

Protein	Amino acid change	Nicotine CLint/μL min−1	Coumarine CLint/nL min−1
CYP2A13.1		6.86 ± 0.10	183 ± 4.91
CYP2A13.4	R101Q	N.D.	N.D.
CYP2A13.5	F453Y	0.93 ± 0.05*** (14%)	33.2 ± 3.28*** (18%)
CYP2A13.6	R494C	1.43 ± 0.14*** (21%)	94.3 ± 8.37*** (52%)
CYP2A13.8	D158E	1.78 ± 0.76* (26%)	34.7 ± 3.66*** (19%)
CYP2A13.9	V323L	0.91 ± 0.25** (13%)	137 ± 10.2*** (75%)

The percentages of intrinsic clearances (CL_int) were compared with CYP2A13.1 CL_int. * p < 0.05, ** p < 0.01, and *** p < 0.005 compared with CYP2A13.1  N.D.: Not Determined.

We conducted molecular dynamics (MD) simulations to explore the structural changes in CYP2A13 variants (CYP2A13.4, CYP2A13.5, CYP2A13.6, CYP2A13.8, and CYP2A13.9). In this study, variants with single mutations were selected to understand the detailed mechanism by which only one mutation reduces the metabolic activities of CYP2A13. The allele frequencies of CYP2A13.4, CYP2A13.5, and CYP2A13.6 were 0.3, 0.3, and 1.0%, respectively, in the Japanese population.¹⁵⁾ CYP2A13.8 and CYP2A13.9 were each present at 1% in French individuals.¹⁶⁾ Three-dimensional (3D) structural information on such variants is crucial for understanding the mechanism of reduced activity in polymorphisms. MD simulations are effective methods for the structural characterization of variant proteins, including CYPs.^19–23) Therefore, MD simulations for the CYP2A13 variants were performed for 3D structural investigation. A simulation of the wild type (CYP2A13.1) was also conducted for comparison.

MATERIALS AND METHODS

Computational Methods

The experimental structure of CYP2A13 was retrieved from the protein data bank (PDB) to construct the initial structure. A crystal structure with high resolution was obtained (PDB ID: 4EJI),¹⁰⁾ and all ligands and water molecules were deleted. A covalent bond was established between the vertical ligand Cys439 and heme iron using the tleap module in AmberTools. To neutralize the charge of the system, chloride ions were placed as counter ions, and the system was solvated with TIP3P water molecules,²⁴⁾ spaced at least 8 Å from the protein surface. Using this as the initial structure, optimization calculation for the water and overall structure was conducted. The minimization for water molecules and counter ions was performed for 1000 steps, and overall structural minimization was performed for 2500 steps. Subsequently, a temperature-increasing simulation of 20 ps was conducted on the optimized structure, raising the temperature from 0 to 300 K. Following this, an equilibrating MD simulation was performed in the NPT ensemble (isothermal–isobaric ensemble) for 2000 ns, with a time step of 2 fs. Amber18 was used for all calculations,²⁵⁾ and the AMBER ff14SB force field was applied for the amino acid parameters.²⁶⁾ The force field parameter developed by Oda et al. was used around the heme iron,²⁷⁾ and that determined by Giammona was used for the other region of the heme molecule.²⁸⁾ The cutoff distance for nonbonding interactions was set to 10 Å, and electrostatic interactions were calculated using the particle mesh Ewald method under periodic boundary conditions.²⁹⁾ The SHAKE algorithm was applied to constrain the lengths of covalent bonds containing hydrogen atoms.³⁰⁾ For the simulations of the CYP2A13 variants, initial structures were constructed by introducing mutations into the wild-type structure after the MD simulation. The side chains of the mutated residues were built using tleap module. For R494C variant, cysteine (Cys) was introduced with a reduced free sulfide. Although the calculation conditions for the variants were the same as those for the wild type, the equilibrium simulation time was 1000 ns. To validate the convergence of the simulations, root mean square deviations (RMSDs) for main-chain atoms were evaluated along the calculated MD trajectories using the 3D structures after optimization as the reference structures. Structural flexibilities were assessed by the root mean square fluctuations (RMSFs) calculated for all Cα atoms in the last 10-ns MD trajectories. For hydrogen bonding criteria, the distance between the two heavy atoms and the angles between the acceptor, hydrogen, and donor atoms were determined using cutoff values of 3 Å and 135°, respectively. Five thousand frames extracted every 2 ps of trajectories during the final 10 ns were used for hydrogen bonding analysis.

RESULTS AND DISCUSSION

The RMSD plots are shown in Fig. 1, and the values converged for all simulations. The RMSD between the experimental structure and the final structure of the simulation for the wild type was 1.322 Å, indicating almost no structural changes. This indicated that no artificial structural changes occurred. The overall structures of all variants were similar to that of the wild type (Fig. 2). Conversely, the RMSF plots showed differences between the wild type and some variants (Fig. 3). There were increases by more than 1 in the RMSF values of the DE loop in CYP2A13.5 and the K”L loop in CYP2A13.8. There was no significant difference in RMSF values between the wild type and CYP2A13.4, CYP2A13.6, and CYP2A13.9. The relationship between these changes and activity will be discussed in each section.

Fig. 1. Root Mean Square Deviation (RMSD) Plots of the Main-Chain Atoms for CYP2A13 Wild Type and Variants

(A) wild type, (B) CYP2A13.4, (C) CYP2A13.5, (D) CYP2A13.6, (E) CYP2A13.8, and (F) CYP2A13.9.

Fig. 2. Overall Structures of the CYP2A13 Wild Type and Variants Obtained by the Final of the Simulations

These figures are superimposed on the wild type (green) and (A) CYP2A13.4 (cyan), (B) CYP2A13.5 (magenta), (C) CYP2A13.6 (orange), (D) CYP2A13.8 (pink), and (E) CYP2A13.9 (light blue).

Fig. 3. Root Mean Square Fluctuation (RMSF) Plots of Cα Atoms in the Last 10 ns of Simulations

(A) wild type, (B) CYP2A13.4, (C) CYP2A13.5, (D) CYP2A13.6, (E) CYP2A13.8, and (F) CYP2A13.9. Plots of differences in RMSF values between the wild type and variants for (G) CYP2A13.4, (H) CYP2A13.5, (I) CYP2A13.6, (J) CYP2A13.8, and (K) CYP2A13.9.

Mutation in CYP2A13.4 Attenuates Interactions with Heme

To investigate the mechanism of enzymatic activity disappearance in CYP2A13.4, the calculated results of this variant were compared with those of the wild type. In this variant, a glutamine (Gln) residue substitutes Arg101 interacting with the carboxyl oxygen of heme. The guanidino group of Arg101 interacted with the carboxyl group of heme in most of the trajectories in the last 10 ns of the simulation (Table 2, Fig. 4A). In contrast, the hydrogen-bond formation rate between Gln101 and heme carboxyl oxygen was only 0.26%. Instead, Arg437 interacted with heme in CYP2A13.4 (Fig. 4B), whereas its frequency was significantly lower than that of Arg101 in the wild type. These results suggest that the interactions with heme are reduced in CYP2A13.4, leading to a decrease in enzymatic activity due to unstable heme binding to the protein. Alternatively, in the wild type, Arg101 formed hydrogen bonds with Ala117, Ala317, and Arg373. Gln101 in CYP2A13.4 exhibited reduced hydrogen bonds with Ala117 and Arg373 and no hydrogen bond with Ala317. Since Ala371 and Arg373 are located on the KK' including Leu370, one of the residues forming the active site.⁹⁾ Therefore, those changes of the hydrogen-bond formation related to the KK' loop are also likely responsible for the decrease in enzymatic activity.

Table 2. Comparison of Hydrogen-Bond Formation between the Wild Type and CYP2A13.4 in the All Trajectories for the Last 10-ns Simulations

Acceptor	Donor	Wild type (%)	2A13.4 (%)
Heme O2A	Arg101 Nη2	100.0	—
Heme O2D	Arg101 Nη2	100.0	—
Heme O2A	Arg101 Nη1	89.98	—
Heme O2D	Arg101 Nε	50.86	—
Heme O1A	Gln101 Nε2	—	0.260
Heme O1A	Arg437 Nη1	0.000	64.14
Heme O1A	Arg437 Nε	0.000	54.48
Heme O2D	Arg437 Nη1	0.000	22.04
Ala117 O	Arg101 Nε	95.92	—
Ala117 O	Gln101 Nε2	—	60.66
Ala371 O	Arg101 Nη2	100.0	—
Arg101/Gln101 O	Arg373 Nη2	48.56	23.20
Gln101 Oε1	Arg373 Nη2	—	23.46

Fig. 4. Structural Comparison between the Wild Type and CYP2A13.4

The (A) wild type and (B) CYP2A13.4 are shown in green and cyan, respectively. Nitrogen, oxygen, and hydrogen are represented in blue, red, and white, respectively, in the stick model. Iron is shown as an orange sphere using a model written as the van der Waals radius. The dotted lines indicate the ionic and hydrogen bonds.

Mutation in CYP2A13.5 Affects Multiple Helices

To investigate the mechanism of enzymatic activity reduction in CYP2A13.5, we compared the calculated structures of the F453Y variant with those of the wild type. Phe453, located on helix L and surrounded by hydrophobic residues (Phe314, Val327, Val355, and Ile356), undergoes substitution with tyrosine, potentially destabilizing hydrophobic cluster formation. In the calculated structure, the orientation of the Tyr453 side chain in CYP2A13.5 significantly deviated from that of Phe453 in the wild type to avoid the hydrophilic effect of the hydroxyl group (Fig. 5). This change allowed the hydroxyl group of the Tyr453 side chain to form hydrogen bonds with Leu310 and Phe314 (Table 3). There was no change in the frequency of hydrogen-bond formation with Met457 and Leu449, which are involved in Helix L formation. Consequently, the F453Y mutation was deemed to impact the assembly of helix K, L, I, and J by altering hydrophobic clusters rather than affecting helix formation. Because helix K is important for interaction with the redox partner,¹⁴⁾ a change in its position causes a decrease in activity. Furthermore, the flexibility in the DE loop were observed in CYP2A13.5 (Fig. 3). Helices D and E were in contact with Helix I; therefore, the structural change due to the subtle shift of Helix I was considered to be transmitted. The maximum distance between the Cα of the Helix I residues in the final structure was 4.7 Å (Gly302). From there, the C-terminal three residues were displaced by about 3 Å. Such subtle changes in each helix could also affect the enzymatic activity of CYP2A13.5.

Fig. 5. Structural Comparison between the Wild Type and CYP2A13.5

The (A) wild type and (B) CYP2A13.5 are presented in green and magenta, respectively. Nitrogen, oxygen, and hydrogen are displayed in blue, red, and white, respectively, in the stick model. Iron is shown as an orange sphere using a model written as the van der Waals radius. The dotted lines indicate the ionic and hydrogen bonds.

Table 3. Comparison of Hydrogen-Bond Formation between the Wild Type and CYP2A13.5 in the All Trajectories for the Last 10-ns Simulations

Acceptor	Donor	Wild type (%)	2A13.5 (%)
Phe453/Tyr453 O	Met457 N	97.80	99.98
Phe453/Tyr453 O	Ile456 N	23.24	18.32
Leu449 O	Phe453/Tyr453 N	84.50	99.44
Phe450 O	Phe453/Tyr453 N	3.620	1.060
Leu310 O	Tyr453 Oη	—	57.28
Phe314 N	Tyr453 Oη	—	10.28

Mutation in CYP2A13.6 Affected the Active Site Structure

To investigate the mechanism of enzymatic activity reduction in CYP2A13.6, we examined structural changes in the R494C variant. Arg494, the C-terminal residue, undergoes a single mutation that reduces V_max by less than half.¹⁸⁾ In the wild type, Arg494 interacted with Asn459, Asp332, Met457, and Gln458 (Table 4). Although Cys494 formed a hydrogen bond with Asn459, no interaction with Asp332 was observed in CYP2A13.6. These changes in hydrogen-bond formation induce a shift in the position of the main chain around Cys494 (the distance of Cα between Cys494 and Arg494 was 2.0 Å in the final structures), resulted in increased hydrogen-bond formation with Lys337 and Arg461 (Figs. 6A, B). The structural changes at Arg461, situated on β3-1, lead to structural changes involving the entire C-terminal loop, including β4. Consequently, hydrogen bonds between Phe480 and Thr212, undetected in the wild type, were detected in 92.04% of the trajectories of CYP2A13.6. These findings indicate that the position of Phe480, a crucial component of the active site, was closer to Helixes F and F′ (Figs. 6C, D). Therefore, the R494C mutation reduces enzymatic activity by causing structural changes in the active site.

Table 4. Comparison of Hydrogen-Bond Formation between the Wild Type and CYP2A13.6 in the All Trajectories for the Last 10-ns Simulations

Acceptor	Donor	Wild type (%)	2A13.6 (%)
Asn459 O	Arg494/Cys494 N	94.66	90.52
Asp332 Oδ1	Arg494 Nη1	85.84	—
Asp332 Oδ2	Arg494 Nη1	79.26	—
Asp332 Oδ2	Arg494 Nη2	77.82	—
Asp332 Oδ1	Arg494 Nη2	13.82	—
Met457 O	Arg494 Nη2	66.88	—
Gln458 O	Arg494 Nη2	25.20	—
Arg494/Cys494 O	Lys337 Nζ	37.50	69.96
Arg494/Cys494 Oxt	Lys337 Nζ	29.76	25.24
Arg494/Cys494 O	Arg461 Nε	0.000	21.02
Arg494/Cys494 O	Arg461 Nη2	0.240	29.90
Thr212 O	Phe480 N	0.000	92.04
Gln218 Oε1	Phe480 N	—	15.52

Fig. 6. Structural Comparison between the Wild Type and CYP2A13.6

The (A) wild type and (B) CYP2A13.6 are shown in green and orange, respectively. Nitrogen, oxygen, hydrogen, and sulfur are visualized in blue, red, white, and yellow, respectively, in the stick model. Iron is shown as an orange sphere using a model written as the van der Waals radius. The dotted lines indicate the ionic and hydrogen bonds.

Mutation in CYP2A13.8 Affects Interaction with Heme

To investigate the mechanism of enzymatic activity reduction in CYP2A13.8, we compared the calculated structures of the D158E variant with those of the wild type. Despite the substitution of an aspartic acid residue with a glutamic acid (Glu) residue, the enzymatic activity of the variant diminished to approximately 20% of that observed in the wild type.¹⁷⁾ D158E mutation did not affect hydrogen-bond formation related to the main-chain atoms, whereas the interactions with Arg161 were increased in the variant (Table 5). The side chain of Phe155 to shifted toward Glu158 because a space was caused by the shift of the negatively charged side chain of Glu158 to interact with Arg161 (Figs. 7A, B). The distance of Cε of Phe155 was 6.5 Å in the final structures. Glu151 interacted with Arg148 in the wild type, whereas in CYP2A13.8, Glu151 also shifted toward Phe155 (the distance of Oε₁ of Glu151 was 4.2 Å in the final structures), and Glu152 interacted with Arg148. These alterations in the orientation of the Glu and arginine (Arg) side chains were accompanied by a significant change in the orientation of Arg143. While Arg143 and Glu146 frequently interacted in the wild type, such interactions were rarely observed in CYP2A13.8; instead, Glu146 and Arg446 interacted frequently. In the wild type, Arg446 formed a hydrogen bond with Leu444 whose side chain interacts hydrophobically with heme. The changes in interaction involving these amino acid residues located on Helix L were further supported by the increase in the RMSF values of the K”L loop (Fig. 3). Therefore, the impact of the mutation on the 3D structure propagated like a domino effect in CYP2A13.8, and the changes in interaction with heme may be responsible for the reduction in enzymatic activity.

Table 5. Comparison of Hydrogen-Bond Formation between the Wild Type and CYP2A13.8 in the All Trajectories for Last 10-ns Simulations

Acceptor	Donor	Wild type (%)	2A13.8 (%)
Asp158/Glu158 O	Gly162 N	94.42	91.44
Gly154 O	Asp158 N	81.46	84.28
Asp158 Oδ1	Arg161 Nη1	13.64	—
Asp158 Oδ1	Arg161 Nη2	13.64	—
Glu158 Oε2	Arg161 Nη1	—	63.48
Glu158 Oε1	Arg161 Nη1	—	52.48
Glu158 Oε2	Arg161 Nε	—	50.56
Glu158 Oε1	Arg161 Nε	—	46.06
Glu151 Oε1	Arg148 Nη2	59.92	0.400
Glu151 Oε2	Arg148 Nη2	53.90	1.480
Glu151 Oε1	Arg148 Nη1	40.64	0.300
Glu151 Oε2	Arg148 Nη1	33.48	0.320
Glu152 Oε1	Arg148 Nη2	0.000	67.42
Glu152 Oε1	Arg148 Nε	0.000	64.02
Glu146 Oε2	Arg143 Nη2	99.82	0.000
Glu146 Oε2	Arg143 Nη2	0.000	0.340
Glu146 Oε1	Arg446 Nη2	0.000	86.10
Glu146 Oε2	Arg446 Nη2	0.000	84.70
Glu146 Oε1	Arg446 Nη1	0.000	59.98
Glu146 Oε2	Arg446 Nη1	0.000	49.90
Leu444 O	Arg446 Nη1	83.26	4.200

Fig. 7. Structural Comparison between the Wild Type and CYP2A13.8

The (A) wild type and (B) CYP2A13.8 are shown in green and pink, respectively. Nitrogen, oxygen, and hydrogen are displayed in blue, red, and white, respectively, in the stick model. Iron is presented as an orange sphere using a model written as the van der Waals radius. The dotted lines indicate the ionic and hydrogen bonds.

Mutation in CYP2A13.9 Affects Interaction with Heme

To investigate the mechanism behind the reduction in enzymatic activity in CYP2A13.9, we scrutinized the structural changes in the calculated structures for the V323L variant. The substitution from valine (Val) to leucine (Leu) introduces a slight increase in bulkiness. The marginally larger bulkiness of Leu compared to Val induces a subtle distortion in helix J and diminished the rate of hydrogen-bond formation between main-chain atoms (Table 6). Although Tyr351 and Val355 are located in proximity to Val323/Leu323 in the 3D structure, no significant differences were observed in the orientation of their side chains. The slight increase in bulkiness caused the position of Helix K to shift, consequently affecting the position of Helix J' as well (Fig. 8A). As a result, the interactions between Arg143 and Glu344 were increased in CYP2A13.9 (Table 6). Moreover, the helix shifts induce a conformational change in the side chain of Phe343 located on Helix J', along with that of Met447 forming CH–π stacking with Phe343 (Figs. 8B, C). The shift of the Met447 side chain (the distance of Cε was 5.9 Å in the final structures) extended to the main chain, influencing the structure of the N-terminal side of Helix L. Specifically, the position and conformation of Ala445 on the same helix and Phe440 on the K”L loop undergo changes. In CYP2A13.9, Phe440 was closer to Helix B. As Ala445, Phe440, and the B'C loop, where structural changes were observed, are in contact with heme, it is suggested that the decreased activity in CYP2A13.9 is attributed to a weakening of the interaction with heme.

Table 6. Comparison of Hydrogen-Bond Formation between the Wild Type and CYP2A13.9 in the All Trajectories for the Last 10-ns Simulations

Acceptor	Donor	Wild type (%)	2A13.9 (%)
His320 O	Val323/Leu322 N	97.22	88.80
Val323/Leu323 O	Val327 N	96.44	79.74
Val323/Leu323 O	Lys326 N	57.42	46.00
Glu344 Oε1	Arg143 Nη2	0.000	79.54
Glu344 Oε2	Arg143 Nη2	0.000	74.76
Glu344 Oε2	Arg143 Nη1	0.000	73.18
Glu344 Oε1	Arg143 Nη1	0.000	70.72

Fig. 8. Structural Comparison between the Wild Type and CYP2A13.9

The (A) wild type and (B) CYP2A13.9 are presented in green and light blue, respectively. Nitrogen, oxygen, hydrogen, and sulfur are displayed in blue, red, white, and yellow, respectively, in the stick model. Iron is shown as an orange sphere using a model written as the van der Waals radius. The dotted lines indicate the ionic and hydrogen bonds. The yellow arrows show the shifts of Helices J’ and K.

CONCLUSION

The conformational changes in CYP2A13 variants were explored through MD simulations, considering structures that incorporate heme. The observed changes of the hydrogen-bond formation are summarized in Table 7. In most variants, there were indications of altered interactions with heme, suggesting that heme misalignment or destabilization of binding could account for the observed reduction in variant activity. In CYP2A13.5, the mutation affected multiple helices, including helix K, implying that diminished interactions with redox partners might underlie the decline in enzymatic activity. Comparatively, the effects of the mutations on 3D structures were deemed minor when contrasted with the polymorphism of CYP2A6.²²⁾ Several cases of undetectable enzymatic activity and reduced substrate affinity have been reported in CYP2A6 polymorphism.³¹⁾ However, in the examined variants of this study, enzymatic activity was detectable in all cases except for CYP2A13.4, and there were no significant differences in K_m values when nicotine was used as the substrate (Table 1). The calculated results, reflecting relatively subtle structural changes, were considered reasonable.

Table 7. Changes of the Hydrogen-Bond Formation in the Variants

Variant	Increased hydrogen bond	Decreased hydrogen bond
CYP2A13.4	Heme-Arg437	Heme-Arg101, Ala371-Arg101/Gln101, Ala117-Arg101/Gln101
CYP2A13.5	Leu310-Tyr453, Phe314-Tyr453	—
CYP2A13.6	Cys494-Lys337, Cys494-Arg461, Ther212-Phe480, Gln218-Phe480	Asp332-Arg494, Met457-Arg494, Gln458-Arg494
CYP2A13.8	Asp158/Glu158-Arg161, Glu152-Arg148, Glu146-Arg446	Glu151-Arg148, Leu444-Arg446, Glu146-Arg143
CYP2A13.9	Glu344-Arg143	—

The reported allele frequencies of CYP2A13.4 and CYP2A13.5 was lower than those of CYP2A13.6, CYP2A13.8, and CYP2A13.9.^15,16) In CYP2A13.4, the mutation affected directly the hydrogen-bond formation with heme. The RMSF changes were observed in CYP2A13.5 and CYP2A13.8. However, no variant investigated in this study indicated drastic structural changes. It would be interesting to find a relationship between the frequency of allele occurrence, structure, and activity in the future.

The MD simulation outcomes from this study will contribute to enhancing our ability to predict structural changes in CYP2A13 variants, offering valuable insights into the genotype-phenotype relationship. It is essential to note that this investigation focused solely on variants with a single mutation; hence, future studies involving variants with multiple mutations are warranted to provide a more comprehensive understanding. In addition, MD simulations for the complexes of the CYP2A13 variants and substrates will be valuable for investigating the effects on the substrate recognition in the variants. The development of biosensors that can rapidly measure drugs and endogenous substances in plasma using CYPs at low cost is expected,³²⁾ and the information obtained by the MD simulations about CYP structural changes in various variants will be useful for this purpose.

Acknowledgments

This work received support from Grants-in-Aid for Scientific Research (Grant Nos. 17K08257 and 21K15244) and the Japan Society for the Promotion of Science.

Conflict of Interest

The authors declare no conflict of interest.

REFERENCES

• 1) Su T, Bao Z, Zhang QY, Smith TJ, Hong JY, Ding X. Human cytochrome P450 CYP2A13: predominant expression in the respiratory tract and its high efficiency metabolic activation of a tobacco-specific carcinogen, 4-(methylnitrosamino)-1-(3-pyridyl)-1-butanone. Cancer Res., 60, 5074–5079 (2000).
• 2) von Weymarn LB, Murphy SE. CYP2A13-catalysed coumarin metabolism: comparison with CYP2A5 and CYP2A6. Xenobiotica, 33, 73–81 (2003).
• 3) Bao Z, He XY, Ding X, Prabhu S, Hong JY. Metabolism of nicotine and cotinine by human cytochrome P450 2A13. Drug Metab. Dispos., 33, 258–261 (2005).
• 4) Murphy SE, Raulinaitis V, Brown KM. Nicotine 5′-oxidation and methyl oxidation by P450 2A enzymes. Drug Metab. Dispos., 33, 1166–1173 (2005).
• 5) Guo Y, Zhu LR, Lu G, Wang H, Hong JY. Selective expression of CYP2A13 in human pancreatic α-islet cells. Drug Metab. Dispos., 40, 1878–1882 (2012).
• 6) Gandini S, Botteri E, Iodice S, Boniol M, Lowenfels AB, Maisonneuve P, Boyle P. Tobacco smoking and cancer: a meta-analysis. Int. J. Cancer, 122, 155–164 (2008).
• 7) Nakajima M, Itoh M, Sakai H, Fukami T, Katoh M, Yamazaki H, Kadlubar FF, Imaoka S, Funae Y, Yokoi T. CYP2A13 expressed in human bladder metabolically activates 4-aminobiphenyl. Int. J. Cancer, 119, 2520–2526 (2006).
• 8) Borlak J, Walles M, Levsen K, Thum T. Verapamil: metabolism in cultures of primary human coronary arterial endothelial cells. Drug Metab. Dispos., 31, 888–891 (2003).
• 9) Smith BD, Sanders JL, Porubsky PR, Lushington GH, Stout CD, Scott EE. Structure of the human lung cytochrome P450 2A13. J. Biol. Chem., 282, 17306–17313 (2007).
• 10) DeVore NM, Scott EE. Nicotine and 4-(methylnitrosamino)-1-(3-pyridyl)-1-butanone binding and access channel in human cytochrome P450 2A6 and 2A13 enzymes. J. Biol. Chem., 287, 26576–26585 (2012).
• 11) Gotoh O. Substrate recognition sites in cytochrome P450 family 2 (CYP2) proteins inferred from comparative analyses of amino acid and coding nucleotide sequences. J. Biol. Chem., 267, 83–90 (1992).
• 12) von Weymarn LB, Zhang QY, Ding X, Hollenberg PF. Effects of 8-methoxypsoralen on cytochrome P450 2A13. Carcinogenesis, 26, 621–629 (2005).
• 13) Vrzal R. Genetic and enzymatic characteristics of CYP2A13 in relation to lung damage. Int. J. Mol. Sci., 22, 12306 (2021).
• 14) Bridges A, Gruenke L, Chang YT, Vakser IA, Loew G, Waskell L. Identification of the binding site on cytochrome P450 2B4 for cytochrome b₅ and cytochrome P450 reductase. J. Biol. Chem., 273, 17036–17049 (1998).
• 15) Fujieda M, Yamazaki H, Kiyotani K, Muroi A, Kunitoh H, Dosaka-Akita H, Sawamura Y, Kamataki T. Eighteen novel polymorphisms of the CYP2A13 gene in Japanese. Drug Metab. Pharmacokinet., 18, 86–90 (2003).
• 16) Cauffiez C, Lo-Guidice JM, Quaranta S, Allorge D, Chevalier D, Cenée S, Hamdan R, Lhermitte M, Lafitte JJ, Libersa C, Colombel JF, Stücker I, Broly F. Genetic polymorphism of the human cytochrome CYP2A13 in a French population: implication in lung cancer susceptibility. Biochem. Biophys. Res. Commun., 317, 662–669 (2004).
• 17) Tamaki Y, Honda M, Muroi Y, Arai T, Sugimura H, Matsubara Y, Kanno S, Ishikawa M, Hirasawa N, Hiratsuka M. Novel single nucleotide polymorphism of the CYP2A13 gene in Japanese individuals. Drug Metab. Pharmacokinet., 26, 544–547 (2011).
• 18) Kumondai M, Hosono H, Maekawa M, Yamaguchi H, Mano N, Oda A, Hirasawa N, Hiratsuka M. Functional characterization of 9 CYP2A13 allelic variants by assessment of nicotine C-oxidation and coumarin 7-hydroxylation. Drug Metab. Pharmacokinet., 33, 82–89 (2018).
• 19) Kobayashi K, Takahashi O, Hiratsuka M, Yamaotsu N, Hirono S, Watanabe Y, Oda A. Evaluation of influence of single nucleotide polymorphisms in cytochrome P450 2B6 on substrate recognition using computational docking and molecular dynamics simulation. PLOS One, 9, e96789 (2014).
• 20) Fukuyoshi S, Kometani M, Watanabe Y, Hiratsuka M, Yamaotsu N, Hirono S, Manabe N, Takahashi O, Oda A. Molecular dynamics simulations to investigate the influences of amino acid mutations on protein three-dimensional structures of cytochrome P450 2D6.1, 2, 10, 14A, 51, and 62. PLOS One, 11, e0152946 (2016).
• 21) Watanabe Y, Fukuyoshi S, Hiratsuka M, Yamaotsu N, Hirono S, Takahashi O, Oda A. Prediction of three-dimensional structures and structural flexibilities of wild-type and mutant cytochrome P450 1A2 using molecular dynamics simulations. J. Mol. Graph. Model., 68, 48–56 (2016).
• 22) Kato K, Nakayoshi T, Nokura R, Hosono H, Hiratsuka M, Ishikawa Y, Kurimoto E, Oda A. Deciphering structural alterations associated with activity reductions of genetic polymorphisms in cytochrome P450 2A6 using molecular dynamics simulations. Int. J. Mol. Sci., 22, 10119 (2021).
• 23) Oda A, Kobayashi K, Takahashi O. Computational study of the three-dimensional structure of N-acetyltransferase 2-acetyl coenzyme a complex. Biol. Pharm. Bull., 33, 1639–1643 (2010).
• 24) Joung IS, Cheatham TE 3rd. III Determination of alkali and halide monovalent ion parameters for use in explicitly solvated biomolecular simulations. J. Phys. Chem. B, 112, 9020–9041 (2008).
• 25) Case DA, Ben-Shalom IY, Brozell SR, et al. AMBER 2018, University of California, San Francisco (2018).
• 26) Maier JA, Martinez C, Kasavajhala K, Wickstrom L, Hauser KE, Simmerling C. ff14SB: improving the accuracy of protein side chain and backbone parameters from ff99SB. J. Chem. Theory Comput., 11, 3696–3713 (2015).
• 27) Oda A, Yamaotsu N, Hirono S. New AMBER force field parameters of heme iron for cytochrome P450s determined by quantum chemical calculations of simplified models. J. Comput. Chem., 26, 818–826 (2005).
• 28) Giammona DA. Thesis, University of California Davis (1984).
• 29) Darden T, York D, Pedersen L. Particle mesh Ewald: an N-log(N) method for Ewald sums in large systems. J. Chem. Phys., 98, 10089–10092 (1993).
• 30) Ryckaert J-P, Ciccotti G, Berendsen HJC. Numerical integration of the cartesian equations of motion of a system with constraints: molecular dynamics of n-alkanes. J. Comput. Phys., 23, 327–341 (1977).
• 31) Hosono H, Kumondai M, Maekawa M, Yamaguchi H, Mano N, Oda A, Hirasawa N, Hiratsuka M. Functional characterization of 34 CYP2A6 allelic variants by assessment of nicotine C-oxidation and coumarin 7-hydroxylation activities. Drug Metab. Dispos., 45, 279–285 (2017).
• 32) Sakaki T. Practical application of cytochrome P450. Biol. Pharm. Bull., 35, 844–849 (2012).