VMD-L Mailing List

From: Michael Robinson (michael.robinson1_at_monash.edu)
Date: Sun Sep 25 2022 - 23:58:49 CDT

Hi Amanda,

Just quickly jumping in here, as I suspect I can help with these questions:
- If you have three distinct, unconnected chains and only wish to perform
analysis on one of them, but the chain naming has them as the same chain,
you can use some of the more advanced selection language. Setting a
selection of something like "same fragment as [a unique identifier that
only occurs in chain C, perhaps a specific resid + resname combination]"
will select the continuous set of connected atoms from that unique
identifier. You can even use this to fix the chain naming, if you want to
do it that way - create that selection in the TCL console and then set the
chain name back to chain C and write updated .psf and .pdb files.
- When it comes to the selection for GLY311 also showing water and ion
atoms, it might be worth double checking how you're selecting the residue.
Residue ID numbering typically restarts for the water and ion chains when
creating your .pdb and .psf files, so if you are selecting for "resid 311",
you'll see the 311th residue of the water and ion chains (though that would
require a lot of ions - not sure how large your system is). One way to
account for this is to specify you only want the protein - "resid 311 and
protein", "resid 311 and not (water or ions)", or "resid 311 and resname
GLY" should all hopefully fix it, if this is the problem. If the ions and
water molecules continue to be visible with these selections, something
much more strange is happening.

Regards,
Michael Robinson

On Mon, 26 Sept 2022 at 14:12, Amanda Jenkins <
amanda.jenkins1379_at_outlook.com> wrote:

> Hi Josh,
>
> Thank you so much! This did work and now the dihedral angle calculations
> appear to be able to calculate values for phi and psi.
>
> My follow-up question is, when I follow your instructions (after using dcd
> and psf as input files), within VMD, how do I delete multiple chain
> selections, save for the chain I want to keep to analyze? I need to do my
> analysis on chain C only but since the psf and dcd files present all three
> chains as chain P, as seen in the Sequence Viewer, all three chain angle
> results are lumped together in each calculated frame that Axel's script
> provides. This doesn't allow me to conveniently discriminate the results
> for chain C from the other chains.
>
> As an aside, I just wanted to mention that after viewing the protein after
> following your instructions, I noticed that there appears to be one
> renegade water molecule and a couple of lone atoms (could be sodiums and
> chlorides) that respond to viewing GLY 311 in isolation, in addition to the
> actual GLY 311, and I don't understand why since sodium, chloride and water
> are not glycines. When I made a PDB of only the protein after following
> your instructions, the resulting PDB file had the same characteristic of
> treating chain C as two fragments and the renegade atoms were no longer
> present in that PDB file. I wasn't sure if these glycine "impostors" are
> crucial to creating a PDB that doesn't split chain C into fragments.
>
> Thank you again Josh for your help, and if there's no convenient way to
> delete the other chains within VMD and run the dihedral angles script
> immediately afterards (that way I only get chain C dihedrals by keeping the
> glycine imposters), then I'll try processing the consolidated dihedrals
> file I mentioned above with scripting.
>
> Amanda J.
> ------------------------------
> *From:* Vermaas, Josh <vermaasj_at_msu.edu>
> *Sent:* Sunday, September 25, 2022 1:29 PM
> *To:* Amanda Jenkins <amanda.jenkins1379_at_outlook.com>; vmd-l_at_ks.uiuc.edu <
> vmd-l_at_ks.uiuc.edu>
> *Subject:* Re: vmd-l: Ramachandran Dihedral Angles Suspiciously Constant
>
>
> Hi Amanda,
>
>
>
> Yeah, that does seem to be the root of the problem. When I loaded the pdb
> you sent along, VMD guesses the bonds present in the system based on a
> distance-based heuristic. As you noted, the long bond between the backbone
> nitrogen and the alpha carbon means that the distance based heuristic did
> not place a bond there. As a result, resid 311 is split between two
> “fragments”. A fragment in VMD is a connected set of atoms that get treated
> as one entity. Since there is no bond between the nitrogen and the carbon,
> the phi/psi calculations act as though Gly311 is not actually a residue
> within a protein, since the nitrogen is in a different fragment from the
> rest of the protein backbone. Since you **did** simulate this in NAMD,
> the solution is simple. Using the same selection you used to create your
> cut-down pdb file, **also** make a cut down psf, and load the psf (which
> has bond information) as well as the pdb file together when running your
> analysis. Alternatively, just analyze the psf/dcd together. Something like
> this is what I’d do:
>
>
>
> mol new step3_charmm2namd.psf
>
> mol addfile snapshot-1.dcd waitfor all
>
> #Once you load these two in, there shouldn’t be a break at residue 311.
>
> #Then you’d use your script as usual.
>
> *set mol [molinfo top]*
>
> *set fp [open "BOB.txt" w]*
>
> *set sel [atomselect $mol "protein and name CA"]*
>
> *set n [molinfo $mol get numframes]*
>
> *for {set i 0} {$i < $n} {incr i} {*
>
> *$sel frame $i*
>
> *puts $fp "\# frame: $i" *
>
> *set a [$sel num] *
>
> *for {set j 0} {$j < $a} {incr j} { *
>
> *puts $fp "[expr $j + 1] [lindex [$sel get {resname phi psi}] $j]" *
>
> *} *
>
> *}*
>
> *$sel delete*
>
> *close $fp*
>
>
>
> -Josh
>
>
>
> *From: *Amanda Jenkins <amanda.jenkins1379_at_outlook.com>
> *Date: *Saturday, September 24, 2022 at 10:22 PM
> *To: *"Vermaas, Josh" <vermaasj_at_msu.edu>, "vmd-l_at_ks.uiuc.edu" <
> vmd-l_at_ks.uiuc.edu>
> *Subject: *Re: vmd-l: Ramachandran Dihedral Angles Suspiciously Constant
>
>
>
> Hi Josh,
>
>
>
> I wanted to let you know that after zooming in on the problematic amino
> acid in isolation in VMD, it looks like that amino acid, for some strange
> reason, has an abnormal distance between its amino nitrogen and alpha
> carbon and VMD is not drawing a bond between these two atoms. I don't know
> if this might be the stumbling block for VMD to calculate the dihedral
> angles. Any thoughts on why this is happening? Thanks.
>
>
>
> Amanda J.
> ------------------------------
>
> *From:* Vermaas, Josh <vermaasj_at_msu.edu>
> *Sent:* Saturday, September 24, 2022 12:28 PM
> *To:* Amanda Jenkins <amanda.jenkins1379_at_outlook.com>; vmd-l_at_ks.uiuc.edu <
> vmd-l_at_ks.uiuc.edu>
> *Subject:* Re: vmd-l: Ramachandran Dihedral Angles Suspiciously Constant
>
>
>
> Huh. Could you send the psf/pdb off list? That is actually what dictates
> what VMD recognizes as protein, and would be more helpful in figuring out
> why your protein isn’t being picked up as protein.
>
>
>
> -Josh
>
>
>
> *From: *Amanda Jenkins <amanda.jenkins1379_at_outlook.com>
> *Date: *Friday, September 23, 2022 at 10:19 PM
> *To: *"Vermaas, Josh" <vermaasj_at_msu.edu>, "vmd-l_at_ks.uiuc.edu" <
> vmd-l_at_ks.uiuc.edu>
> *Subject: *Re: vmd-l: Ramachandran Dihedral Angles Suspiciously Constant
>
>
>
> Hi Josh,
>
>
>
> Yes, this script was indeed originally authored by Axel Kohlmeyer.
>
>
>
> The problematic amino acid in question is not a terminal amino acid but it
> is actually located in the middle of the protein which is what surprised me
> with the 0.0 at all snapshots. I do not have any other molecules or
> ligands in the system, just my protein alone.
>
>
> I have attached to this e-mail a representative conf file with the
> settings I used. Hopefully it can help diagnose the issue. Thanks for
> your help.
>
>
>
> Amanda J
> ------------------------------
>
> *From:* Josh Vermaas <vermaasj_at_msu.edu>
> *Sent:* Friday, September 23, 2022 3:36 PM
> *To:* Amanda Jenkins <amanda.jenkins1379_at_outlook.com>; vmd-l_at_ks.uiuc.edu <
> vmd-l_at_ks.uiuc.edu>
> *Subject:* Re: vmd-l: Ramachandran Dihedral Angles Suspiciously Constant
>
>
>
> Hi Amanda,
>
> This is a fun one! You did indeed do exactly what I would have done, and
> use the "$sel frame $i" to set the correct frame. The update is strictly
> speaking not needed, unless the selection changes as a function of time.
> This looks like what Axel suggested back in 2005 (
> https://www.ks.uiuc.edu/Research/vmd/mailing_list/vmd-l/4885.html
> <https://urldefense.com/v3/__https:/www.ks.uiuc.edu/Research/vmd/mailing_list/vmd-l/4885.html__;!!HXCxUKc!1qS9Pd488sYskLzSfGlq8feeJPzLSN_KMuLzPMOR773KOUpWvsruYpkBn1z31HRm4T6cLaAXTYc3V6oM0Bz7I4X_WHvP$>),
> so this looks reasonable.
>
> But I think I now understand the problem. Is the amino acid that is
> reporting consistently 0 either at the beginning or end of the chain? The
> first amino acid doesn't have a phi angle, and the last amino acid doesn't
> have a psi angle. Likewise, if you had something that wasn't a protein in
> your system, but *did* have a CA atom, that might also not have a phi or
> psi angle. In VMD, not having a phi or psi angle manifests as *always*
> reporting 0. What I did to test this was to take your script and run it
> against the ubiquitin in a box simulation from the NAMD tutorial. That made
> it clear that your script was actually fine, and that you might be running
> into something weird with how your system is constructed.
>
> -Josh
>
> On 9/23/22 2:30 PM, Amanda Jenkins wrote:
>
> Hello Josh,
>
>
>
> I am calculating the dihedral angles (phi and psi) using a Tcl TkCon
> script which is below:
>
>
>
> *set mol [molinfo top]*
>
> *set fp [open "BOB.txt" w]*
>
> *set sel [atomselect $mol "all and name CA"]*
>
> *set n [molinfo $mol get numframes]*
>
> *for {set i 0} {$i < $n} {incr i} {*
>
> *$sel frame $i*
>
> *$sel update*
>
> *puts $fp "\# frame: $i" *
>
> *set a [$sel num] *
>
> *for {set j 0} {$j < $a} {incr j} { *
>
> *puts $fp "[expr $j + 1] [lindex [$sel get {resname phi psi}] $j]" *
>
> *} *
>
> *}*
>
> *$sel delete*
>
> *close $fp*
>
>
>
> Is the line *$sel update* (bolded) insufficient to update the
> trajectory? Thanks.
>
>
>
> Amanda J.
> ------------------------------
>
> *From:* Josh Vermaas <vermaasj_at_msu.edu> <vermaasj_at_msu.edu>
> *Sent:* Friday, September 23, 2022 7:54 AM
> *To:* Amanda Jenkins <amanda.jenkins1379_at_outlook.com>
> <amanda.jenkins1379_at_outlook.com>; vmd-l_at_ks.uiuc.edu <vmd-l_at_ks.uiuc.edu>
> <vmd-l_at_ks.uiuc.edu>
> *Subject:* Re: vmd-l: Ramachandran Dihedral Angles Suspiciously Constant
>
>
>
> Hi Amanda,
>
> How are you calculating phi/psi angles? I'm betting that your script isn't
> updating the selection over your trajectory. A way to confirm this visually
> is to label the specific dihedral of interest (on the VMD main window
> Mouse->Label->Dihedral), and animate through your trajectory. I suspect
> you'll see non-zero values. It is basically impossible for a trajectory to
> give you exactly a 0.0 dihedral for all frames, but it is relatively easy
> (and I've done it) to forget to update your selection within Tcl so that
> VMD is calculating on the right set of coordinates.
>
> -Josh
>
> On 9/23/22 6:14 AM, Amanda Jenkins wrote:
>
> Dear community,
>
>
>
> I have results/PDB trajectories from several MD simulations/repeats of a
> protein composed of several chains. For one of the chains, one of the
> amino acids, always has values of 0.0 for both phi and psi angles
> throughout the trajectory simulation. This is observed throughout each
> different set/repeat of simulation results and at the same amino acid in
> the same chain. I have RMSD aligned each individual trajectory/simulation
> result set, so my initial theory was that maybe this amino acid in question
> is serving as the "anchor" point upon which all the other snapshots are
> superimposed for a given simulation set, but it still seems a bit odd to me
> that both phi and psi angles would be perfectly 0.0 at every single
> trajectory snapshot. I have other amino acids that have small/fractional
> phi and psi angles so I am suspicious that this perfect set of 0.0 phi and
> psi angles might be an artifact of some sort and not an actual physical
> phenomenon.
>
>
>
> Would someone please tell me what could be causing this situation? Is it
> an artifact?
>
>
>
> For additional information, my total simulation time is 200 nanoseconds
> with a snapshot being printed every 1 nanosecond, and this is done for all
> of my distinct simulation result sets/repeats.
>
>
>
> Thanks for any insight anyone can provide.
>
>
>
> Amanda J.
>
>
>
> --
>
> Josh Vermaas
>
>
>
> vermaasj_at_msu.edu
>
> Assistant Professor, Plant Research Laboratory and Biochemistry and Molecular Biology
>
> Michigan State University
>
> vermaaslab.github.io
>
>
>
> --
>
> Josh Vermaas
>
>
>
> vermaasj_at_msu.edu
>
> Assistant Professor, Plant Research Laboratory and Biochemistry and Molecular Biology
>
> Michigan State University
>
> vermaaslab.github.io
>
>