Computational-Assisted Development of Molecularly Imprinted Polymers for Synthetic Cannabinoid Recognition

TL;DR

This paper uses computational methods to design polymers that can selectively detect synthetic cannabinoids, which are harmful drugs spreading globally.

Contribution

The novel use of computational modeling to optimize molecularly imprinted polymers for synthetic cannabinoid recognition.

Findings

TFAA and BA were found to form the most stable complexes with synthetic cannabinoids due to acidity and aromatic interactions.

Six solvents were evaluated for solvation energy, identifying suitable candidates for polymerization.

Computational predictions guide efficient experimental validation of MIPs for SC extraction.

Abstract

Synthetic cannabinoids (SCs), a prominent class of new
psychoactive
substances, pose growing challenges to public health due to their
severe toxic effects and widespread global presence. In this study,
we employed computational methods to develop molecularly imprinted
polymers (MIPs) for the selective recognition of seven SCs, chosen
based on seizure reports from the Narcotics Examination Unit of the
Scientific Police of the State of São Paulo. Density functional
theory and extended tight binding for geometry, frequency, and noncovalent
model 2 (GFN2-xTB) calculations were used to optimize the molecular
geometries and predict ideal monomer–solvent combinations for
MIP synthesis. We assessed six solventsacetone, acetonitrile,
dichloromethane, chloroform, diethyl ether, and dimethyl sulfoxidebased
on their solvation energy, identifying suitable candidates for the
polymerization step. Hydrogen bonding interaction sites were mapped,
guiding the selection of functional monomers such as acrylic acid
(AA), 4-vinylbenzoic acid (BA), 2-(trifluoromethyl)­acrylic acid (TFAA),
and methacrylic acid. Our findings suggest that TFAA and BA offer
the most stable complexation with SCs, influenced by their acidity
and aromatic interactions. These computational predictions pave the
way for resource-efficient experimental validation and enhance the
development of MIPs as tools for the extraction of SCs in complex
matrices, contributing to efforts to combat the global SC epidemic.

Figures

General structural framework of the seven SCs studied in this work, with their particular substituents and functional groups.

Scaffold of SCs with four molecular regions highlighting the tail (in blue), core (in red), linker (in orange), and linked group (in yellow).

Total solvation energy for the SCs in the solvents acetone (Ace), acetonitrile (ACN), chloroform (Chl), dichloromethane (DCM), diethyl ether (EtOEt), and dimethyl sulfoxide (DMSO), calculated with the SMD model.

ESP for molecule 6, where negative density regions are represented by red color, positive by blue color, and intermediate by green color.

Pattern for HB formation between a carboxylic acid and an MDMB group (a) and ADB group (b).

FMs applied in this work, with their respective names and abbreviations.

Complexation energy (E complexation) for SCs and FMs.

Representatives of each interaction group with their respective strongest interaction energy for a given FM.

Complexation energy for SCs in protics solvents, water (H2O), methanol (MeOH), ethanol (EtOH), and isopropyl alcohol (i-PrOH).

Full text

Introduction

1

Since the early 2000s, synthetic cannabinoid (SC) receptor agonists, or more simply SCs, have been present in the global scene of new psychoactive substances, posing an increasingly significant challenge for security and public health agencies worldwide.? SCs constitute the largest class of new psychoactive substances, with 353 substances already notified to the United Nations Office of Drugs and Crime (UNODC) by 84 countries through the Early Warning Advisory system.

Mirroring the effects of THC (i.e., (−)trans-Δ^9^-tetrahydrocannabinol) on the endocannabinoid system’s CB1 and CB2 receptors, many SCs available on the illicit market today appeared from research aimed at the potential use of THC analogs or substances capable of activating cannabinoid receptors for therapeutic use, such as pain reduction, antiepileptic effect, anti-inflammatory, antiemetic, and antineoplastic activity, without causing dependence. ?,? However, many adverse effects are associated with SCs, such as hallucinations, anxiety, suicidal thoughts, rhabdomyolysis, cardiovascular effects, as well as severe intoxications leading to cerebral and cardiac ischemia, heart attack, and even death.?

In Brazil, in 2023, the number of intoxication cases involving mixtures of more than one SC increased significantly. According to data from the Municipality of São Paulo, in the first report compiled by the Surveillance Coordination of the Municipal Health Secretariat, there were 493 notifications of suspected cases of SC intoxication in 2023.? In the neighboring state of Rio de Janeiro, SCs have also been identified within prisons, showing similar patterns of substances found in São Paulo, with a greater emphasis on MDMB-4en-Pinaca, ADB-Butinaca, and 5F-MDMB-Binaca.?

Given the risk to public health and the imminent global epidemic related to SCs, it becomes imperative to develop methods that assist in the detection of these compounds. One strategy described in the literature is the use of molecularly imprinted polymers (MIPs), which have the ability to selectively extract a specific molecule from a complex matrix.? This is a significant advantage, as these SCs are often found in matrices such as cellulose or vegetables, mixed with other controlled or uncontrolled substances.?

MIPs are based on the formation of specific cavities on the surface of a polymer, which mimic the shape of a template molecule.? The application process begins with the interaction of monomers with the target molecule, followed by polymerization of the monomers. Finally, the target molecules are removed, leaving only the desired cavity. In developing an MIP system for a specific molecule, it is essential to determine the optimal combination of monomers and solvents for each stage of the production process. ?,?

For example, different strategies are being developed to propose MIPs for recognizing certain natural and SCs. These strategies are based on the ability of MIPs to extract specific molecules, followed by trace level detection analysis. In the literature, the most commonly used technique for this detection is liquid chromatography tandem mass spectrometry (LC–MS/MS). ?−? ? ? However, there are also studies that employ Raman spectroscopy in solution? and voltammetry,? demonstrating the versatility of MIP applications for cannabinoid preconcentration. In addition, Yang et al., in 2023, conducted a study on cannabidiol extraction using MIP-coated magnetic nanoparticles, aiming to improve extraction yield to meet the market demand for therapeutic cannabidiol.?

The experimental proposition of MIPs can be resource-intensive, involving significant amounts of reagents and solvents. An alternative to mitigate these costs is the application of computational methods to predict a better combination of monomers and solvents for MIP production. ?,? Different computational strategies can be employed to evaluate the steps to produce efficient MIPs, as well described by Mizaikoff et al.? In general, when the template is a macromolecule, as proteins, and to study the polymerization, which involves a large number of atoms, molecular dynamics simulations are the most appropriate strategy.? On the other hand, quantum mechanics-based methods are computationally expensive but accurate in calculating the interactions among the molecules and the possibility of (de)­protonation, for example. Density functional theory (DFT) is a quantum mechanics method with a good cost-accuracy ratio. So, for studying systems with a large number of atoms, DFT methods became very useful, particularly during the prepolymerization steps, since the interactions between the template and the functional monomers (FMs) are key to this stage.?

In a previous study,? our group employed a computational approach for the development of MIPs targeting Δ^9^-THC and THC–COOH. To broaden the selection of FMs, the complexation energy between these target molecules and seven different FMs was estimated by using DFT calculations. Furthermore, the effect of four different solvents was incorporated through an implicit solvent model.? This approach is adopted to identify the most suitable FM for a given target molecule and has been successfully applied to other compounds, such as thiamethoxam,? atenolol,? and carvedilol.?

Computational methods play a crucial role in these studies by enabling the evaluation of various monomers and solvents and optimizing conditions for subsequent experimental assays. However, these often focus on a single target molecule rather than addressing a broader family of structurally related compounds, such as natural cannabinoids, cathinones, or SCs, where slight structural modifications can lead to significant differences in properties. Furthermore, the implicit solvent model has limitations, especially when used for protic solvents, as it may not accurately capture key solvent–molecule interactions.

So, while MIPs hold immense promise as a cost-effective recognition method, their experimental development often demands significant resources due to the myriad variables at play. Our computational approach contributes to modifying this landscape by assessing polymer composition, solvent selection, and monomers pK a. The last is an often overlooked but critical factor in such studies, since pK a evaluation may ensure that the hydrogen bonds formed in the complexes would be consistent with experimental data. We went beyond the one-molecule-at-a-time approach prevalent in this field,? by applying a combination of the xTB and DFT methods, which allowed for studying systems with a higher number of atoms at a lower computational cost than methods that rely solely on DFT. This approach enabled the addition of underexplored aspects, such as the combination of 7 SCs (templates) with 4 FMs and 4 explicit solvents, resulting in 56 complexes that were studied. In this way, we systematically investigated the optimal conditions to produce MIPs for seven different recently seized SCs reported by the Narcotics Examination Unit of the Scientific Police of the State of São Paulo (NEE-SP).

Methodology

2

The selection of the seven SCs studied in this work was carried out in partnership with the Narcotics Examination Unit of the Scientific Police of the State of São Paulo (NEE-SP). The NEE-SP conducted a reference study in 2021 on SCs that were submitted to the institution in cellulose matrices,? and from this report, the three most seized SCs were selected, 5F-MDMB-Pica (1), 4F-MDMB-Butinaca (2), and MDMB-4en-Pinaca (3, Figure). The remaining four SCs were selected based on the total number of seizures in the period from July 2022 to December 2023, considering all matrices: MDMB-Butinaca (4), ADB-Butinaca (5), ADB-4en-Pinaca (6), and ADB-Fubinaca (7, Figure).

General structural framework of the seven SCs studied in this work, with their particular substituents and functional groups.

Geometry optimization calculations were performed using GFN2-xTB (version 6.5.0)? at the most stringent convergence level, without applying solvent models. Only geometries without imaginary frequencies were considered, indicating minimum energy structures. Total energy calculations were carried out using DFT with the ωB97XD/6–31+G­(d,p) level of theory, ?,? as implemented in the Gaussian 09 package.?

The study of the interactions is a crucial aspect that requires a thorough evaluation. Initially, SCs are extracted from the matrices where they were commercialized, a process that typically involves the use of aprotic solvents. So, in this case, a continuous solvation model can be employed; here, the solvent effects were incorporated through the solvent model density (SMD) within the DFT framework.? The next step consists of evaluating the more adequate monomer to be used. Atomic charges and the electrostatic surface potential (ESP) were computed for the SCs using the ChelpG scheme,? allowing us to determine the regions for specific interactions with four different monomers. Molecular surfaces and intermolecular interactions were visualized with Binana 2.2 and Jmol 2.1.? The last step consists of washing the MIPs to remove the template (SC). In this case, strong hydrogen bonding may be established among the washing solvent and SCs, so four protic solvents were evaluated for each template. As hydrogen bonding must be evaluated, explicit solvent molecules were included at each proton donor or acceptor group in the calculations.

For the theoretical determination of the pK a value, all geometry optimizations were performed at the theoretical level ωB97XD/aug-cc-TZVP, using the implicit solvation model SMD, along with the addition of two explicit water molecules. ?,? For more details on this specific determination, see Supporting Information.

Results and Discussion

3

Selection of Synthetic Cannabinoids

3.1

The SCs studied in this work (Figure) have a standard scaffold, which can be divided into four molecular regions: tail, core, linker, and linked group (Figure).? For each of these regions, there is a wide variety of modifications; thus, the total number of synthetically possible cannabinoids is significantly large. However, there is a predominance, not only in Brazil but also worldwide, of a small recurrent number of derivatives. In other words, the SCs with the highest number of seizures also appear among the most prevalent SCs worldwide.?

Scaffold of SCs with four molecular regions highlighting the tail (in blue), core (in red), linker (in orange), and linked group (in yellow).

Therefore, the various substitutions recurrently found in SCs and that were studied in this work (Figure) within the context of MIP production are, specifically, tail substitutions (R_1_, including butyl, pent-4-enyl, 4-fluorobenzyl, 4-fluorobutanyl, and 5-fluoropentanyl groups) and two core structures, indole (R_2_ = CH) and indazole (R_2_ = N). The linker utilized in our study was carboxamide (i.e., CONH_2_), and the groups attached to this linker were tert-leucinamide and its methyl ester.

Solvent Selection for the Extraction of SC

from the Matrices

3.2

Strategically, to promote the association between a given SC and the monomer of an MIP through hydrogen bonding interactions, the polymerization process should occur in an organic, aprotic, and polar solvents. This allows for the solvation of the reagents of interest without competition with hydrogen bonding interactions that could arise from the solvent. Six solvents were selected for this purpose: acetone (Ace), acetonitrile (ACN), chloroform (Chl), dichloromethane (DCM), diethyl ether (EtOEt), and dimethyl sulfoxide (DMSO). To obtain the total solvation energy of the molecule (E sol(SMD)), the difference between the total energies of the molecule in vacuum (E SC) and in the presence of solvent (E SMD) was determined, following the SMD formalism, according to eq.

$$[eqn]$$

In general, from our data, Ace, ACN, and DCM emerge as the best candidates for the polymerization stage as well as for the extraction of the SC from the seizure matrices, given the higher interaction (i.e., more negative E sol(SMD)) with the SCs. Thus, we can conclude that solvents with moderate dipole moment values are preferable. Both the solvents that are the most and least polar, respectively, DMSO and EtOEt, did not present satisfactory solvation properties in the context of the analysis above, considering their lower interaction (i.e., less negative E sol(SMD)) values (Figureall values are also presented in Table S3).

Total solvation energy for the SCs in the solvents acetone (Ace), acetonitrile (ACN), chloroform (Chl), dichloromethane (DCM), diethyl ether (EtOEt), and dimethyl sulfoxide (DMSO), calculated with the SMD model.

Selection of Reaction Sites

3.3

The cannabinoid molecules were optimized using the semiempirical GFN2-xTB method, followed by determination of total energy in vacuum, with the theoretical level ωB97XD/6–31+G­(d,p). ChelpG maps (within a range of −0.1 to 0.1 eÅ^–3^) were plotted to better identify the reaction sites (Figure). In the Supporting Information, we provide all the representations of the electrostatic potential surfaces presented with the 2D molecular structures (Table S4).

ESP for molecule 6, where negative density regions are represented by red color, positive by blue color, and intermediate by green color.

Cannabinoids 1, 2, and 7 have three regions with a negative charge density, while the other compounds have only two regions with a negative charge density. After the regions of intramolecular interaction were disclosed, three patterns of HB formation were identified. The first pattern occurs only in SCs with a fluorine atom at the end of the tail group. The second pattern occurs in SCs of the MDMB class, which have an amide group and an ester group, with an arrangement that allows interaction with a HB donor group and another HB acceptor group in one face of the molecule and interaction with a HB donor group in the other face of the molecule, as shown in Figurea; such MDMB SCs are molecules 1, 2, 3, and 4. The third pattern occurs in the ADB class, which is composed of two amide groups; thus, on both faces of the molecule, there will be the possibility of interaction with an HB donor group and another HB acceptor group on each side, as shown in Figureb. All the aforementioned interaction patterns are suitable for carboxylic acid monomers, as this functional group can participate in hydrogen bonding both as a donor and an acceptor.

Pattern for HB formation between a carboxylic acid and an MDMB group (a) and ADB group (b).

Selection of the Functional Monomer

3.4

With the identification of potential intramolecular interaction points of HB nature and the selection of the carboxylic acid functional group as the most suitable within this context, among the selected FMs, there are acrylic acid (AA), methacrylic acid (MA), 2-(trifluoromethyl)­acrylic acid (TFAA), and 4-vinylbenzoic acid (BA), Figure. After the positioning of the monomers, considering the interactions and their quantities relative to each SC as described above, geometry optimization was performed by using the GFN2-xTB method. Upon obtaining the geometry at a minimum in the potential energy surface, the total energy was calculated using ωB97XD/6–31+G­(d,p).

FMs applied in this work, with their respective names and abbreviations.

The total complexation energy (E complexation) was determined using the total energies of the generated complex (E SC–FM), the isolated SC in vacuum (E SC), the isolated FM in vacuum (E FM), and the number of monomers present in each complex (n), according to eq.

$$[eqn]$$

Among the evaluated FMs, TFAA and BA exhibited the strongest interaction, as indicated by the E complexation values, meaning they were the FMs that most stabilize the SCs, as shown in Figure (and Table S5). This behavior parallels the trend observed in the pK a values of these acids, which are organized in increasing order, as follows: TFAA, AA, BA, and MA (see Table S3 for theoretical and experimental pK a values). ?,?,?,? The acidity of the monomers decreases in this order; thus, the ability to donate a hydrogen-bonding to a Lewis base (i.e., the SC framework) can be seen as increasing; indeed, TFAA presents the highest complexation energy (Figure) and lowest pK a value among the evaluated FMs.

Complexation energy (E complexation) for SCs and FMs.

Despite BA having a pK a value higher than that of TFAA, BA possesses an additional π-stacking interaction. In most complexes formed with BA, its aromatic ring is positioned parallel to the aromatic ring of the SC’s core, especially in the case of 7, which features an aromatic ring in its apolar extension, favoring this additional interaction.

Observing the range of complexation energies for each SC with a given FM (Figure), we can divide the SC into four groups. SCs with a fluorine atom in the tail moiety and MDMB as the linked group: 1 and 2; SCs formed with a tail moiety consisting only of an aliphatic chain and MDMB as the linked group: 3 and 4; SCs without the fluorine atom in the tail moiety and ADB as the linked group: 5 and 6; and a last group consisting only of the SC 7, with an aromatic ring and a fluorine atom in the tail moiety, and the presence of the ADB group as the linked group.

In Figure, all interactions between a representative of the four groups mentioned above and the FM resulting in the strongest interaction complexation energy are observed, and the other SCs and their interactions are presented in Table S6. For 2 and 6 (Figurea,c), four possible hydrogen bonds (black arrows, Figure) are formed between the ADB and the carboxylic acid groups, in addition to the hydrophobic interaction (gray circles, Figure) between the tail moiety and the aromatic ring present in the BA monomer. Justifying this set as being the strongest interaction, 2 was the only SC that exhibited a stronger interaction with FM TFAA (Figurea). Despite TFAA performing a hydrophobic interaction of lesser intensity compared to BA, it has a greater hydrogen bonding capacity than BA, and thus, the existence of a HB with the fluorine atom present in the tail moiety provides a higher binding energy with TFAA of −44.1 kcal/mol, as shown in Figure.

Representatives of each interaction group with their respective strongest interaction energy for a given FM.

Among all the calculated combinations, SC 7 and FM BA provided the strongest interaction of −46.42 kcal/mol. In this case, there is also a π–π stacking interaction between the aromatic ring of BA and the aromatic ring present in the tail moiety (Figured). This interaction occurs when there is an arrangement between the aromatic rings in such a way that they are face-to-face. This new interaction provides great stability to this set. Thus, in a mixture where all SCs are combined, the polymerization process with BA would be preferred for 7, given its strongest interaction.

Solvent for Washing the MIP

3.5

After polymerization, it is necessary to remove the SCs from the surface of the obtained polymer. Another instance where this removal is required is during the application of the MIP. After extracting the intact SC from the seized matrix, it is necessary to remove the SC for analysis, followed by cleaning the MIP for proper reuse.?

The most efficient cleaning method involves using a solvent that competes with the MIP for hydrogen bonding interactions with the SC.? Thus, a polar protic solvent is desirable. In this regard, the following solvents were selected: water (H_2_O), methanol (MeOH), ethanol (EtOH), and isopropyl alcohol (i-PrOH).

In this case, the use of a continuous solvation model, as applied before for the aprotic solvents, would not be suitable for determining the involved complexation energy because this type of model cannot include hydrogen bonding in its formalism.? Therefore, to construct a more representative system, it is necessary to form a complex between the SC and the solvent, analogous to the treatment applied to the systems of SCs and monomers. Such a complex is shown in Table S7.

Thus, it was observed that, in general, water did not yield satisfactory results (Figure and Table S8): its highly polar nature is incompatible with the structure of the SCs, which, despite being capable of hydrogen bonding, also possess nonpolar regions. On the other hand, isopropanol showed excellent results for all SCs, positioning itself as the most suitable solvent or very close to that position.

Complexation energy for SCs in protics solvents, water (H2O), methanol (MeOH), ethanol (EtOH), and isopropyl alcohol (i-PrOH).

Conclusion

4

This study provides a comprehensive theoretical framework for designing and optimizing MIPs tailored specifically for new SCs. The analysis identifies the most effective solvents and FMs based on their ability to form stable complexes with SCs, with special emphasis on hydrogen bonding interactions. Geometry optimization calculations were performed using GFN2-xTB, followed by total energy calculations at the ωB97XD/6–31+G­(d,p) level of theory, ensuring accurate structural characterization. Acetone, acetonitrile, and dichloromethane emerged as the most suitable solvents for the extraction and polymerization processes, as evaluated using the SMD solvation model within the DFT framework. At the same time, isopropyl alcohol was identified as the optimal solvent for MIP cleaning, a conclusion supported by explicit solvation calculations, which accounted for strong hydrogen bonding interactions. Among the evaluated FMs, TFAA and BA offered the best stabilization of SCs, with BA showing additional π-stacking interactions that further enhanced its efficacy. This work not only highlights the critical factors influencing MIP performance, such as the structural characteristics of SCs and their interaction patterns with monomers, but also provides practical guidance for the preparation of MIPs, integrating a stepwise computational approach to solvent and monomer selection.

Supplementary Material