lorenzo-rovigatti / oxDNA

A new version of the code to simulate the oxDNA/oxRNA models, now equipped with Python bindings
https://dna.physics.ox.ac.uk/
GNU General Public License v3.0
38 stars 26 forks source link

Can a backbone connection be approximated with mutual harmonic traps? #87

Closed alec-cook closed 5 months ago

alec-cook commented 6 months ago

Hello,

Is it possible to create a mutual trap which approximates the strength and equilibrium distance of a backbone connection? My goal is to design simulations where the backbone force is only applied in one direction (i.e. one particle exerts a force on its neighbor without feeling a force in return). I decided to nick the backbone between two nucleotides and apply mutual traps to simulate the backbone force, first recapitulating the behavior of a normally bonded strand (with equal traps on either side) before then attempting to change the forces so that they are not equal. However, even in the case with symmetric traps, the simulations always throw the "distance between bonded neighbors exceeds acceptable values" error.

My dna structure contains a 15nt single-strand flanked by two 18bp double-strands, one on each side. I intend to connect the double strands to the single strand with these traps. I am running simple MD simulations at 20C without any additional forces beyond those described.

For the rest distance r0, I have tried 0.7525 (in oxDNA units, the FENE R0 from the oxDNA code) and 0.575 (roughly the average center-of-mass distance between adjacent nucleotides I observe in my simulations). I have varied the stiffness from a default of 1 up to 2.8, a roughly physical spring constant for DNA taken from Epstein et al (2012): https://web.mit.edu/cepstein/Public/CEpsteinDNA.pdf

I have attempted to resolve this by reducing the step size, but the exceeded distance error still persists

Do there exist harmonic parameters that would approximate the behavior of the oxDNA backbone? Do the modifications to the FENE potential (described in section 3.2 here https://arxiv.org/ftp/arxiv/papers/2004/2004.05052.pdf) make that impossible? Is there a simpler way to implement this functionality that I'm not seeing? The forces don't need to exactly match the behavior of a standard backbone connection, but I would at least like to apply a similarly strong force that keeps the unbound, adjacent nucleotides close together without causing the simulation to crash.

Thank you so much for your help, Alec Cook

ErikPoppleton commented 5 months ago

You can use asymmetric harmonic potentials to mimic a backbone, however there are a couple of unphysical details to this approach:

To get a more physical description, you could implement a FENE external force relatively easily, however this would not fix the center of mass problem. If, however, all you care is that the strands are held in proximity, then the harmonic potential with a high stiffness is sufficient.

The "distance between bonded neighbors..." error means that some pair of bonded particles in your starting configuration have a distance outside the values allowed by the backbone FENE potential. To fix this you need to run a short simulation with the modified FENE potential you mention at the end (you do this by setting the max_backbone_force parameter in the input file, which will instruct the model to use the modified potential and bypass the distance check. Documentation here). Make sure you remove this option when running your production run.

alec-cook commented 5 months ago

Thank you for the info! There is nothing wrong with my starting configuration. The simulation will begin to run without issue, but eventually, the distance between bonded neighbors error will crop up. This does not happen if the mutual traps are not added, or if the nucleotides are connected with a backbone instead.

ErikPoppleton commented 5 months ago

Ah, I see the DNA interaction actually throws that error during simulations (I mostly use RNA these days). In that case, it sounds like your asymmetric force is pulling too hard on the nucleotide. Normally I would suggest decreasing the dt, but you already did that. How does the distance between the nucleotides interacting with the harmonic force compare to the ones using normal backbones? How different are the distributions?

lorenzo-rovigatti commented 5 months ago

Do you mind sharing the files so that I can test what's happening?

alec-cook commented 5 months ago

I think my problem may have been a typo in one of the force statements. I am running some additional stimulations now, and the fact that they haven't thrown an error yet is encouraging. I won't get the chance to share my files for a day or two, but I think the bonded neighbor error may be resolved for now.

My focus now is trying to make the harmonic trap as similar as possible to the backbone force, but for the reasons mentioned above, I understand they won't be identical. Is there a stiffness I could use to closely approximate the strength of the backbone at small displacements?

Thank you so much, Alec

Get Outlook for Androidhttps://aka.ms/AAb9ysg


From: Lorenzo Rovigatti @.> Sent: Saturday, January 6, 2024 10:17:48 AM To: lorenzo-rovigatti/oxDNA @.> Cc: Cook, Alexander W @.>; Author @.> Subject: Re: [lorenzo-rovigatti/oxDNA] Can a backbone connection be approximated with mutual harmonic traps? (Issue #87)

Do you mind sharing the files so that I can test what's happening?

— Reply to this email directly, view it on GitHubhttps://github.com/lorenzo-rovigatti/oxDNA/issues/87#issuecomment-1879726136, or unsubscribehttps://github.com/notifications/unsubscribe-auth/AQYPWDQJJNK7TARYWPWT64DYNFTJZAVCNFSM6AAAAABBNNNAGOVHI2DSMVQWIX3LMV43OSLTON2WKQ3PNVWWK3TUHMYTQNZZG4ZDMMJTGY. You are receiving this because you authored the thread.Message ID: @.***>

lorenzo-rovigatti commented 5 months ago

As Erik said, the problem is that external forces are applied to the centres of mass, whereas the FENE interaction acts between the backbone sites, and therefore it's hard to find a good approximant of the latter using the former. However, you could do the following: simulate a single strand made of a pair of T at ~high temperature (50-60 °C), so that stacking doesn't play a big role and use the distance observable to sample the distribution between the centre of mass distance between the two nucleotides, $P(R)$. Now, using a potential of mean-force approach, your effective centre-to-centre interaction is just $U(R) \approx -k_B T \log(P(R))$. Fit a harmonic function to what you find and you should be good to go.

alec-cook commented 5 months ago

Thank you! You may mark this thread as closed; my simulations are working, and I understand the model better than I did before. Thank you for your help!