How to Dock a Ligand to DNA or RNA (rDock and Vina)

Docking a ligand to DNA or RNA follows the same three steps as protein docking, but each one changes. Keep the structural waters and divalent ions the pocket depends on, define the site from a reference ligand instead of a blind cavity search, and score with rDock, which was built for nucleic acids.

Most docking tutorials assume the target is a protein. Then a project hands you a DNA duplex and a set of phytochemicals, or an RNA riboswitch and a fragment library, and the familiar recipe quietly stops working. The cavity detector returns nothing useful. The scores come out oddly uniform. Nobody explains which part broke. This guide is the missing part, written by the StemSkills Lab team from 10+ years in structural bioinformatics, drug design and multiscale molecular modeling.

What actually changes when the receptor is DNA or RNA instead of a protein?

Four things change, and every downstream problem traces back to one of them.

• The site is open, not enclosed. A protein pocket is usually a buried cavity with a clear mouth. A DNA minor groove is a long continuous channel, and an intercalation site does not exist until the base pairs separate. Cavity detectors that look for enclosure find either nothing or the entire helix.
• The surface is charged and polar. Every phosphate carries a negative charge, so electrostatics and desolvation dominate where a protein pocket would offer a hydrophobic wall.
• Waters and ions are part of the structure. Ordered waters in the groove and divalent cations bridging phosphates hold the fold in place. Deleting every non-standard record removes real chemistry.
• The scoring functions were fitted mostly on proteins. Public protein-ligand complexes outnumber nucleic-acid-ligand complexes by orders of magnitude, so most scoring functions have simply seen far fewer of your kind of target.

The sequence below assumes you already know what a grid box and a pose are. If not, start from our molecular docking pillar guide first.

How do you prepare a clean DNA or RNA receptor?

Work from a real structure and keep a record of every atom you remove. Three examples you can download today and practise on:

• PDB 121D, the A-tract dodecamer d(CGCAAATTTGCG) with the minor groove binder netropsin, solved at 2.2 Å. A clean groove-binding reference.
• PDB 1D10, adriamycin and daunomycin bound to DNA at 1.5 Å. An intercalation reference with the base pairs already prised apart.
• PDB 2GDI, the thiamine pyrophosphate riboswitch with TPP bound, solved at 2.05 Å. This entry also contains magnesium, potassium and sodium ions, which is exactly the situation that catches students out.

Then work through this order:

1. Strip the protein, keep the nucleic acid. Many nucleic acid ligand complexes are deposited inside a protein assembly. Pull out only the chains you want before anything else, using the same chain-selection approach as our guide on extracting a single chain from a PDB file. A leftover protein chain will silently become part of your receptor surface.
2. Decide about waters deliberately, not by habit. The default advice to delete all waters comes from protein docking. In a groove complex, one or two bridging waters can be the difference between a sensible pose and a floating ligand. Our guide on removing water molecules from a PDB file shows how to keep a named subset instead of all or nothing.
3. Keep the ions that organise the pocket. In 2GDI the magnesium ions coordinate the pyrophosphate group of the bound ligand. Remove them and you have changed the problem. The handling rules are the same ones in our guide on metal ions in docking and MD.
4. Fix the protonation of the ligand, not the backbone. Phosphates are deprotonated at physiological pH, so the receptor side is easy. The ligand is where mistakes happen, especially amines and carboxylates on natural products.
5. Convert to the format your docking program demands. rDock requires a mol2 receptor. AutoDock requires PDBQT. Neither accepts a raw PDB. Use the conversions in our guide to converting molecular file formats with Open Babel and check the output by eye rather than trusting the converter.

On the AutoDock side, one piece of good news: receptor preparation now understands nucleic acids properly. Meeko’s default residue template set is built from Amber24 OFF libraries and ships RNA.lib, DNA.OL15.lib, parmBSC1.lib and all_modrna08.lib, stored in Meeko/meeko/data/residue_chem_templates.json. Modified nucleotides are covered by a separate curated set based on the Nucleic Acid Knowledgebase, loaded with --add_templates. The details are in the Meeko template documentation.

Where do you put the box: groove binding, intercalation, or a defined RNA pocket?

This is where nucleic acid docking most often fails silently, because a box that is technically valid can still be scientifically meaningless.

• A defined pocket, such as a riboswitch aptamer. Easiest case, and it behaves like a protein. If a crystal ligand exists, build the box around it. With Meeko you do not have to read coordinates by hand: --box_enveloping takes the reference ligand file and --padding adds the margin.
• Groove binding. There is no single pocket. Restrict the box to the sequence stretch your hypothesis is about, usually four to six base pairs of an AT-rich tract for a minor groove binder, and say so in your methods. A box around the whole duplex is a blind search, not a hypothesis.
• Intercalation. Rigid-receptor docking cannot open a base-pair step, so you cannot dock an intercalator into a closed B-DNA model and expect a correct pose. Either start from a structure already in the intercalated state, such as 1D10, or treat the result as a groove-surface pose and say so.

The mechanics of centre and size are unchanged from the protein case, so our guide on setting the AutoDock Vina grid box still applies. What changes is the justification you owe your reader.

Which docking program should you use for a nucleic acid target?

The answer depends on what you are docking, not on which program you already know. A small molecule into RNA and a protein onto DNA are different computational problems handled by different software.

rDock is the one most students have never heard of, and the one built for this job. It is published as rDock: A Fast, Versatile and Open Source Program for Docking Ligands to Proteins and Nucleic Acids by Ruiz-Carmona and colleagues in PLOS Computational Biology, 2014. It is not a new experiment: the code descends from RiboDock, developed at RiboTargets, later Vernalis (R&D) Ltd, between 1998 and 2006, with RNA targets as a first-class use case rather than an afterthought.

How do you run rDock on an RNA or DNA target?

rDock is a three-stage command-line workflow. Prepare a system definition parameter file that points at your mol2 receptor, map the cavity, then dock.

Stage 1, cavity mapping. For an open nucleic acid surface, use the reference ligand mapper rather than a blind two-sphere search. In the SECTION MAPPER of your parameter file, set the site mapper to RbtLigandSiteMapper and give it your reference ligand with REF_MOL, then tune RADIUS, SMALL_SPHERE, MIN_VOLUME, MAX_CAVITIES, VOL_INCR and GRIDSTEP. Then run:

rbcavity -was -d -r system.prm

The -d flag writes a grid file you can load next to the receptor in PyMOL or ChimeraX. Look at it before you dock. If the mapped cavity covers the whole helix, your RADIUS or MIN_VOLUME is wrong, and no amount of sampling rescues the run.

Stage 2, docking. The minimum rDock invocation names the input and output ligand SD files, the system definition file and the docking protocol:

rbdock -i ligands.sd -o docked -r system.prm -p dock.prm -n 50

Here -n 50 requests 50 docking runs per ligand. The rDock docking strategies documentation covers the parts that matter specifically for nucleic acids. Explicit structural waters load from a separate PDB file through a SECTION SOLVENT block, and each water can be set independently to FIXED, TETHERED or FREE for translation and rotation, which gives nine flexibility modes in total. That is the feature to reach for when a groove water matters but you do not want to nail it in place. Pharmacophore restraints go in a SECTION PHARMA block, and tethered scaffold docking needs the ligand SD file annotated first with rbtether.

Stage 3, the AutoDock cross-check. If the site is a defined pocket, run the same system through Vina and compare:

mk_prepare_receptor.py --read_pdb rna_clean.pdb -o rna_receptor -p -v \
  --box_enveloping reference_ligand.sdf --padding 4
mk_prepare_ligand.py -i ligand.sdf -o ligand.pdbqt
vina --receptor rna_receptor.pdbqt --ligand ligand.pdbqt \
  --config rna_receptor.box.txt --exhaustiveness=32 --out docked.pdbqt

Vina’s default exhaustiveness is 8. The official AutoDock Vina basic docking tutorial uses 32, and on an open, featureless surface the higher setting is worth the time. Two programs agreeing on a binding mode is a far stronger claim than one program producing a number.

Why do the scores look worse than on a protein target?

Because they are worse, and this is measured rather than anecdotal. The rDock paper validated the same program on two sets. On 85 diverse protein-ligand complexes from the CCDC/Astex set, the success rate was about 80% at an RMSD threshold of 2 Å. On 56 RNA-ligand complexes, the success rate was 54 plus or minus 3% at a looser threshold of 2.5 Å. The authors are explicit about why RNA is harder, describing RNA sites as having “less closed cavities, less hydrophobic, featureless” character compared with protein pockets, with larger and more flexible ligands on top of that.

Two consequences follow. Roughly a coin flip on a program designed for nucleic acids is the honest baseline, so treat a single pose from a single program as a hypothesis, not a result. And the threshold moves: the RNA benchmark used 2.5 Å, not the 2.0 Å convention from protein redocking, so validate against the threshold your target class actually uses.

What if you are docking a protein to DNA or RNA instead?

Then none of the above applies, because you have a macromolecular docking problem. Two free servers handle it well.

HDOCK (Yan, Zhang, Zhou, Li and Huang, Nucleic Acids Research 45(W1):W365-W373, 2017) combines template-based modelling with free FFT-based docking and an iterative knowledge-based scoring function, and handles protein to protein, protein to DNA and protein to RNA. It does not dock small molecules, so it is the wrong tool for a phytochemical into a groove. NPDock (Tuszynska and colleagues, Nucleic Acids Research, 2015) is built specifically for protein-nucleic acid complexes and runs a pipeline of global docking with GRAMM, scoring with a statistical potential, clustering of the best-scored models, and local refinement of the most promising ones.

If macromolecular docking is where your project is heading, our guide on protein-protein docking with ClusPro, HDOCK and HADDOCK covers the server workflow and how to compare their outputs, and the same reasoning about flexibility carries over from docking a peptide to a protein.

How do you troubleshoot the errors this workflow throws?

• Receptor preparation reports an unknown or unmatched residue. Usually a modified nucleotide, a pseudouridine or a 2′-O-methyl sugar. Load the curated modified-nucleotide templates with --add_templates, or map the residue explicitly with -n / --set_template.
• rDock refuses your receptor. It wants mol2, not PDB. Convert it, then confirm the atom and bond counts survived the conversion.
• rbcavity returns no cavity, or one cavity the size of the helix. The two-sphere method needs enclosure that a groove does not have. Switch to RbtLigandSiteMapper with a REF_MOL, and raise MIN_VOLUME to discard the trivial surface dimples.
• Every pose scores within a few tenths of a kcal/mol. Expected on a featureless surface. Cluster the poses by geometry and look at which contacts recur, rather than ranking numbers inside the noise.
• The ligand drifts along the groove between runs. Either the box is too long for the hypothesis, or a bridging water is missing. Narrow the box to the base-pair stretch you care about, and reintroduce the specific waters through SECTION SOLVENT as TETHERED.
• Alternate conformations confuse the receptor. High-resolution nucleic acid structures often carry altlocs. Pick one with --default_altloc or --wanted_altloc instead of letting the tool guess.
• The intercalator sits on the outside of the duplex. Not a bug. A rigid closed duplex has nowhere to intercalate. Start from an intercalated structure instead.

How do you report a nucleic acid docking result so a reviewer believes it?

The methods paragraph carries the burden, because the reader cannot see what you deleted. State the PDB ID and resolution, which chains you kept, which waters and ions you kept and why, the program and version, how the cavity or box was defined, and the number of runs per ligand. If a crystal ligand exists, redock it and report the RMSD against the 2.5 Å threshold. Then describe the binding mode structurally, naming the groove, the base-pair step or the residues contacted, rather than only quoting a score.

For the full reporting template, our guide on writing the methods section of a docking study gives the paragraph structure, and if your ligands are natural products the selection logic in molecular docking of phytochemicals applies unchanged. To see where this sits in a full skill set, follow the computational biology skills roadmap.

Frequently asked questions

Can AutoDock Vina dock a ligand to DNA at all?

Yes. Receptor preparation with Meeko ships DNA and RNA residue templates, so the PDBQT file builds correctly. The limitation is the scoring function, which was developed for protein-ligand complexes, not the file handling. Use Vina on a well-defined pocket and cross-check with rDock.

Should I keep magnesium ions in an RNA receptor?

Keep the ones that are structurally involved, such as ions coordinating the ligand or bridging phosphates in the pocket. Remove distant surface ions. Record the decision in your methods either way, because it changes the result.

Can I use HDOCK to dock a small molecule into an RNA pocket?

No. HDOCK docks macromolecules: protein to protein, protein to DNA and protein to RNA. For a small molecule into RNA, use rDock or Vina.

What RMSD threshold counts as a successful redocking on a nucleic acid?

The rDock RNA-ligand benchmark used 2.5 Å, looser than the 2.0 Å convention for protein redocking, because nucleic acid sites are more open and the ligands are more flexible. Report which threshold you used.

Do I need molecular dynamics after docking to a nucleic acid?

If the claim matters, yes. Nucleic acids are more deformable than proteins, so a rigid-receptor pose is a weaker piece of evidence here than it is on a protein. A short MD run that keeps the ligand in place is the cheapest way to show the pose survives motion.