INTRODUCTION

TOP ABSTRACT INTRODUCTION METHODS RESULTS AND DISCUSSION CONCLUSIONS REFERENCES

For many biochemical reactions, such as initiation of transcription and recombination, two distant DNA sites should be brought into close contact (see, for example, Rippe et al., 1995; Wasserman and Cozzarelli, 1986). This site juxtaposition is a random event caused by thermal fluctuations of the DNA global conformation. The purpose of this paper is to model this process numerically for a supercoiled plasmid in a dilute solution.

Consider an irreversible chemical reaction between two reactive groups attached to different sites of a supercoiled DNA molecule. We assume that the activity of the reactive groups is "switched on" at the initial time instant t = 0, with the DNA chain being in statistical equilibrium. Due to the internal diffusion, the reactive groups repeatedly collide with each other until the reaction takes place. The reaction radius R is the distance between the centers of the reactive groups at the instant of collision. The kinetics of the reaction are characterized by the probability distribution P(t) of the reaction time, i.e., the time corresponding to the final, "successful," collision.

In the present study we will consider only diffusion-controlled reactions. In such reactions the mean time between the first and the final collisions is negligible in comparison with the mean time of the first collision.

The general analytical theory of intrachain reactions of polymers was developed by Wilemski and Fixman (1974a). These and other authors applied it to various polymer systems (Wilemski and Fixman, 1974b; Doi, 1975a, b; de Gennes, 1982a, b; Noolandi et al., 1984), including non-supercoiled DNA (Berg, 1984). For supercoiled DNA, however, no fruitful analytical approach has been found yet. Several simple theoretical models were proposed to account for the quasi-one-dimensional "slithering" motion of the strands forming a superhelix (Sessions et al., 1997; Marko, 1997; Wedemann et al., 1998), but these models completely ignore the real three-dimensional organization of the DNA molecule. The most adequate tool for modeling intrachain reactions of supercoiled DNA in a general sense is a numerical simulation using, e.g., the Brownian dynamics (BD) technique.

Recently, a number of computer algorithms have been developed for simulations of supercoiled DNA by the BD method (Chirico and Langowski, 1994, 1996; Heath et al., 1996; Klenin et al., 1998; Jian et al., 1998). Jian et al. (1998) were the first to apply this approach for calculation of the mean time of the first collision, _F, for various DNA lengths L, superhelical densities , and distances S between the reactive groups along the DNA contour. The typical _F value was 1 ms (L = 1500 bp,  = 0.05, S = 600 bp). However, these simulations were performed only for an NaCl concentration of 0.01 M and a reaction radius of 10 nm. It should be noted that, at the given ionic conditions, the effective diameter of the double helix is d_eff = 15 nm (Stigter, 1977). Thus, R was considerably smaller than d_eff. In this case, the first collision time should be strongly affected by the necessity of diffusion against the electrostatic potential. This restriction vanishes at higher ionic strengths, which are more typical for biochemical reactions.

In the present study we use our earlier BD program (Klenin et al., 1998) to calculate the probability distribution P_F(t) of the first collision time for supercoiled DNA in a 0.1 M NaCl solution. We concentrate our attention on the general form of the function P_F(t) and the dependence of the mean time of the first collision, _F, on the reaction radius R. In our calculations, L = 2500 bp, S = 470 bp, and  = 0.05.

It is well known that the intrachain interactions in DNA can efficiently be mediated by a DNA-bending protein or by a curved sequence of the basepairs (see, for example, Rippe et al., 1995). To model this effect quantitatively for diffusion-controlled reactions, we repeat our calculations for a DNA chain with a permanent bend of 100° in the middle position between the reactive groups along the DNA contour. As expected, the presence of the permanent bend dramatically reduces the mean time of the first collision, _F, and changes the form of the distribution P_F(t).

	METHODS

TOP ABSTRACT INTRODUCTION METHODS RESULTS AND DISCUSSION CONCLUSIONS REFERENCES

The probability distribution P_F(t) of the first collision time for two sites of a supercoiled DNA was calculated by the BD method.

Our BD algorithm is described in detail elsewhere (Klenin et al., 1998). Here, we present only the principal points. The DNA molecule was modeled by a closed chain of beads connected through straight elastic segments. To each segment, a local reference frame was attached. The energy of the system was given by harmonic potentials with respect to 1) the angles between consecutive segments, 2) the twist angles between consecutive reference frames, and 3) the deviations from the equilibrium segment length. Electrostatic interactions were taken into account through a Debye-Hückel potential. The linear charge density of DNA was renormalized as described by Stigter (1977). Each step in a BD trajectory was performed according to the algorithm of Ermak and McCammon (1978), with second-order corrections. The hydrodynamic interaction between the beads was given by the Rotne-Prager tensor (Rotne and Prager, 1969). Initial chain conformations were obtained by a Monte Carlo (MC) method.

In the MC simulations, we followed the algorithm proposed by Vologodskii et al. (1992), with the modifications described in Klenin et al. (1995) and Klenin and Langowski (2000). The geometry and the energy of the MC model were the same as those for the BD model, with the only exception that the segment length was fixed. In addition, knotted conformations were explicitly forbidden. Two types of MC steps were used: 1) a pivoting step, by which a subchain bounded by two randomly chosen beads was rotated about its end-to-end vector by a random angle, and 2) a reptational step, by which a randomly chosen subchain of n + 1 segments was exchanged with a randomly chosen subchain of n segments, the end-to-end lengths being properly readjusted. (In the present study, n = 3.) The pivoting steps were accepted or rejected in accordance with the classical Metropolis criterion. For the reptational steps, the criterion was modified as described in Klenin and Langowski (2000).

The following set of parameters was used. The molecule length was L = 850 nm (2500 bp). The hydrodynamic radius of DNA was r_h = 1.2 nm (Hagerman and Zimm, 1981). The superhelical density was  = 0.05, the persistence length L_p = 50 nm, the torsion rigidity C = 2.5 × 10¹⁹ erg cm, the NaCl concentration I = 0.1 M, the temperature T = 293 K. The elastic stretch modulus was "softened" to the value _S = 63 pN to provide reasonable computational time. The equilibrium segment length was l₀ = 10 nm, the BD time step t = 1.9 ns. The parameters correspond to a bead radius r_b = 2.3 nm. The distance between the reactive groups along the DNA contour was S = 160 nm (470 bp).

The scheme of the simulations was as follows. First, we generated a sufficiently large number N_start of independent chain conformations (N_start  100) by an MC method. From each such conformation, a BD trajectory was initiated. The segment lengths were first allowed to relax for 100 BD steps, after which the simulation was started. For each subchain of the length S, the end-to-end distance was monitored as a function of time. These functions were then used to obtain histograms of the first collision time for a set of reaction radii in the range between 6 and 20 nm.

In total, we performed two series of simulations with the set of parameters given above. In the second series, the chain contained a permanent bend of 100° in the center of the shortest stretch between the reactive groups. Consequently, only one pair of beads was regarded as reactive groups. N_start was 2000. The permanent bend was realized as a sequence of three "bent" joints, with the equilibrium bending angles lying in the same plane (Klenin et al., 1995, 1998). At each "bent" joint, the equilibrium angle between the adjacent segments was equal to 33.3°.

The following limitations of our model should be mentioned. First, the beads representing the reactive groups have the same hydrodynamic radius as all the other beads. Second, we have not achieved a plateau in the dependence of the results on the segment length l₀. Estimations show that a twofold decrease in l₀ leads to an 10% decrease in the mean time of the first collision for R = 6 nm. However, elimination of these limitations would cost too much CPU time and is beyond the scope of the present study.

Overall, our simulations took 1 year CPU time on the HP-S2000 installation of the DKFZ. We have verified that our algorithm can reproduce the results reported by Jian et al. (1998).

	RESULTS AND DISCUSSION

TOP ABSTRACT INTRODUCTION METHODS RESULTS AND DISCUSSION CONCLUSIONS REFERENCES

Our calculations are performed for supercoiled DNA of 2500 bp length in aqueous solution at a NaCl concentration of 0.1 M. The reactive groups are separated by 470 bp along the DNA contour. It should be noted that the effective diameter of the double helix at the given ionic strength is 5.6 nm (Stigter, 1977). The mean superhelix diameter under these conditions, as calculated by MC simulations and measured by small angle neutron scattering, is 9 nm (Hammermann et al., 1998).

At first, we present the results for an isotropic chain model. Fig. 1 shows the probability distribution P_F(t) of the first collision time for different reaction radii R. The function P_F(t) is essentially non-exponential. Within the initial time interval, P_F(t) scales approximately as t^1/2. During this period the most part of the first collisions take place. This is demonstrated in Fig. 2, where the probability (t) that the first collision occurs before the time t is shown:
Fig. 3 presents the mean time of the first collision, _F, as a function of R. The value of _F is inversely proportional to R. Note that we consider only those R values that are larger than the effective diameter of the double helix (5.6 nm).

View larger version (27K): [in this window] [in a new window]   FIGURE 1   Distribution function PF(t) of the first collision time for various reaction radii R. Isotropic chain. The interpolation curves were calculated according to Eq. 9 with the parameter F as presented in Fig. 3. For convenience, the data are multiplied by an arbitrary factor z.

View larger version (23K): [in this window] [in a new window]   FIGURE 2   Probability (t) that the first collision occurs before the time t for various reaction radii R. Isotropic chain. The interpolation curves were obtained by integration of the corresponding curves from Fig. 1.

View larger version (12K): [in this window] [in a new window]   FIGURE 3   Mean time of the first collision, F, as a function of the reaction radius R. Isotropic chain. The interpolation curve is the best fit of the form F ~ R1.

The function P_F(t) changes dramatically when a permanent bend of 100° is inserted between the reaction groups (Figs. 4-6). The time dependence of P_F is significantly stronger: for R  10 nm, P_F(t) scales almost as t¹ (10⁶ s < t < 10³ s). The mean time of the first collision, _F, is approximately by one order of magnitude smaller than in the case of the isotropic chain. It should be noted that the effect of the permanent bend would be opposite if the reactive groups were located non-symmetrically relative to the bend.

View larger version (33K): [in this window] [in a new window]   FIGURE 4   Distribution function PF(t) of the first collision time for various reaction radii R. Chain with a permanent bend of 100° between the reactive groups. For convenience, the data are multiplied by an arbitrary factor z.

View larger version (26K): [in this window] [in a new window]   FIGURE 5   Probability (t) that the first collision occurs before time t for various reaction radii R. Chain with a permanent bend of 100° between the reactive groups.

View larger version (12K): [in this window] [in a new window]   FIGURE 6   Mean time of the first collision, F, as a function of the reaction radius R. Chain with a permanent bend of 100° between the reactive groups.

In a supercoiled DNA molecule, there are two possible mechanisms of collisions between reactive groups, as illustrated in Fig. 7. One can distinguish the collisions by relative reptation of the two strands forming the superhelix and the collisions by bending of the superhelix axis (Marko and Siggia, 1995; Marko, 1997). The process of reptation is known to be very slow. According to theoretical estimations (Marko, 1997), the mean time of the first collision resulting from reptation, _F^(r), depends on the contour separation of the reactive groups, S, as _F^(r) ~ S², whereas, for the bending motion of an unbranched molecule, _F^(b) ~ S^3/2 (Berg, 1984). Branching reduces the mean distance between the reactive groups along the superhelix axis and makes the dependence of _F^(b) on S still weaker. For sufficiently large S, one can therefore expect that _F^(r) _F^(b). Even though in our BD simulations the value of S is relatively small (470 bp), visual inspection of the chain conformations shows that both mechanisms contribute to the simulated function P_F(t). In the case of the isotropic chain, the part of the collisions by bending is ~20% for R = 6 nm, and 65% for R = 20 nm. The permanent bend freezes the reptational motion in a favorable position and reduces the part of the collisions by bending to 3% for R = 6 nm, and 10% for R = 20 nm.

View larger version (24K): [in this window] [in a new window]   FIGURE 7   Various types of collisions. (A) Collision by mutual reptation of the DNA strands; (B) collision by bending of the superhelix axis; (C) collision by bending of the superhelix axis in the region of an end loop. The positions of the reactive groups are indicated by the arrows. The chain conformations are taken from the BD simulations. The angles between segments are smoothed. The chain thickness corresponds to the effective diameter of DNA (5.6 nm).

The reptational motion of a short unbranched DNA molecule can be described by a simple analytical model. One can imagine the superhelix end loops as two "particles" diffusing along the DNA contour in a cyclic one-dimensional space. The distance a between the particles is fixed and is equal to the half DNA length, L/2. A collision takes place whenever one of the particles passes through the point x = 0 lying exactly in the middle between the reactive groups. This system is equivalent to a single particle that diffuses along an infinite straight line and, at the initial time instant, is uniformly distributed in the interval (0, a). The boundaries are transparent. Passing through a boundary corresponds to a collision. The problem is to find the distribution function P_F(t) of the first passage time.

Let P_F(t|x₀) be the distribution of the first passage time for a particle starting from the point x₀ (0 < x₀ < a). Then the probability that such a particle is within the distance dx/2 to one of the boundaries at the time t is given by

(1)

where the function G_B(t|0)dx under the integration sign is the probability that, at the time t, the particle can be found within the distance dx/2 to one of the boundaries, provided it has started from one of the boundaries. Note that G_B(t|0) = G_B(t|a). Averaging Eq. 1 over x₀, we get

(2)

This equation can be solved with respect to P_F(t) by means of the Laplace transform, which for an arbitrary function F(t) is defined by

(3)

The Laplace transform converts the convolution integral in the right side of Eq. 2 to a product, _B(s|0)_F(s), and Eq. 2 can be rewritten as

(4)

For the function G_B(t|x₀) we have

(5)

where

(6)

is the usual diffusion propagator, with D being the mean diffusion coefficient of an end loop. After the substitution of Eqs. 5 and 6 into Eq. 4, the latter becomes

(7)

with the mean time of the first collision

(8)

The representation in the form of an infinite series in Eq. 7 makes it possible to take the inverse Laplace transform of _F(s):

(9)

For the times of interest, the series in Eq. 9 converges very rapidly. One needs, practically, to take into account only a few first terms.

In our BD simulations, the DNA is relatively short: L = 2500 bp. Approximately 80% of the molecules are not branched. As shown in Fig. 1, Eq. 9 predicts the form of the simulated function P_F(t) for the isotropic chain reasonably well. According to Eq. 8, the mean diffusion coefficient of an end loop, D, is of the order 10⁷ cm²/s or 100 kbp²/s, in agreement with earlier direct calculations (Chirico and Langowski, 1996; Wedemann et al., 1998).

We have, however, to answer the question: why does Eq. 9 still hold for R = 20 nm, when the collisions by bending prevail? By inspecting the chain conformations, we noted that approximately one-half of the collisions by bending occur near the end loops (Fig. 7 c), where we observed that the superhelix is most flexible. We conclude that the probabilities of collisions of both types are strongly correlated. This fact explains why the form of the function P_F(t) can be predicted by a one-dimensional model, whereas the dependence of _F on R remains essentially "three-dimensional."

Another question (to which, at present, we have no answer), is whether the motion of an end loop along the DNA contour can be described by the usual diffusion propagator (Eq. 6). For this motion, the term "slithering" is commonly used; however, slithering, in the literal meaning of this word, was never observed in BD simulations. One reports rather global reshaping of the molecule through the formation and disappearance of the end loops (Chirico and Langowski, 1996; Jian et al., 1998). We made the similar observation: any substantial displacement of the two opposite loops is always accompanied by the "birth" and the "death" of a superhelix branch. Branching is necessary for the loop migration, although the percentage of the branched conformations might be relatively low. In short DNA molecules (L  1500 bp), Wedemann et al. (1998) observed processes which they called "slithering"; however, upon closer inspection, also in that case the reshaping takes place by the creation and disappearance of small loop-like deformations. The end loops seem to interact with each other not by slithering, but by exchanging small loops that do not form real branches (Fig. 8). The efficiency of such interactions should rapidly decrease with the distance along the superhelix axis. If the distance exceeds a certain threshold, it is more probable for two small loops to collide and form a separate branch than to migrate one after another in the same direction. Except for a very short DNA, this threshold is apparently smaller than the mean "end-to-end distance." The applicability of the usual diffusion propagator (Eq. 6) for such a system remains to be investigated.

View larger version (19K): [in this window] [in a new window]   FIGURE 8   Scheme of interaction between two end loops. The "quantum" of the interaction is a small loop, which detaches itself from one of the large end loops, propagates along the superhelix axis, and joins the other large end loop. No slithering is required for this process. This kind of interaction is only possible when the distance between the end loops along the superhelix axis is sufficiently small. Otherwise, the "quanta" collide with each other and form a separate branch.

	CONCLUSIONS

TOP ABSTRACT INTRODUCTION METHODS RESULTS AND DISCUSSION CONCLUSIONS REFERENCES

For moderate distances between the reactive groups, kinetics of a diffusion-controlled intrachain reaction of supercoiled DNA combines the features of one- and three-dimensional diffusion processes. On the one hand, the probability distribution P_F(t) of the first collision time scales as P_F ~ t^1/2 in a large initial time interval; on the other hand, the mean time of the first collision is inversely proportional to the reaction radius. The collisions are coursed by the two types of internal motion: 1) quasi-one-dimensional mutual reptation of the DNA strands, and 2) three-dimensional bending of the superhelix axis. The two types of motion are strongly correlated, because bending is most probable in the regions of the end loops.

A permanent bend of 100° in the middle position between the reactive groups dramatically changes the form of the distribution P_F(t) and reduces the mean time of the first collision by approximately one order of magnitude.