Question: Use for protein-dna binda energy

[EvanKomp]

Hello!

Personal context: I have run some basic MD simulations, and familiar with high level concepts, but not fine tuned considerations on parameterization, etc. I have recently been reading up on AFE methods, particularly using this resource: doi:10.33011/livecoms.2.1.18378. I have a high level understanding of the steps in your example scripts by referencing the API

Use case: I need a relatively high throughput way to estimate the change in BFE of a variant protein to a DNA compared to a WT protein. I have a high resolution structure of the WT protein-dna complex, and I have some code that aligns the variant protein structure into the complex replacing the wt structure. I would like to do this algorithmically from within python for an upcoming publication.

Questions:

• Is the library capable of dealing with the protein-DNA use case?
• If so, I assume I should minimize the complex before running AFE sims?
• How much difference and which parts of the example script would require a lot of change? I am particularly worried about the atom matching step eg.

  #alignment restraints between lig1 and lig2
  lig1_ref_atoms  = [ refatoms_lig1[i]+lig1_atoms[0] for i in range(3)]
  lig2_ref_atoms  = [ refatoms_lig2[i]+lig2_atoms[0] for i in range(3)]
  atm_utils.addAlignmentForce(liga_ref_particles = lig1_ref_atoms,
                          ligb_ref_particles = lig2_ref_atoms,
                          kfdispl = 2.5 * kilocalorie_per_mole/angstrom**2,
                          ktheta =  10.0 * kilocalorie_per_mole,
                          kpsi =  10.0 * kilocalorie_per_mole,
                          offset = lig2_restr_offset)

Thanks for any input, Evan

[egallicc]

This package is aimed at modeling the binding of small molecules. It cannot reliably estimate the binding free energy of a protein to a nucleic acid nor the change in binding free energy upon mutation, even if the mutation is small. This limitation stems from using a dual-topology approach.

[EvanKomp]

I see. Thank you for the prompt response.

Though this is not the ideal forum, from a field expert's perspective, what are your thoughts on the following:

Digging into it more, I see some recent work that can build the topology considering some number of amino acid mutations (the software pmx). For my use case, I have potentially deletions and insertions, which that software does not support. This seems like an extra formidable challenge, but still possible? Something like:

1. Align the two proteins using a global sequence alignment and a scoring matrix in order to identify "equivalent" residues
2. Treat equivalent residues as pmx does, eg. build the hybrid topology.
3. For residues with no equivalent, have ghost atoms in the topology
4. These residues fade in and out over lambda similar to how the ghost atoms would for a ligand with varying functional groups

This would be extra challenging for the simulation if multiple residues in a row are appearing or disappearing but seems doable?

Thanks for your contributions to the community!

[egallicc]

I am far from an expert, but, yes, my sense is that implementations exist that should do what you describe out of the box (more or less). I describe a NAMD implementation here: Methods in Molecular Biology, Computational Peptide Science: Methods and Protocols, Springer Nature, edited by Thomas Simonson (2021). See Section 3.2. In general, any MD engine that supports alchemical transformations that modify a portion of a molecule should, in principle, be applicable to your problem.