Molecular Dynamics in the Study and Development of Molecu- larly Imprinted Materials – Status Quo, Quo Vadis?

The past two decades have witnessed the introduction of and then a steady increase in the use of computational techniques in the study and development of molecularly imprinted polymers (MIPs). Molecular dynamics (MD) based studies have had a significant role in this development as they can provide insights concerning the mechanisms governing the molecular level events underlying MIP synthesis and MIP-ligand interactions and can be used for the identification of preferred monomer compositions and for the prediction of MIP properties. We here review the role that MD has played in the development of molecular imprinting and examine the different types of MD strategies that have been used, including their advantages and challenges. Recent trends in the application of MD to the study of MIPs are presented, along with a perspective on the future importance of MD-based studies for the development of molecular imprinting science and technology.

While the molecular imprinting concept is in principle simple, Figure 1, the underlying molecular level events are many and complex [26]. Attempts to elucidate the molecular level mechanisms for MIP synthesis and MIP-ligand recognition, and to rationally design MIPs, have been significant for the development of the field [4,5]. These efforts have been facilitated by the introduction of a range of methods for the study and design of MIPs and for the development of novel polymer compositions and formats. Prior to the introduction of computational techniques into the field, researchers had explored the use of various theoretical treatments [27][28][29], combinatorial polymer synthesis approaches [30,31] and various spectroscopy-based studies [32][33][34][35]. It was, however, with the introduction of computational chemistry into the molecular imprinting field, particularly when used together with experimental studies, that significant interest in rational MIP design and understanding of molecular level details of MIP systems began to evolve, as can be followed through a series of reviews of the area over the period [36][37][38][39][40][41][42][43][44][45][46][47][48][49]. The increased use of computational chemistry for better understanding the various aspects of the molecular imprinting process and for the design of novel MIPs has been driven by the establishment of more powerful computers, better theoretical descriptors, and purpose-oriented software. Collectively, these developments have, and continue to, pave the way for more precise and efficient modeling of molecular level structure and interactions in molecular imprinting systems, as reflected in the growing numbers of publications deploying computational chemistry [48]. Three broad classes of computational methods have been applied to the study and development of aspects of the molecular imprinting process and for understanding specific MIP systems: quantum mechanicsbased studies, statistical studies (most notably multivariate analysis) and molecular mechanics-based molecular dynamics studies [48]. With respect to MIP design and a broader understanding of the molecular-level mechanisms underlying the various stages of the molecular imprinting process, it is molecular dynamics-based studies that can be argued to have had, to this point in time, the greatest impact on molecular imprinting, and is a trend that continues, Figure 2. This development has arisen due to the capability to model complete systems with multiple examples of components in experimentally relevant stoichiometries, together with a documented capacity for validation or comparison with experimental studies and has resulted in a variety of studies for prognostic or diagnostic purposes [48].

Figure 2.
Number of published papers where MD simulations have been applied to molecular imprinting systems in order to either select MIP components (fragment screening), perform comprehensive studies of pre-polymerization interactions (full system) or create polymer models (structure and function).
Here we present a review of the establishment, development, current status and projected future (status quo, quo vadis?) of MD-studies in the molecular imprinting field. This review covers, and illustrates with examples from seminal studies, from the initial simulated annealing-based studies pioneered by the Piletsky group [50] and the broader use of fragment-screening-like [51,52] studies for ranking candidate monomers through to the introduction of full system all-atom studies for diagnostic and prognostic purposes pioneered by the Mizaikoff and Nicholls groups [53,54], all with a primary focus in the interaction of a number of pre-polymerization mixture components. The subsequent application of MD-studies for understanding or predicting molecular recognition phenomena in MIP systems through, most notably, site modeling and theoretical docking of ligands with in silico MIP-recognition sites is presented, together with a projection of the importance of MD-based techniques for the future of the molecular imprinting field.

Fragment screening-like simulations of template-monomer interactions
Since the ability of a MIP to selectively recognize and bind an analyte is derived from the complexes formed between template and functional monomers in the pre-polymerization solution, researchers have simulated this process in attempts to select the most appropriate monomer for a given template. In 2001 Piletsky and co-workers described a method for screening the template ephedrine against a library of 20 functional monomers [50]. The method has since been successfully applied in a range of studies . In this protocol [94], an initial scoring function based on a Leapfrog-algorithm is followed by an MD-based annealing simulation, Figure 3. In the scoring step, interaction points on the template are identified based on electronic, steric, and lipophilic properties (colored crosses in Figure 3a) and the monomers in the library are then rotated around these points and interaction energies are calculated. Next, the template is placed in a virtual box saturated with the monomer(s) identified in the scoring step followed by an annealing simulation using MD. Finally, the system is energy minimized and the complexes formed are examined to determine a suitable monomer-template ratio, Figure 3b.  after an annealing simulation and energy minimization. Reproduced from [94].
The above type of MD-based studies of isolated MIP components or systems of limited size require relatively little computational cost and can provide valuable guidance to the composition of polymer systems. However, to obtain a more holistic view of the multitude of equilibria existing in pre-polymerization mixtures, and to reflect the polyclonality of non-covalent MIPs, simulations require the presence of multiple copies of template together with monomers (both functional and cross-linking) and solvent in relevant stoichiometries. These more comprehensive, all component MD simulations of MIP prepolymerization mixtures are highlighted in the following section.

Simulations of all-component systems
A combination of the inherent complexity of molecular imprinting systems and the limitations of computer power and software together long impeded the application of MD-based techniques to systems reflecting those used in the laboratory for polymer synthesis, i.e. with multiple copies of templates and corresponding stoichiometries for monomers, explicit solvent and initiator. The parallel developments in software and access to computing power made the first all-atom full system studies possible [41]. Today, the simulation of all-component systems with stoichiometries corresponding to experimental pre-polymerization mixtures is both possible, and quite readily accessible [159].
Essentially, MD-simulations of MIP pre-polymerization mixtures representative of synthetic protocols are undertaken after several equilibration steps, and the results of the production runs, effectively a temporal array of molecular trajectories, are subsequently analyzed, Figure 4a. Detailed information regarding the nature and extent of interactions between the molecular species can then be extracted and correlated with MIP recognition performance.
(a) (b) Figure 4. Examples of all-atom, all-component simulation systems: (a) The components of a phenylalanine anilide imprinted polymer are randomly placed in a simulation box (top). After several equilibration steps (bottom left) the system is subjected to MD simulations, whereafter the resulting molecular trajectories are analyzed using e.g. radial distribution functions, hydrogen bond-and grid density analyses (bottom right). Reproduced from [41] with permission from Springer; (b) Snapshot from simulations of miniemulsion polymerization systems. Adapted from [160].
Several studies with the aim to explain the molecular basis of recognition have demonstrated the significance of including all the polymer building blocks in the simulation by illuminating the nature of complex formation and underlying mechanisms [53,54, [160][161][162][163][164][165][166][167][168][169]. Some strategies involved the evaluation of the choice of suitable polymer components such as the type of functional monomer, crosslinker, template or porogen [170][171][172][173][174][175][176] and/or the investigation of the relative stoichiometries between the monomers [172,[175][176][177][178][179][180][181][182][183][184][185][186][187][188] suitable for use in imprinting protocols. A recent example by Rebelo et al. [183] reported the development of a MIP-based electrochemical sensor, comprised of a carbon paste electrode modified with MIP microparticles and multi-walled carbon nanotubes, for detection of the antibiotic furazolidone in environmental waters. The group employed MD simulation of six all-atom, all-component pre-polymerization mixtures to determine the most favorable template to monomer molar ratio to be used in the preparation of the sensor. Another comprehensive theoretical study was recently performed by Garcia and colleagues [184] using acrylate-based monomers in the design of nano-MIPs selective for L-fucose and D-mannose as prototype disease biomarkers. Here, MD simulations of fully atomistic pre-polymerization mixtures displaying various molar ratios between the monomers were performed to determine the most favorable composition to maximize recognition potential towards the templates In a series of reports, Golker et al. [175,176,[185][186][187] explored the relationships between MIP composition, final MIP morphology and recognition. By employing extensive MD simulations of pre-polymerization mixtures displaying a broad range of variations in composition and/or relative stoichiometry, correlations between the nature and extent of interactions observed in the pre-polymerization phase with results from physical characterization studies and equilibrium re-binding studies were described. Similar relationships between MIP composition, morphology and recognition were also investigated by Olsson et al. in a study of bisphenol A-imprinted nanoparticles synthesized by miniemulsion polymerization, Figure 4b [160].
It is commonly accepted that monomer-template interactions in the pre-polymerization phase, if preserved during polymerization, represent the basis for recognition in MIPs, and as evident from the published literature reviewed here, MD simulations representing fully atomistic all-component MIP pre-polymerization mixtures have paved the way for a deeper understanding of the molecular events underlaying recognition. Nevertheless, multifaceted MD-based simulation approaches relating to other aspects of the imprinting process have been reported and are described in the following section.

Simulations of MIP structure and function
Since force-field based methods are unable to account for electron transfer necessary for bond formation, to date, naturally, most efforts on utilizing MD simulations in MIP development have been focused on modeling solution aspects prior to polymer gelation. Yet, a number of studies have attempted to model the polymerization reaction [189][190][191][192][193][194]. In a paper by Cowen et al. [193], a novel algorithm capable of mimicking atomistic radical polymerization of all-component pre-polymerization systems was presented. The usefulness and accuracy of the algorithm was demonstrated by simulating an ephedrine imprinted nanoparticle, Figure. 5a, and comparing the simulation data in terms of relative affinity and selectivity, prior to and after polymerization, with theoretical and empirical results derived from earlier studies [50]. Schauperl and Lewis [192] presented a polymer model to investigate aspects of binding heterogeneity of xanthine imprinted polymers. A growing chain of functional-and crosslinking monomers was formed by successively adding the monomers to the model system. As the monomer chain was allowed to grow around one or more template molecules up to a pre-defined density, MD simulations were utilized to study template-polymer interaction.
(a) (b) Figure 5. (a) Ephedrine imprinted nanoparticle resulting from the polymerization algorithm by Cowen et al. [193]. Reproduced with permission from John Wiley and Sons; (b) Construction of 17β-estradiol imprinted tubes for virtual chromatography (top). Inside and lateral views of the nonimprinted (NVIP) and virtually imprinted (VIP) tubes (bottom). Adapted from [195] with permission from the authors.
Attempts of using the MD method to model binding [195][196][197][198][199][200][201][202][203][204][205][206][207] and/or structural aspects [208][209][210] of the resulting polymer have been reported. In a study by Mizaikoff and co-workers [195], hollow tubes were created by sequential stacking of rings of functional monomers, Figure 5b. The inside void was then filled with solvent and template and the whole assembly subjected to MD simulations, allowing the monomers to move along the Z-axis. The resulting imprinted tubes were used for "virtual chromatography", again using MD simulations, by pushing solvent and template or analogues through the tubes.
With the aim of identifying factors underlying MIP affinity and selectivity, Sobiech and coworkers [198][199][200]203] employed a MD-based strategy to mimic binding site models. The method involved MD simulation of pre-arranged complexes between template, functional monomer and crosslinker whereafter the template was removed, and the polymerizing step mimicked by creation of new bonds between the vinyl groups of monomers and crosslinkers. The resulting cavity was then MD-simulated together with the analyte to model binding. A similar approach was used by Gajda et al. [201] to simulate a MIP binding site for aripiprazole, by Janczura et al. [204] to probe the adsorption behavior and selectivity of a MIP model for determination of nitroxidative stress products in human urine samples and by Zeng et al. [207] to simulate the adsorption process of a pHresponsive switchable MIP cavity with affinity for tylosin, a diol-containing macrolide antibiotic. Mazouz and co-workers [202] utilized MD simulations to develop a binding model of a polypyrrol-based MIP sensor for detection of prostate specific antigen. In the study, monomer-template interactions in terms of adsorption energies were characterized and related to the ability of the polymer to generate imprinted binding sites for the protein. A coarse-grained simulation approach of protein imprinted hydrogels was presented by Zadok and Srebnik [210] where all-component MD simulations of pre-polymerization mixtures in combination with Monte Carlo methods were used to describe the effect of composition on template complexation and material characteristics Even though the calculations describing the intermolecular forces are based on classical Newtonian mechanics, hence providing simplified approximations of the molecular level events, the diversity of studies exploring the MD method to model the various stages of the imprinting process reflect the capability of the method to provide a holistic picture that can aid in the design of complex molecularly imprinted systems.

Conclusions
MD-based studies are having an increasing impact on the development of molecular imprinting science and technology. Since the introduction of MD-based fragment-screening-like studies of template-monomer interactions for monomer selection and subsequent full system simulations, MD-based studies, after successful benchmarking with experimental techniques, have found legitimacy as prognostic and diagnostic tools for the assessment of MIP systems. More recent years have even seen the application of MD for the prediction of MIP recognition and morphological properties and even MIP function. The ongoing improvements in computer power, software capabilities and understanding of the molecular imprinting process should reinforce the importance of in silico studies of MIP-systems, and in particular MD-based studies, over the coming years.