Interactions of Nucleotide Bases with Decorated Si Surfaces from

ARTICLE pubs.acs.org/JPCC

Interactions of Nucleotide Bases with Decorated Si Surfaces from Molecular Dynamics Simulations Susanna Monti,*,† Giacomo Prampolini,‡ and Vincenzo Barone‡ † ‡

Istituto di Chimica dei Composti Organo Metallici (ICCOM-CNR), Area della Ricerca, via G. Moruzzi 1, I-56124 Pisa, Italy Scuola Normale Superiore, Piazza dei Cavalieri 7, I-56126 Pisa, Italy

bS Supporting Information ABSTRACT: All-atom molecular dynamics simulations and potential of mean force calculations have been carried out to define a complete picture of the adsorption properties of the four nucleotide bases 9-methyladenine, 9-methylguanine, 1-methylthymine, and 1-methylcytosine on a Si(111) surface functionalized with alkyl-amine molecules in aqueous solution. A detailed knowledge of the interactions between the free or bound nucleotide groups and the decorated surfaces is important to design efficient and selective DNA-based sensors, which could be used as powerful biotechnological detection systems. Indeed, the conformation of a DNA strand tethered to an electrode (probe) depends on the strength of its interactions with the layer, which influences the probe capability to capture the target oligonucleotides. However, the arrangement of each single DNA base in close proximity to a substrate and how it may be controlled to improve sensor performance have not received much attention. The results of this investigation showed that all of the bases had favorable adsorption properties on the amine layer, and these interactions could be sometimes mediated by interposing water molecules. The plane of the nucleotide rings was inclined and oriented so as to maximize the number of hydrogen bonds with the surface and all the molecules adsorbed with a similar strength, but a slight preference for pyrimidine bases was observed.

’ INTRODUCTION Biosensing devices such as DNA-based sensors,1
9 which rely on the ability of single-stranded DNA probes, bound to an appropriately functionalized substrate, to recognize and capture the complementary DNA targets which move in solution, have become increasingly important in designing new biomedical detection systems10
13 and finely tuned therapies.14
16 These kinds of biosensors are essentially composed of a recognition and a transduction unit. The recognition unit usually consists of a solid-state inorganic support (i.e., gold, oxidized silicon, crystalline silicon, etc.) which is coated with a stable organic film (self-assembled monolayer, SAM) in order to tailor in a convenient way its properties and thus the characteristics of the interface. The organic monolayer can be further modified for DNA grafting.17
19 The physical and chemical properties of the SAM covered surface are closely related to the structure of the SAM and to its specific attributes. In fact, thermal and chemical stability as well as surface polarity and hydrophobicity can have a strong impact upon nonspecific binding properties and could prevent DNA hybridization. Thus, an accurate investigation of the configuration and dynamics of the adsorbed chains might be important to understand the mechanisms which take place at the deposition interfaces and would help develop robust strategies for creating efficient and selective DNA-based devices. Notwithstanding monolayers have been characterized from an experimental point of view using different techniques including, for r 2011 American Chemical Society

example, Fourier transform infrared (FT-IR) spectroscopy,20 scanning tunneling microscopy (STM),21,22 X-ray reflectivity, and atomic force microscopy (AFM);23
25 only a small number of computational studies describing these systems at the atomic level have been reported to date.22,26
30 Even fewer are the investigations which have used a classical all-atom molecular dynamics approach to estimate the free energy change upon adsorption of DNA strands at the aqueous SAM interface.31
34 To the best of our knowledge, none of the cited works has examined in detail the adsorption and binding strength of single nucleotide bases on a decorated interface. Previous studies35
38 showed that both single and double stranded DNA structures tethered to a functionalized substrate could adopt an upright orientation relative to the surface and swing randomly in the bulk region of the solvent or bend, twist, and tilt toward the interface, forming relatively stable hydrogen bonding interactions with the head groups of the SAM. As a matter of fact, the permanence of DNA bases near the substrate and the frequency of their contacts with the SAM are closely related to the DNA base sequence and the affinity of each single nucleotide base for the surface head groups. Indeed, once the bases reach the Received: February 2, 2011 Revised: April 4, 2011 Published: April 19, 2011 9146

dx.doi.org/10.1021/jp201103u | J. Phys. Chem. C 2011, 115, 9146–9156

The Journal of Physical Chemistry C

ARTICLE

surface, they can spend some time there resting or moving (sliding) on top of the layer. In this work, the potential of mean force (PMF) molecular dynamics (MD) simulations are used to determine the free energy change upon adsorption of the four different nucleotide bases that form DNA strands, namely, adenine, thymine, guanine, and cytosine. As done by other authors (see, for example, refs 39
41), N-methylated nucleic acid bases 9-methyladenine (ADE, A), 1-methylthymine (THY, T), 9-methylguanine (GUA, G), and 1-methylcytosine (CYT, C) are considered in this study. A procedure similar to the one adopted in a previous investigation, concerning the adsorption free energy of aminoacid analogues at the aqueous titania interface,42 has been followed in order to obtain the binding trend of each nucleotide base as a function of the distance from the amine functionalized surface. The investigation of the individual binding properties of DNA bases will complement the studies already performed36
38 and will be able to explain particular binding behaviors observed experimentally.

’ MODELS AND METHODS Starting Models. The substrate model consisted of a H-passivated Si(111) surface functionalized with n-propylamine chains. Force field (FF) parameters for the functionalized silicon layer were previously30 derived from quantum mechanical calculations by means of the Joyce software.43 A H-passivated Si(111) unit cell (3.84  3.84 Å2) was extended along both the x and y directions to generate the final simulation substrate (≈ 19.8  22.8 Å2); a two-layer system made of 72 silicon atoms each was generated, and some of the terminal hydrogens of the Si(111) slab were replaced, randomly, with 18 propylamine chains. This substitution, which corresponds to a 50% coverage, turned out to be one of the most favorable arrangements of the chains in agreement with theoretical and experimental data, measured on similar systems26,27 The starting slab geometry was taken from a previous work30 and 720 TIP3P44 water molecules were added on top of the substrate (z dimension, ≈ 60 Å). A full atomistic model was employed except for hydrogen atoms bonded to the bottom silicon layer, which were grouped in a united atom representation. The model was relaxed from the starting conformation through a series of energy minimizations and constant volume molecular dynamics runs (at T = 300 K). The temperature of the system was maintained to 300 K using a Nose-Hoover thermostat.45,46 The silicon atoms of the bottom layer were frozen at their initial positions. Long range electrostatic interactions were treated with the particle mesh Ewald (PME) method47,48 using a spacing of 1.2 Å and a 4th order spline interpolation, while the short-range interactions were truncated at Rc = 9.5 Å. The equations of motion were integrated using a time step of 1 fs and the configurations were saved every 250 steps (0.25 ps). All simulations were performed with the GROMACS 4.0.7 package.49 Nucleic acid bases are attached to the DNA backbone through an “N-glycoside” bond between the sugar C10 and N9 of the purine (or N1 of the pyrimidine). We chose to replace the backbone portion with a methyl group (N-methylated nucleic acid bases) in order to prevent hydrogen bond formation of this region of the rings with the amine headgroups. The four N-methylated nucleic acid bases (Figure 1) parameters were

Figure 1. Stick representation of the four nucleotide bases used in this investigation.

extracted from previous simulations.36 Each nucleotide base was inserted in the pre-equilibrated periodic box, described above, in a location far from the layer and all the water molecules falling within 1.5 Å of the nucleotide base (five) were removed . Equilibration was repeated fixing the position of the nucleotide base. With this arrangement, a series of MD simulations were performed (see next subsection). Potential of Mean Force Calculations. The adsorption of each DNA base on the functionalized silicon surface was investigated by means of free energy profiles obtained through the potential of mean constraint force approach.42,50
53 The free energy difference ΔG(z) = G(z)
G(z0) between the nucleotide base at height z and in the bulk (z0) was calculated by integrating the component of the force exerted on the molecule in the direction of the constrained z, averaged over the constant ensemble (ÆFz(z)æ).42,53 dΔGðzÞ ¼
ÆFz ðzÞæ dz

ð1Þ

The position-dependent PMF could be reconstructed through a series of dynamic runs where the center of mass of the nucleotide fragment was positioned at fixed heights on the frozen silicon layer. Indeed, the reaction coordinate was the vertical separation between the bottom layer of the silicon slab, which was frozen during all the simulations, and the center of mass the nucleotide base. Despite the force fluctuated heavily, long simulations provided accurate averages.54 During each MD run, the molecule could diffuse freely in a plane parallel to the layer (xy plane), rotate in all directions around its center of 9147

dx.doi.org/10.1021/jp201103u |J. Phys. Chem. C 2011, 115, 9146–9156

The Journal of Physical Chemistry C

Figure 2. PMF free energy curve of ADE adsorption at the aqueous propylamine interface. The reaction coordinate is the distance between the ADE center of mass and the plane of the frozen silicon atoms. The adsorption geometry corresponding to the most favorable arrangement is displayed below. Silicon atoms are represented by light-yellow spheres. Possible hydrogen bonds are evidenced by dashed black lines.

mass, and explore different regions of the potential energy surface. About 20 runs per adsorbate were performed. Each run was 3.5 ns long, with the first 0.5 ns used for equilibration (this corresponded to a total simulation time of about 70 ns per adsorbate). The accuracy was estimated from separate averaging over blocks of 100 ps.42,55,53 After having identified probable binding arrangements of all the bases in close proximity to the amine layer, these configurations were used as starting structures for 5.5 ns unrestrained MD simulations. A new series of 10 ns MD runs was then carried out starting from the final geometries of the systems extracted from the unrestrained trajectories in order to further sample the configurational space. The trajectories explored the region close to the surface without any imposed restriction on the molecular coordinates, and the nucleotide base could also migrate far from the layer. These unrestrained dynamics could be very useful on the one hand to check the ability of the adopted PMF methodology to identify the preferred binding modes and on the other

ARTICLE

Figure 3. PMF free energy curve of GUA adsorption at the aqueous propylamine interface. The reaction coordinate is the distance between the GUA center of mass and the plane of the frozen silicon atoms. The adsorption geometry corresponding to the most favorable arrangement is displayed below. Silicon atoms are represented by light-yellow spheres. Possible hydrogen bonds are evidenced by dashed black lines.

hand to define an effective investigation strategy that could be applied to larger adsorbates.

’ RESULTS AND DISCUSSION Analysis and Discussion of the Free Energy Profiles and Comparison with the Results of the Short Unrestrained MD Simulations. The PMF profiles along the reaction pathways of

the molecules onto the SAM are shown in Figures 2
5 together with two different views of the configuration of the molecule representing a possible binding arrangement. Density profiles obtained from the short unrestrained MD simulations are also shown (Figures 2
5, blue lines). Even though the total simulation time (5.5 ns) was perhaps not sufficient to explore the whole configurational space of the free molecules, it was long enough to sample a quite large portion of the potential energy surface comprising different minimum energy regions. As it can be noticed in the plots, these regions turned out to be located at the same distance identified through 9148

dx.doi.org/10.1021/jp201103u |J. Phys. Chem. C 2011, 115, 9146–9156

The Journal of Physical Chemistry C

ARTICLE

Figure 4. PMF free energy curve of THY adsorption at the aqueous propylamine interface. The reaction coordinate is the distance between the THY center of mass and the plane of the frozen silicon atoms. The adsorption geometry corresponding to the most favorable arrangement is displayed below. Silicon atoms are represented by light-yellow spheres. Possible hydrogen bonds are evidenced by dashed black lines.

Figure 5. PMF free energy curve of CYT adsorption at the aqueous propylamine interface. The reaction coordinate is the distance between the CYT center of mass and the plane of the frozen silicon atoms. The adsorption geometry corresponding to the most favorable arrangement is displayed below. Silicon atoms are represented by light-yellow spheres. Possible hydrogen bonds are evidenced by dashed black lines.

the PMF profiles. Thus, the density checks might be considered supporting evidence of the reliability of the chosen methodology. Examination of Figure 6, where the locations on the layer explored by the molecules are displayed, reveals that the whole interface was spanned. In the background of each figure, the motion of the surface chains is rendered by means of contour plots colored using fading tones of blue. Denser regions are evidenced by light-cyan areas while black sections characterize holes where no chains are present. These small holes are only a few and are filled with water molecules. The position of the molecule is described, instead, by yellow/brown areas and black/ white contour lines. Dark brown patches indicate more populated regions. Considering that the interface consists of amine moieties only and the chains are quite packed, but mobile, it could be inferred that most of the time the adsorbing molecule sees a very similar landscape wherever it moves. As a consequence, in principle, no specific sector of the simulation cell should be preferred. Moreover, the surface chains reorient

themselves when the molecule reaches the layer and binds on top. The extension of the molecular contour plots clearly indicates that more than half the cell was visited by the DNA bases during the simulation time and different binding sites and positions were sampled. Well-defined minima, which correspond to a position of the center of mass of the fragment at a distance of about 11 Å (z ≈ 22.8 Å) from the frozen layer (z ≈ 11.8 Å), were obtained for all the molecules. The reaction coordinate could have been converted to the distance of the center of mass of the molecule from the average plane of the amine head groups but, considering that the position of the amine plane was continuously changing, the adopted reference was maintained (the position of the minimum corresponds to a separation between the center of mass of the molecule and the average plane defined by the nitrogen atoms of about 2.5 Å). Representative snapshots of the minima were characterized by an orientation of the base rings almost parallel to the surface plane. The most probable inclination was about 15 in all cases 9149

dx.doi.org/10.1021/jp201103u |J. Phys. Chem. C 2011, 115, 9146–9156

The Journal of Physical Chemistry C

ARTICLE

Figure 6. Density maps. In the background, the motion of the surface chains is rendered by means of contour plots colored using fading tones of blue. Denser regions: light-cyan areas, black areas represent holes (no chains). The position of the molecule is depicted by means of yellow/brown areas and black/white contour lines. Brown patches indicate more populated regions: (a) ADE, (b) THY, (c) GUA, and (d) CYT.

(see Figure 7), but the angle distribution was broader for THY and CYT, which could also adopt values in the range 30
60. For these cases, frequent tilting was observed during the whole simulation time, whereas the pyrimidine bases were more stably adsorbed and their distributions were sharper. As suggested by the shapes of the angle distributions, GUA had the tendency to adopt a more definite orientation than ADE probably because more active hydrogen bonding sites were available for the interaction with the surface groups. The evolution of angles and distances as a function of the simulation time is displayed in Figure S1
S4 of the Supporting Information for the four cases. A comparison with the results obtained from the short unrestrained MD simulations is also shown. The data evidence that when the molecules were in close proximity to the surface, both methods (restrained and unrestrained dynamics) found very similar values. The amine moieties of the bases were very often hydrogen bonded to the surface nitrogen atoms, and the N
N distance distributions were centered at about 3.5 Å (Figure 7). Instead, distance analysis for THY, which describes the behavior of the thymine oxygen atom, revealed a preferential location of the carbonyl group at 2.8 Å from the nitrogens of the amine head-groups. As far as the free energy profiles are concerned, they show a very similar trend having no significant barriers along the approaching path of the molecule toward the amine layer and very favorable binding energies. The energy values are very close to each other being about
20 kcal/mol for both the pyrimidine bases and
11 and
14 for thymine and cytosine, respectively. These findings are in line with the number and type of hydrogen bond donor/acceptors and hydrophobic groups of the fragments which determine the characteristic features of the interaction.

Figure 7. Distributions of the angle between the ring of the bases and the surface plane (defined considering the frozen silicon atoms) for all the nucleotide bases. Distributions of the minimum distance between the amine group of the bases and the surface chains head-groups (nitrogen atoms). For THY, one of the oxygens has been considered.

In all of the stable complexes, the amine moiety of the bases was oriented most of the time toward the layer and interacted strongly with the NH2 groups of the SAM forming a network of hydrogen bonds. These interactions effectively tethered the nucleotide bases to the surface and were responsible for the long residence times in that region. Average residence times measured during the free MD trajectories, where desorption events were also observed, were in the ranges 140
280 and 50
70 ps for stronger and weaker adsorption behaviors, respectively. These times could be further increased when, instead of isolated nucleotide bases, a whole DNA chain, covalently linked to a substrate, is considered (Figure 8).36 Indeed, a previous investigation, where the behavior of a short DNA chain tethered to an amine functionalized silicon substrate was investigated, it was observed that, after bending, the DNA filament extended its bases toward the layer and the concerted adsorption of most of the ring hydrogen bond donors and acceptors anchored the whole sequence to the surface. During the adsorption, conformational rearrangements of the various portions of the backbone guided most of the bases toward the surface and direct interactions between the bases and the amine headgroups were established and maintained for tenths of nanoseconds due to cooperative effects.36 9150

dx.doi.org/10.1021/jp201103u |J. Phys. Chem. C 2011, 115, 9146–9156

The Journal of Physical Chemistry C

ARTICLE

Figure 8. DNA fragment (probe) attached to a functionalized surface through its bases. Desorption of a THY nucleotide base during unrestrained 5.5 ns MD simulations at T = 300 K.

Table 1. Interaction Energies (kcal/mol) Calculated during the Constrained MD Simulations Where the Molecule Was in Close Proximity to the Minima of the PMF Plotsa base

surface Elec σ

L-J

σ

Table 2. Interaction Energies (kcal/mol) Calculated during the Short Unrestrained MD Simulations Starting from the Structures Corresponding to the Minima of the PMF Plotsa

water Tot

σ

Elec

σ

L-J

base

σ Total σ

surface Elec σ

L-J

σ

water Tot

σ

Elec

σ

L-J

σ

Tot

σ

A


4.4 3.0
6.4 1.3
10.8 2.4
32.6 6.1
5.7 3.0
38.3 4.6

A


6.5 5.5
6.1 2.2
12.7 5.6
29.8

7.4
5.8 2.9
35.6 6.7

G T


3.6 2.7
7.8 1.1
11.3 2.5
40.5 6.3
5.9 2.9
46.5 5.1
4.5 3.9
6.8 1.3
11.3 3.3
24.2 6.0
6.4 2.4
30.5 4.9

G
5.2 4.6
5.9 2.9
11.1 4.3
38.6 T
3.1 4.2
3.5 2.2
6.6 5.4
27.0

7.2
6.8 3.4
45.4 6.2 6.5
8.4 2.9
35.4 6.2

C


3.6 3.4
5.2 1.2


8.7 3.1
33.0 6.3
5.0 2.7
38.0 5.2

C
1.3 2.8
2.4 2.1
3.7 4.0
35.7
6.1
6.4 3.1
42.1 5.4

a

a

The driving force was essentially of electrostatic nature, and the interaction did not involve the inner portion of the amine chain but just the top head groups. However, as can be noticed in Tables 1 and 2, both the electrostatic and van der Waals components were favorable and well-balanced in all cases. On the contrary, as expected, the interaction with water was always dominated by the electrostatic term. Long Unrestrained MD. Density Profile of Water and Distribution of the Nucleobases near the SAM Interface. Figure 9 shows water and base densities as a function of the z coordinate. The passivating hydrogens, bound to the top silicon layer, are

located at z ≈ 16.7 Å, whereas the amine headgroups span the region extending from z = 20.8 to z = 21.5 Å. As it was expected, a clear water layering took place in proximity to the amine moieties, and two distinct high density zones were identified (first and second peak in the density profiles). The first peak corresponds to those water molecules which formed direct and indirect hydrogen bonds with the SAM (molecules between z = 20.5 Å and z = 23.9 Å). These waters originated from a thick layer which could be only slightly perturbed by the presence of the nucleotide base, as testified by the height variation observed only in the case of cytosine.

The total interaction energy between the adsorbed molecule and the surface chains and its Coulombic and Lennard-Jones components are reported. The nucleotide
solvent interaction energy is also shown.

The total interaction energy between the adsorbed molecule and the surface chains and its Coulombic and Lennard-Jones components are reported. The nucleotide
solvent interaction energy is also shown.

9151

dx.doi.org/10.1021/jp201103u |J. Phys. Chem. C 2011, 115, 9146–9156

The Journal of Physical Chemistry C

ARTICLE

Figure 9. Water density and distribution of the nucleotide bases as a function of the z coordinate.

Indeed, cytosine spent more time than the other nucleotide bases far from the surface, and as a consequence the water layer was slightly denser. The trend of the water density curves between z = 16.7 and z = 20.5 Å suggests also that a few solvent molecules filled the cavities of the SAM and their permanence is confirmed by the intense yellow spots visible in all of the density maps displayed in Figure 10. The first peak is followed by a wider region of lower water density, which extends from z = 23.9 to z = 26.9 Å before converging to the value of the bulk after z = 30.0 Å. As far as the position of the nucleotide base is concerned, broad density peaks are located in the region corresponding to the first water layer, suggesting that all the DNA bases favorably interacted with the surface chains and had the tendency to insert their groups deeply into the functionalized layer if no constraints were imposed on their motion. The four nucleotide bases showed a quite different adsorption behavior which seemed to depend on the number, orientation, and location of the hydrophilic groups. However, base
surface interaction/insertion was not so strong that the hydrogen bonding network between the SAM and water was disrupted. This was confirmed by analyzing the residence time of water in the first layer above the amine chains which was nearly constant in all simulations (22.8 ( 0.5 ps). In fact, the residence time of water was maintained even when the nucleotide base settled on the surface as in the case of ADE and GUA, which showed the longest residence times (in the range 140
300 ps). The shortest residence times were observed in the case of THY (≈ 30 ps) and CYT (≈ 70 ps). Moreover, water residence times agree well with the value obtained from the simulation without the nucleotide bases.56 Surface Morphology. The formation of several dynamic cavities in different regions of the SAM was observed during the whole simulation runs. This was due to the motion and flexibility of the chains and to the presence of the H-passivated silicon sites. The characteristic features of these clefts (which were most of the time filled with water, being too small to host larger molecules) strongly depended on the spatial arrangement and orientation of the chains, which determined not only their size but also the extension of their hydrophobic and hydrophilic portions. However, the dimensions of one of these pockets were such that it could accommodate the nucleotide base inside. Even so, steric criteria are not sufficient to determine the stability of a possible insertion, since also the electrostatic properties of the ligand

Figure 10. Two-dimensional density map from the simulation of ADE (a), GUA (b), THY (c), and CYT (d). In the background, different tones of blue (from navy blue to light cyan) are employed to depict the amine chains motion. Denser regions are represented by lighter colors. Water molecules, characterized by long permanence times, are evidenced by yellow/orange areas, whereas the nucleotide bases are identified by multicolored patches (blue-light blue-cyan-green-yelloworange-red-burgundy-white, increasing density). Low density lines are not displayed for clarity. Persistent contacts are orange-red-burgundywhite, while weaker interactions are blue.

should complement those of the receptor. The inner part of the cavity was formed by the aliphatic sections of the propylamine chains and thus was prevalently hydrophobic, whereas the mouth of the pocket was bordered by the NH2 moieties and as a consequence was more hydrophilic. These considerations prompted us to analyze in deeper detail the molecular electrostatic potential in the space surrounding the four nucleotide bases. Molecular Electrostatic Potential (MEP). The molecular electrostatic potential of each nucleotide base is visualized on its solvent accessible surface and is displayed in Figure 11. Negative and positive potential regions are red and blue, respectively, while fading colors are used to represent decreasing values. The electrostatic effects observed in the region surrounding the molecule were responsible for its recognition and binding to the surface sites. The comparison of the MEP produced by the partial charges of the four nucleotide bases and the examination of the same property inside the big binding pocket (identified during the dynamic runs) revealed three contiguous traits with alternate positive and negative potentials (Figure 11, black arrows). These features were recognized and attracted by the complementary regions located in the interior of the binding cleft. Instead, the different disposition and character of the cytosine MEP features were responsible for its shorter permanence in the cavity and its alternative binding orientation. Cytosine position seemed to be more superficial, and only the NH2 group was inserted deeply into the layer. GUA and CYT oxygen and nitrogen lone pairs produce a noticeable negative potential (red mesh) which is more extended than in the other molecules due to the closeness of the atoms. 9152

dx.doi.org/10.1021/jp201103u |J. Phys. Chem. C 2011, 115, 9146–9156

The Journal of Physical Chemistry C

ARTICLE

Figure 11. Molecular electrostatic potential of ADE, GUA, CYT, and THY visualized on the solvent accessible surface of each nucleotide base using color codes: positive regions (blue), negative regions (red). Black arrows indicate groups which were found inside the binding pocket.

The effect of the lone pair of the ring nitrogens in ADE is similar to those observed in the other cases, but the presence of the positive potential regions between the negative regions somewhat confines the negative portions and prevent their combination. On the contrary, all the hydrogen atoms give rise to positive areas which prevail over the negative ones in THY and CYT because of the cooperative effects. Two-Dimensional Density Maps. Two-dimensional surface density maps, describing the shifts of the amine chains and the motion and location of the nucleotide base and solvent molecules on the SAM were computed and analyzed in detail (Figure 10). From the examination of the maps, it is evident that the chain motion was essentially a rotation of the top head groups around their mean position (Si
C axis) and the surface configuration is very similar and almost conserved in all the simulations. Water molecules were adsorbed and distributed on the surface quite uniformly, and only a few well-localized regions (already observed in the short unrestrained MD), evidenced by the yellow/orange areas, showed higher density. These zones correspond to the holes filled with water. The position of the nucleotide base on the surface was mainly localized inside the cavity, ”appeared” during the dynamics, or near its opening (in the case of CYT, Figure 10d), and was witnessed by the multicolored (blue-light blue-cyan-green-yellow-orange-red-burgundywhite, increasing density) area at the bottom of each map. The presence of more populated regions (orange-red-burgundywhite areas) in the case of ADE and GUA (Figure 10a,b) denotes more persistent contact points and deeper insertion of the nucleotide base in the layer, whereas the reduced area of the high density regions in the case of CYT indicates that lateral diffusion on the surface, exchange between sites, and desorption from the layer were more probable and took place. This behavior was observed for the other cases as well but to a lower extent, as suggested by the examination of the evolution of the z coordinate of selected groups of the bases and the minimum distance of the

Figure 12. Evolution of the minimum distance between the amine nitrogen of the nucleotide bases (or the oxygen atom for THY) and the nitrogen atoms of the SAM. Distance distributions and representative structures.

center of mass of the bases from the surface nitrogens as a function of the simulation time. Base
Surface Distances and Orientation of the Nucleotide Bases. In order to localize the molecule in relation to the surface during the unrestrained dynamics, the minimum distance between the amine nitrogen of the nucleotide base (or the oxygen atom for THY) and the nitrogen atoms of the SAM was calculated and plotted as a function of the simulation time (Figure 12). Moreover, the orientation of the ring plane was checked computing the angle between the normal vectors of the base and of the surface (Figure 13). As it can be observed in Figure 12, ADE distance oscillates around an average value of about 3.94 Å (σ = 0.67 Å), which suggests that ADE maintained this position until the end of the simulation, whereas the other base distances show larger fluctuations implying that desorption and adsorption events took place. The average distances when the molecules were adsorbed on the layer were on average 4.1 Å (σ = 0.7 Å), 3.8 Å (σ = 0.6 Å), and 4.8 Å (σ = 0.8 Å) for CYT, GUA, and THY, respectively, supporting the fact that these groups are most of the time involved in hydrogen bonding interactions with the surface atoms. 9153

dx.doi.org/10.1021/jp201103u |J. Phys. Chem. C 2011, 115, 9146–9156

The Journal of Physical Chemistry C

ARTICLE

Table 3. Interaction Energies (kcal/mol) Calculated during the Long Unrestrained MD Simulationsa base

surface Elec σ

A

L-J

σ

water Tot

σ

Elec

σ

L-J

σ

Tot

σ


5.3 3.3
10.6 1.8
15.9 3.2
22.1 5.8
3.8 2.5
25.9 5.0

G
4.6 3.3
11.6 1.5
16.3 3.4
32.2 5.9
3.0 2.6
35.2 5.0 T C


4.2 2.5
8.3 1.4
12.5 2.5
21.4 5.4
5.2 2.3
26.6 4.3
5.3 1.9
10.5 4.4
15.8 4.5
23.6 7.0
5.4 2.1
30.0 6.8

a

The total interaction energy between the adsorbed molecule and the surface chains and its Coulombic and Lennard-Jones components are reported. The nucleotide
solvent interaction energy is also shown.

Figure 13. Evolution of the angle between the normal vector of the base and the z axis: angle distributions.

To give a more complete picture of the motion of these moieties, on the layer, the plot of the evolution of their z coordinate during the simulation was added to the Supporting Information (Figure S5). The z positions around 21
24 Å represent hydrogen bonds with the amine headgroups, whereas lower values indicate insertion into the SAM (the green line marks the average position of the passivating hydrogens). As it can be noticed, both GUA and CYT, once in contact with the surface chains, had the tendency to insert their NH2 groups inside the layer (into the cavity already described). However, GUA is more stably adsorbed than CYT, and even if it switched the position of the NH2 group from inside to outside the cavity after about 3.8 ns, it remained bound to the SAM through hydrogen bonding interactions until the end of the simulation. On the contrary, CYT remained partially inside the layer for about 2.5 ns, then shifted to a position on top of the headgroups. There it did not find a stable arrangement, as confirmed by the frequent desorption/adsorption events observed in the plot. Indeed, after about 9 ns it migrated to the bulk region of the solvent. As far as the orientation is concerned, ADE and GUA rings formed an angle of about 60 (or 120) with the xy plane, whereas CYT and THY were more mobile and changed their orintations more frequently. Notwithstanding this, some

preferential arrangements were identified with the centers of the broad peaks of the distributions which were located at about 126 and 83 for THY and at 163, 115, 54, and 13 for CYT. Intermolecular Interaction Energy. The interaction energies between the different bases and the surrounding environment were also analyzed and compared with the previous data (Table 3, data refer to positions of the bases close to the interface). As expected, the major contribution to the total base
solvent interaction energy came from the electrostatic term. When the nucleotide bases were surrounded by water (far from the SAM), the total interaction energy was more favorable for GUA (≈ 64 kcal/mol, σ ≈ 6 kcal/mol) followed by ADE (≈ 56 kcal/ mol, σ ≈ 6 kcal/mol), CYT (≈ 45 kcal/mol, σ ≈ 6 kcal/mol), and THY (≈ 38 kcal/mol, σ ≈ 6 kcal/mol). Upon binding, a reduction of the interaction with water ranging from 30% to 50% (molecule inside the cavity) was observed. Unlike the other simulations, the van der Waals component dominates the interaction of the nucleotide base with the layer, in all cases. The presence of a strong L-J interaction could be ascribed to the different binding mode identified during these longer simulations where the molecule was prevalently located inside the quite large hydrophobic binding pocket, which was potentially present in the SAM at the beginning of the simulation but was not evident and inaccessible. From the examination of the sampled data, and in particular the two-dimensional density maps where ligand isodensity lines are found on the cavity position, it could be speculated that during both the restrained and short unrestrained MD runs the nucleotide bases explored the region of the pocket only superficially and did not spend enough time there to gain access to the inner portion. The pocket (central bottom black region of each map) was partially covered by three amine head groups which prevented the insertion of external moieties. Apparently, during the longer unrestrained MD runs a combination of events such as closeness of specific groups of the ligand and concerted movements of the nearby chains caused opening of the “lid” and intercalation of the nucleotide base. As already mentioned before, the interior surfaces of the cleft were essentially formed by the CH2 groups of the propylamines, thus being hydrophobic, whereas the top section, in contact with the solvent, was more hydrophilic because of the neighboring NH2 groups. As a consequence, the CH3 groups of the ligands entered the core of the cleft whereas its hydrogen bond donors and acceptors remained in contact with the hydrophilic moieties at the opening and water molecules. 9154

dx.doi.org/10.1021/jp201103u |J. Phys. Chem. C 2011, 115, 9146–9156

The Journal of Physical Chemistry C

’ CONCLUSIONS A specific and accurate force field, parametrized in a previous work30 to describe the most probable packing/coverage of solvated alkyl-amine chains attached to a model silicon surface (Si(111)), was employed here to investigate the binding properties of a similar layer interacting with nucleotide residues. Both restrained and unrestrained all-atom molecular dynamics simulations were carried out to characterize the adsorption of the fragments. The free energy change upon adsorption was calculated by means of the potential of the mean constraint force method and several probable binding configurations were identified. All the PMF profiles were attractive and lacked barriers along the approaching paths. From the selected structures, unrestrained MD simulations were run to explore the regions of the minima. The position of the DNA bases and their local interactions were analyzed and characterized in detail. All the fragments had favorable free energy of binding to the amine layer. The strength of the adsorption was determined by the number of contacts between the functional groups of the molecules and the head-groups of the chains. Sometimes, stable arrangements were mediated by solvent molecules close to the interface. The orientation of the rings was such as to maximize the number of hydrogen bonds, to insert a portion of the nuocleobase inside the layer, and to minimize unfavorable interactions. These data were consistent and satisfactorily agreed with the results of investigations on tethered DNA dynamics36 where the bases were bonded to the DNA backbone through their N1 or N9 atoms. In this study, the DNA backbone was replaced with a methyl group and no restrictions were applied either to the motion or position of the base. Thus, the molecules could adopt unrealistic orientations which were hindered in the more complex cases analyzed previously. Surface induced attraction arose from both propylamine headgroups and the stabilization of water layers close to the interface. The long simulation runs revealed that the amine layer is a highly dynamic substrate, which can rearrange its overall configuration and present a new face to the ligand. These morphological changes could be induced by environmental effects, by the solvent influence, and by the action of the nucleotide base itself, but they could not be easily predicted due to the complex nature of the systems, to the huge number of their components, and to the inability to obtain reliable and accurate experimental measurements at the atomic/molecular level. The tendency of the nucleotide bases to reduce the degree of exposure of their nonpolar groups to water was evidenced by the preferential location of these moieties inside the hydrophobic sections of surface cavities. However, these arrangements were highly improbable in linked systems. The properties of the interface and, at the same time, the complexity of the ligand render the use of PMF methodology, based on the definition of appropriate reaction coordinates, quite difficult and elaborate, nevertheless, a useful tool to understand, in combination with other techniques, the principles that control molecular adsorption on functionalized surfaces. The force driving the fragment to the layer was found to be of electrostatic nature, but the binding at the interface was a balanced combination of both electrostatic and van der Waals components. All the molecules were attached with a similar strength even though a slight preference for pyrimidine bases was observed. ’ ASSOCIATED CONTENT

bS

Supporting Information. Evolution of the angle between the ring plane and the plane defined by the frozen silicon atoms

ARTICLE

of the layer; evolution of the minimum distance between selected atoms of the four bases and the nitrogens of the amine chains; evolution of the z coordinate of the nitrogen atoms of the amine groups of ADE, CYT, and GUA and the oxygen atom of THY during the simulation. This material is available free of charge via the Internet at http://pubs.acs.org.

’ AUTHOR INFORMATION Corresponding Author

*Phone: þ39-050-3152520. Fax: þ39-050-3152555. E-mail: [email protected].

’ ACKNOWLEDGMENT The authors gratefully acknowledge the computing facilities of CASPUR (Inter-University Consortium for the Application of Super-Computing for Universities and Research-Rome). ’ REFERENCES (1) Healey, B. G.; Matson, R. S.; Walt, D. R. Anal. Biochem. 1997, 251, 270–279. (2) Palecek, E.; Fojta, M.; Tomschik, M.; Wang, J. Biosens. Bioelectron. 1998, 13, 621–628. (3) Liu, J.; Cao, Z.; Lu, Y. Chem. Rev. 2009, 109, 1948–1998. (4) Zheng, G.; Patolsky, F.; Cui, Y.; Wang, W.; Lieber, C. M. Nat. Biotechnol. 2005, 23, 1294–1301. (5) Drummond, T. G.; Hill, M. G.; Barton, J. K. Nat. Biotechnol. 2003 21, 1192–1199. (6) Yang, C.; Zhong, Z.; Lieber, C. M. Science 2005, 310, 1304–1307. (7) Levicky, R.; Horgan, A. Trends Biotechnol. 2005, 23, 143–149. (8) Binder, H. J. Phys.: Condens. Matter 2006, 18, S491–S523. (9) Tarlov, M. J.; Steel, A. B. DNA-based Sensors. In Biomolecular Films: Design, Functions and Applications; Marcel Dekker: New York, 2003. (10) Ramsey, G. Nat. Biotechnol. 1998, 16, 40–44. (11) Kitano, H. Science 2002, 295, 1662–1664. (12) Ziegler, C. Anal. Bioanal. Chem. 2004, 379, 946. (13) Patolosky, F.; Lieber, C. M. Mater. Today 2005, 8, 20. (14) Hood, L.; Heath, J. R.; Phelps, M. E.; Lin, B. Science 2004, 306, 640–643. (15) Draghici, S.; Khatrri, P.; Eklund, A. C.; Szallasi, Z. Trends Genet. 2006, 22, 101–109. (16) Cummings, C. A.; Relman, D. A. Emerg. Infect. Dis. 2000, 6, 513. (17) Ruiz-Hitzky, E. Chem. Rec. 2003, 3, 88–100. (18) Strother, T.; Cai, W.; Zhao, X.; Hamers, R. J.; Smith, L. M. J. Am. Chem. Soc. 2000, 122, 1205–1209. (19) Hansen, M. N.; Farjami, E.; Kristiansen, M.; Clima, L.; Pedersen, S. U.; Daasbjerg, K.; Ferapontova, E. E.; Gothelf, K. V. J. Org. Chem. 2010, 75, 2474–2481. (20) Johansson, E.; Hurley, P. T.; Brunschwig, B. S.; Lewis, N. S. J. Phys. Chem. C 2009, 113, 15239–15245. (21) Yu, H.; Webb, L. J.; Heath, J. R.; Lewis, N. S. Appl. Phys. Lett. 2006, 88, 252111. (22) Yu, H.; Webb, L. J.; Solares, S. D.; Cao, P.; Goddard, W. A.; Heath, J. R.; Lewis, N. S. J. Phys. Chem. B 2006, 110, 23898–23903. (23) Linford, M. R.; Fender, P.; Eisenberger, P. M.; Chisey, C. E. D. J. Am. Chem. Soc. 1995, 117, 3145–3155. (24) Sieval, A. B.; Vleeming, V.; Zuilhof, H.; Sudh€olter, E. J. R. Langmuir 1999, 15, 8288–8291. (25) Zhu, X. Y.; Boiadjiev, V.; Mulder, J. A.; Hsung, R. P.; Major, R. C. Langmuir 2000, 16, 6766–6772. (26) Sieval, A. B.; van der Hout, B.; Zuilhof, H.; Sudh€olter, E. J. R. Langmuir 2000, 16, 2987–2990. (27) Sieval, A. B.; van der Hout, B.; Zuilhof, H.; Sudh€olter, E. J. R. Langmuir 2001, 17, 2172–2181. 9155

dx.doi.org/10.1021/jp201103u |J. Phys. Chem. C 2011, 115, 9146–9156

The Journal of Physical Chemistry C

ARTICLE

(28) Zhang, L.; Wesley, K.; Jiang, S. Langmuir 2001, 17, 6275–6281. (29) Peig, Y.; Ma, J.; Jiang, Y. Langmuir 2003, 19, 7652–7661. (30) Barone, V.; Cacelli, I.; Ferretti, A.; Monti, S.; Prampolini, G. Phys. Chem. Chem. Phys. 2010, 12, 4201–4209. (31) Vainrub, A.; Pettitt, B. M. Biopolymers 2003, 68, 265. (32) Wong, K. Y.; Vainrub, A.; Powdrill, T.; Hogan, M.; Pettitt, B. M. Mol. Simul. 2004, 30, 121–129. (33) Ambia-Garrido, J.; Pettitt, B. M. Commun. Comput. Phys. 2008, 3, 1117–1131. (34) Wong, K. Y.; Vainrub, A.; Powdrill, T.; Hogan, M.; Pettitt, B. M. (35) Wong, K. Y.; Pettitt, B. M. Biopolymers 2004, 73, 570–578. (36) Barone, V.; Cacelli, I.; Ferretti, A.; Monti, S.; Prampolini, G. Phys. Chem. Chem. Phys. 2009, 11, 10644–10656. (37) Monti, S.; Cacelli, I.; Ferretti, A.; Prampolini, G.; Barone, V. J. Phys. Chem. B 2010, 114, 8341–8349. (38) Monti, S.; Cacelli, I.; Ferretti, A.; Prampolini, G.; Barone, V. Mol. BioSyst. 2010, 6, 2230–2240. (39) Pranata, J.; Wierschke, S. G.; Jorgensen, W. L. J. Am. Chem. Soc. 1991, 113, 2810–2819. (40) Miller, J. L.; Kollman, P. A. J. Phys. Chem. 1996, 100, 8587–8594. (41) Jurecka, P.; Hobza, P. J. Am. Chem. Soc. 2003, 125, 15608–15613. (42) Monti, S.; Walsh, T. R. J. Phys. Chem. C 2010, 114, 22197–22206. (43) Cacelli, I.; Prampolini, G. J. Chem. Theory Comput. 2007, 3, 1803. (44) Jorgensen, W. L. J. Am. Chem. Soc. 1981, 103, 335–350. (45) Nose, S. Mol. Phys. 1984, 52, 255. (46) Hoover, W. G. Phys. Rev. A 1985, 31, 1695. (47) Darden, T.; York, D.; Pedersen, L. J. Chem. Phys. 1993, 98, 10089. (48) Essmann, U.; Perera, L.; Berkowitz, M.; Darden, A.; Lee, H.; Pedersen, L. J. Chem. Phys. 1995, 103, 8577. (49) Hess, B.; Kutzner, C.; van der Spoel, D.; Lindahl, E. J. Chem. Theory Comput. 2008, 4, 435. (50) Skelton, A. A.; Walsh, T. R. Mol. Simul. 2007, 33, 379. (51) Skelton, A. A. Modelling the Interaction Between Titania Surfaces and Strong-Binding Peptides. Thesis, University of Warwick, Coventry, U.K., 2008. (52) Skelton, A. A.; Liang, T.; Walsh, T. R. ACS Appl. Mater. Interfaces 2009, 1, 1482. (53) Kerisit, S.; Parker, S. C. J. Am. Chem. Soc. 2004, 126, 10152–10161. (54) Marrink, S.-J.; Berendsen, H. J. C. J. Phys. Chem. 1994, 98, 4155– 4168. (55) Hess, B. J. Chem. Phys. 2002, 116, 209–217. (56) Barone, V.; Cacelli, I.; Ferretti, A.; Monti, S.; Prampolini, G. J. Phys. Chem. C 2011, 115, 4145–4154.

9156

dx.doi.org/10.1021/jp201103u |J. Phys. Chem. C 2011, 115, 9146–9156