To calculate global elasticity, dnaMD globalElasticity command can be used. However, it takes HDF5 file as input.
Following steps can be performed to generate HDF5 files. The tutorial file can be downloaded here. We will prepare HDF5
file for both free and bound DNA.
Both free and bound DNA is superimposed on to the same DNA structure. Careful that bending calculation is fitting dependent. Therefore, at first we aligned both free and bound DNA to a common DNA structure as follows
In [ ]:
%%bash
# Align bound DNA
echo 6 0 | gmx trjconv -f inputs/F1_complex_DNA.xtc -s inputs/dna.tpr -n inputs/dna.ndx -o complex_dna_aligned.xtc -fit rot+trans
# Align free DNA
echo 6 0 | gmx trjconv -f inputs/F1_free_DNA.xtc -s inputs/dna.tpr -n inputs/dna.ndx -o free_dna_aligned.xtc -fit rot+trans
1. Run do_x3dna on DNA trajectory, -nofit is used because DNA is already superimposed to a common DNA structure using trjconv.
In [ ]:
%%bash
# For free DNA
echo 2 | do_x3dna -f free_dna_aligned.xtc -s inputs/dna.tpr -n inputs/dna.ndx -ref -noavg -nofit -name free
mv *_free.dat outputs/.
# For bound DNA
echo 2 | do_x3dna -f complex_dna_aligned.xtc -s inputs/dna.tpr -n inputs/dna.ndx -ref -noavg -nofit -name bound
mv *_bound.dat outputs/.
2. Run dnaMD to extract the parameters from do_x3dna output files and save as HDF5 file. Also calculate global axis, curvature and tangents.
In [8]:
%%bash
# For free DNA
dnaMD saveAsH5 -tbp 27 -i outputs/L-BPS_free.dat,outputs/L-BPH_free.dat,outputs/HelAxis_free.dat -o free_dna.h5
dnaMD axisCurv -tbp 27 -bs 2 -be 25 -ctan -scp 100 -s 1000 -cta 30 -io free_dna.h5 -ap free_dna_axis.pdb
# For bound DNA
dnaMD saveAsH5 -tbp 27 -i outputs/L-BPS_bound.dat,outputs/L-BPH_bound.dat,outputs/HelAxis_bound.dat -o bound_dna.h5
dnaMD axisCurv -tbp 27 -bs 2 -be 25 -ctan -scp 100 -s 1000 -cta 30 -io bound_dna.h5 -ap bound_dna_axis.pdb
Now, we have HDF5 files of both free and bounds DNA. It can be used for the calculation of elastic properties. These files
can be used with either dnaMD Python module or dnaMD globalElasticity.
Following command calculate Bending Stretching Twist modulus matrix. Output matrix will be stored in csv file. Elastic modulus matrix is printed as output and average values of contour length and cumulative twist angle is also printed.
In [1]:
%%bash
dnaMD globalElasticity -i free_dna.h5 -tbp 27 -bs 4 -be 20 -estype BST -paxis X -o elastic_modulus_BST.csv
The above modulus matrix is in this form:
$$\text{modulus matrix} = \begin{bmatrix} M_{Bx} & M_{Bx,By} & M_{Bx,S} & M_{Bx,T} \\ M_{Bx,By} & M_{By} & M_{By,S} & M_{By,T} \\ M_{Bx,S} & M_{By,S} & M_{S} & M_{S,T} \\ M_{Bx,T} & M_{Bx,T} & M_{S,T} & M_{T} \end{bmatrix} $$Where:
Following command calculate Bending Stretching Twist modulus matrix. Output matrix will be stored in csv file. Elastic modulus matrix is printed as output and average values of contour length and cumulative twist angle is also printed.
In [2]:
%%bash
dnaMD globalElasticity -i free_dna.h5 -tbp 27 -bs 4 -be 20 -estype ST -o elastic_modulus_ST.csv
The above modulus matrix is in this form:
$$\text{modulus matrix} = \begin{bmatrix} M_{S} & M_{S,T} \\ M_{S,T} & M_{T} \end{bmatrix} $$Where:
Same command can be used to calculate elasticity as a function of time with option
-ot/--output-time and save it in csv format file. This result can beused to check
their convergence.
-fgap/--frame-gap is an essential option.NOTE: Elastic properties cannot be calculated using a single frame because fluctuations are required. Therefore, here time means trajectory between zero time to given time.
In [3]:
%%bash
dnaMD globalElasticity -i free_dna.h5 -tbp 27 -bs 4 -be 20 -estype BST -paxis X -fgap 100 -em block -ot modulus_time_BST.csv
# Print first and last 3 line of output file
echo "====================================="
echo "Elastic modulus as a function of time"
echo "====================================="
head -4 modulus_time_BST.csv
printf ".\n.\n.\n"
tail -3 modulus_time_BST.csv
Some of the plots from above data can be found here
Same as above but only for stretching and twisting motions.
In [4]:
%%bash
dnaMD globalElasticity -i free_dna.h5 -tbp 27 -bs 4 -be 20 -estype ST -paxis X -fgap 100 -gt "gmx analyze" -em "block" -ot modulus_time_ST.csv
# Print first and last 3 line of output file
echo "====================================="
echo "Elastic modulus as a function of time"
echo "====================================="
head -4 modulus_time_ST.csv
printf ".\n.\n.\n"
tail -3 modulus_time_ST.csv
To caluclate global deformation enrgy, dnaMD globalEnergy can be used. At first, elastic matrix from reference
DNA (most often free or unbound DNA) is calculated and subsequently this matrix is used to calculate deformation free
energy of probe DNA (most often bound DNA).
The deformation free energy is calculated using elastic matrix as follows
$$G = \frac{1}{2L_0}\mathbf{xKx^T}$$$$\mathbf{x} = \begin{bmatrix} (\theta^{x} - \theta^{x}_0) & (\theta^{y} - \theta^{y}_0) & (L - L_0) & (\phi - \phi_0) \end{bmatrix}$$Where, $\mathbf{K}$, $\theta^{x}_0$, $\theta^{y}_0$, $L_0$ and $\phi_0$ is calculated from reference DNA while $\theta^{x}$, $\theta^{y}$, $L$ and $\phi$ is calculated for probe DNA from each frame.
This command gives output energy as a function of time in csv file and also average energies with error.
In [5]:
%%bash
dnaMD globalEnergy -ir free_dna.h5 -ip bound_dna.h5 -tbp 27 -bs 4 -be 20 -estype BST -et "all" -paxis X -gt "gmx analyze" -em "block" -o energy_all_BST.csv
# Print first and last 3 line of output file
echo "==========================================================="
echo "Deformation free energy of bound DNA as a function of time"
echo "==========================================================="
head -4 energy_all_BST.csv
printf ".\n.\n.\n"
tail -3 energy_all_BST.csv
Some of the plots from above data can be found here
Deformation free energy can be calculated as the following terms:
full : Use entire elastic matrix -- all motions with their couplingdiag : Use diagonal of elastic matrix -- all motions but no couplingb1 : Only bending-1 motionb2 : Only bending-2 motionstretch : Only stretching motiontwist : Only Twisting motionsst_coupling : Only stretch-twist coupling motionbs_coupling : Only Bending and stretching couplingbt_coupling : Only Bending and Twisting couplingbb_coupling : Only bending-1 and bending-2 couplingbend : Both bending motions with their couplingst : Stretching and twisting motions with their couplingbs : Bending (b1, b2) and stretching motions with their couplingbt : Bending (b1, b2) and twisting motions with their couplingWhen all is used, all above terms were calculated.
In [6]:
%%bash
dnaMD globalEnergy -ir free_dna.h5 -ip bound_dna.h5 -tbp 27 -bs 4 -be 24 -estype ST -et "all" -gt "gmx analyze" -em "block" -o energy_all_ST.csv
# Print first and last 3 line of output file
echo "==========================================================="
echo "Deformation free energy of bound DNA as a function of time"
echo "==========================================================="
head -4 energy_all_ST.csv
printf ".\n.\n.\n"
tail -3 energy_all_ST.csv
Same as above but only for stretching and twisting motions.
Local elastic properties can be caluclated using either local base-step parameters or local helical base-step parameters.
In case of base-step parameters: Shift ($Dx$), Slide ($Dy$), Rise ($Dz$), Tilt ($\tau$), Roll ($\rho$) and Twist ($\omega$), following elastic matrix is calculated.
$$ \mathbf{K}_{base-step} = \begin{bmatrix} K_{Dx} & K_{Dx,Dy} & K_{Dx,Dz} & K_{Dx,\tau} & K_{Dx,\rho} & K_{Dx,\omega} \\ K_{Dx,Dy} & K_{Dy} & K_{Dy,Dz} & K_{Dy,\tau} & K_{Dy,\rho} & K_{Dy,\omega} \\ K_{Dx,Dz} & K_{Dy,Dz} & K_{Dz} & K_{Dz,\tau} & K_{Dz,\rho} & K_{Dz,\omega} \\ K_{Dx,\tau} & K_{Dy,\tau} & K_{Dz,\tau} & K_{\tau} & K_{\tau, \rho} & K_{\tau,\omega} \\ K_{Dx,\rho} & K_{Dy,\rho} & K_{Dz,\rho} & K_{\tau, \rho} & K_{\rho} & K_{\rho,\omega} \\ K_{Dx,\omega} & K_{Dy,\omega} & K_{Dz,\omega} & K_{\tau, \omega} & K_{\rho, \omega} & K_{\omega} \\ \end{bmatrix} $$In case of helical-base-step parameters: x-displacement ($dx$), y-displacement ($dy$), h-rise ($h$), inclination ($\eta$), tip ($\theta$) and twist ($\Omega$), following elastic matrix is calculated.
$$ \mathbf{K}_{helical-base-step} = \begin{bmatrix} K_{dx} & K_{dx,dy} & K_{dx,h} & K_{dx,\eta} & K_{dx,\theta} & K_{dx,\Omega} \\ K_{dx,dy} & K_{dy} & K_{dy,h} & K_{dy,\eta} & K_{dy,\theta} & K_{dy,\Omega} \\ K_{dx,h} & K_{dy,h} & K_{h} & K_{h,\eta} & K_{h,\theta} & K_{h,\Omega} \\ K_{dx,\eta} & K_{dy,\eta} & K_{h,\eta} & K_{\eta} & K_{\eta, \theta} & K_{\eta,\Omega} \\ K_{dx,\theta} & K_{dy,\theta} & K_{h,\theta} & K_{\eta, \theta} & K_{\theta} & K_{\theta,\Omega} \\ K_{dx,\Omega} & K_{dy,\Omega} & K_{h,\Omega} & K_{\eta, \Omega} & K_{\theta, \Omega} & K_{\Omega} \\ \end{bmatrix} $$
In [7]:
%%bash
dnaMD localElasticity -i free_dna.h5 -tbp 27 -bs 4 -be 7 -o local_elasticity_4-7bps.csv
Same command can be used to calculate elasticity as a function of time with option
-ot/--output-time and save it in csv format file. This result can be used to check
their convergence.
-fgap/--frame-gap is an essential option.NOTE: Elastic properties cannot be calculated using a single frame because fluctuations are required. Therefore, here time means trajectory between zero time to given time.
In [8]:
%%bash
dnaMD localElasticity -i free_dna.h5 -tbp 27 -bs 4 -be 7 -fgap 200 -ot local_elasticity_time_4-7bps.csv
# Print first 3 and last 3 line of first seven columns from output file
echo "====================================="
echo "Elastic matrix as a function of time"
echo "====================================="
head -4 local_elasticity_time_4-8bps.csv | awk '{print $1, $2, $3, $4, $5, $6 ,$7, "... ... ..."}'
printf ".\n.\n.\n"
tail -3 local_elasticity_time_4-8bps.csv | awk '{print $1, $2, $3, $4, $5, $6 ,$7, "... ... ..."}'
Above method gives local elasticities of a small local segment of the DNA. However, we mostly interested in large segment of the DNA. This large segment can be further divided into smaller local segments. For these smaller segments local elasticities can be calculated. Here these segments overlapped with each other.
Same command can be used to calculate elasticity as a function of time with option
-os/--output-segments and save it in csv format file.
-fgap/--frame-gap is an essential option.
In [9]:
%%bash
dnaMD localElasticity -i free_dna.h5 -tbp 27 -bs 4 -be 20 -span 4 -fgap 200 -gt "gmx analyze" -em "acf" -os local_elasticity_segments.csv
# Print first 3 and last 3 line of first seven columns from output file
echo "========================================"
echo "Elasticity as a function of DNA segments"
echo "========================================"
head -4 local_elasticity_segments.csv | awk '{print $1, $2, $3, $4, $5, $6 ,$7, "... ... ..."}'
printf ".\n.\n.\n"
tail -3 local_elasticity_segments.csv | awk '{print $1, $2, $3, $4, $5, $6 ,$7, "... ... ..."}'
At first, elastic matrix from reference DNA (most often free or unbound DNA) is calculated and subsequently this matrix is used to calculate deformation free energy of probe DNA (most often bound DNA).
$$G = \frac{1}{2}\mathbf{xKx^T}$$When helical='False'
When helical='True'
In [10]:
%%bash
dnaMD localEnergy -ir free_dna.h5 -ip bound_dna.h5 -tbp 27 -bs 10 -be 13 -et all -gt "gmx analyze" -em "block" -o local_energy_time_4-7bps.csv
# Print first 3 and last 3 line from output file
echo "==============================================="
echo "Local deformation energy as a function of time"
echo "==============================================="
head -4 local_energy_time_4-7bps.csv
printf ".\n.\n.\n"
tail -3 local_energy_time_4-7bps.csv
Some of the plots from above data can be found here
Same as the above but energy is calculated using helical base-step parameters
In [11]:
%%bash
dnaMD localEnergy -ir free_dna.h5 -ip bound_dna.h5 -tbp 27 -bs 10 -be 13 -et all -helical -gt "gmx analyze" -em "block" -o local_helical_energy_time_4-7bps.csv
# Print first 3 and last 3 line from output file
echo "======================================================"
echo "Local helical deformation energy as a function of time"
echo "======================================================"
head -4 local_energy_time_4-7bps.csv
printf ".\n.\n.\n"
tail -3 local_energy_time_4-7bps.csv
Above method gives energy of a small local segment of the DNA. However, we mostly interested in large segment of the DNA. This large segment can be further divided into smaller local segments. For these smaller segments local deformation energy can be calculated. Here these segments overlapped with each other.
In [12]:
%%bash
dnaMD localEnergy -ir free_dna.h5 -ip bound_dna.h5 -tbp 27 -bs 4 -be 20 -span 4 -et all -gt "gmx analyze" -em "block" -os local_energy_segments.csv
# Print first 3 and last 3 line from output file
echo "=================================================="
echo "Local deformation energy as a function of segments"
echo "=================================================="
head -4 local_energy_segments.csv
printf ".\n.\n.\n"
tail -3 local_energy_segments.csv
Some of the plots from above data can be found here
Same as the above but energy is calculated using helical base-step parameters
In [13]:
%%bash
dnaMD localEnergy -ir free_dna.h5 -ip bound_dna.h5 -tbp 27 -bs 4 -be 20 -span 4 -et all -helical -gt "gmx analyze" -em "block" -os local_helical_energy_segments.csv
# Print first 3 and last 3 line from output file
echo "=========================================================="
echo "Local helical deformation energy as a function of segments"
echo "=========================================================="
head -4 local_helical_energy_segments.csv
printf ".\n.\n.\n"
tail -3 local_helical_energy_segments.csv