[gmx-users] Percentage of H-bonds

Justin A. Lemkul jalemkul at vt.edu
Wed Oct 27 13:12:26 CEST 2010



Carla Jamous wrote:
> Hi Justin,
> 
> I guess the problem is in the way I concatenate my trajectory.
> So my question may seem stupid but I can't seem to understand this part 
> of the manual:
> 
> c (continue) - The start time is taken from the end
> of the previous file. Use it when your continuation run
> restarts with t=0.
> 
> l (last) - The time in this file will be changed the
> same amount as in the previous. Use it when the time in the
> new run continues from the end of the previous one,
> since this takes possible overlap into account.
> 

Where is this in the manual?

> When I use trjcat to concatenate three different trajectories:
> 
> trjcat -f traj1.xtc traj2.xtc traj3.xtc -o traj_cat.xtc -settime
> 
> Should I use "c" or "l"?
> 

I see no way to incorporate any such option, if it even exists - what version of 
Gromacs are you using?  In any case, the use of -settime makes all of this 
irrelevant.  You're specifically setting the start time of each of the 
trajectories.  Use gmxcheck to make sure you got what you think you did.

What is the reason for -settime?  Are the trajectories non-continuous?

I'm CC'ing this back to the list.  I only offered to take a look at your files 
off-list.  Further questions and feedback need to stay on gmx-users.

-Justin

> Thanks,
> Carla
> 
> On Tue, Oct 26, 2010 at 2:39 PM, Justin A. Lemkul <jalemkul at vt.edu 
> <mailto:jalemkul at vt.edu>> wrote:
> 
> 
> 
>     Carla Jamous wrote:
> 
> 
>         Yes, I'm always looking at the same one, it means the
>         corresponding one in the other trajectory.
> 
> 
>     If you want me to try to figure out what's going on, send me your
>     .xpm and .ndx files (off-list), as well as a description of what I'm
>     looking for and why you think it shouldn't be there.
> 
>     -Justin
> 
> 
>         Carla
> 
> 
>         On Tue, Oct 26, 2010 at 2:17 PM, Justin A. Lemkul
>         <jalemkul at vt.edu <mailto:jalemkul at vt.edu>
>         <mailto:jalemkul at vt.edu <mailto:jalemkul at vt.edu>>> wrote:
> 
> 
> 
>            Carla Jamous wrote:
> 
> 
>                Please I have a problem:
> 
>                I concatenated 3 trajectories and got one long trajectory on
>                which I ran the perl script it gives me lets say 60% for
>         hbond 1.
> 
>                But when I look in each of the 3 trajectories alone, this
>         hbond
>                doesn't even exist in 2 of them.
>                Do you have an idea where the problem might be because I
>         checked
>                if it's well concatenated and it is. And I can't figure
>         out why
>                I have this error.
> 
> 
>            The hydrogen bond indices are regenerated in each trajectory,
>         since
>            not all H-bonds might be present in the same order.  Are you sure
>            you're continually looking at the same one?
> 
>            -Justin
> 
>                Thanks,
>                Carla
> 
> 
>                On Tue, Oct 26, 2010 at 2:08 PM, Justin A. Lemkul
>                <jalemkul at vt.edu <mailto:jalemkul at vt.edu>
>         <mailto:jalemkul at vt.edu <mailto:jalemkul at vt.edu>>
>                <mailto:jalemkul at vt.edu <mailto:jalemkul at vt.edu>
>         <mailto:jalemkul at vt.edu <mailto:jalemkul at vt.edu>>>> wrote:
> 
> 
> 
>                   Carla Jamous wrote:
> 
>                       Hi Justin,
> 
>                       Please can you explain this particular comment in
>         your perl
>                       script: "# There should now be $nres lines left in
>         the file
>                       # The HB map for the last index is written first
>         (top-down in
>                       .xpm file)
>                       #   * Element 0 is index $nres, element 1 is
>         $nres-1, etc."
> 
>                       Is the perl script reading the index file
>         beginning from the
>                       bottom? If so, why is that? Because in the manual, it
>                says that
>                       in the matrix, the ordering is identical to that
>         in the
>                index file.
> 
> 
>                   The order is the same, and the indices are mapped
>         exactly as you
>                   would expect them.  The .xpm file is written, however,
>         from top
>                   down.  So the final index in hbond.ndx is the top line
>         in the
>                   matrix.  As such, even though they are in the same
>         "order,"
>                per se,
>                   the two files (.xpm and .ndx) have to be read in opposite
>                directions.
> 
>                   -Justin
> 
>                       Thank you.
> 
>                       Carla
> 
> 
>                       On Mon, Oct 11, 2010 at 4:31 PM, Justin A. Lemkul
>                       <jalemkul at vt.edu <mailto:jalemkul at vt.edu>
>         <mailto:jalemkul at vt.edu <mailto:jalemkul at vt.edu>>
>                <mailto:jalemkul at vt.edu <mailto:jalemkul at vt.edu>
>         <mailto:jalemkul at vt.edu <mailto:jalemkul at vt.edu>>>
>                       <mailto:jalemkul at vt.edu <mailto:jalemkul at vt.edu>
>         <mailto:jalemkul at vt.edu <mailto:jalemkul at vt.edu>>
>                <mailto:jalemkul at vt.edu <mailto:jalemkul at vt.edu>
>         <mailto:jalemkul at vt.edu <mailto:jalemkul at vt.edu>>>>> wrote:
> 
> 
>                          I wrote a Perl script to do a similar task
>         (appended
>                below).
>                           Perhaps it will be useful to you.  I hope it
>         works; I
>                had to
>                       hack
>                          out some things that were specific to my needs and
>                have done only
>                          limited testing.
> 
>                          -Justin
> 
> 
> 
>                          #!/usr/bin/perl
> 
>                          #
>                          # plot_hbmap.pl <http://plot_hbmap.pl>
>         <http://plot_hbmap.pl>
>                <http://plot_hbmap.pl> <http://plot_hbmap.pl>
> 
>                       - plot the probability of
> 
>                          finding a particular hydrogen bond
>                          # based on several input files:
>                          #   1. coordinate file (for atom naming) - MUST
>         be a .pdb
>                       file with
>                          NO CHAIN IDENTIFIERS
>                          #   2. hbmap.xpm
>                          #   3. hbond.ndx (modified to contain only the atom
>                numbers
>                       in the
>                          [hbonds...] section, nothing else)
>                          #
> 
>                          use strict;
> 
>                          unless(@ARGV) {
>                             die "Usage: perl $0 -s structure.pdb -map
>         hbmap.xpm
>                -index
>                          hbond.ndx\n";
>                          }
> 
>                          # define input hash
>                          my %args = @ARGV;
> 
>                          # input variables
>                          my $map_in;
>                          my $struct;
>                          my $ndx;
> 
>                          if (exists($args{"-map"})) {
>                             $map_in = $args{"-map"};
>                          } else {
>                             die "No -map specified!\n";
>                          }
> 
>                          if (exists($args{"-s"})) {
>                             $struct = $args{"-s"};
>                          } else {
>                             die "No -s specified!\n";
>                          }
> 
>                          if (exists($args{"-index"})) {
>                             $ndx = $args{"-index"};
>                          } else {
>                             die "No -index specified!\n";
>                          }
> 
>                          # open the input
>                          open(MAP, "<$map_in") || die "Cannot open input map
>                file!\n";
>                          my @map = <MAP>;
>                          close(MAP);
> 
>                          open(STRUCT, "<$struct") || die "Cannot open input
>                coordinate
>                       file!\n";
>                          my @coord = <STRUCT>;
>                          close(STRUCT);
> 
>                          open(NDX, "<$ndx") || die "Cannot open input index
>                file!\n";
>                          my @index = <NDX>;
>                          close(NDX);
> 
>                          # determine number of HB indices and frames
>                          my $nres = 0;
>                          my $nframes = 0;
>                          for (my $i=0; $i<scalar(@map); $i++) {
>                             if ($map[$i] =~ /static char/) {
>                                 my $res_line = $map[$i+1];
>                                 my @info = split(" ", $res_line);
> 
>                                 $nframes = $info[0];
>                                 my @nframes = split('', $nframes);
>                                 shift(@nframes);    # get rid of the "
>                                 $nframes = join('', @nframes);
> 
>                                 $nres = $info[1];
>                             }
>                          }
> 
>                          print "Processing the map file...\n";
>                          print "There are $nres HB indices.\n";
>                          print "There are $nframes frames.\n";
> 
>                          # initialize hashes for later output writing
>                          # counter $a holds the HB index from hbond.ndx
>                          my %hbonds;
>                          for (my $a=0; $a<$nres; $a++) {
>                             $hbonds{$a+1} = 0;
>                          }
> 
>                          # donor/acceptor hashes for bookkeeping purposes
>                          my %donors;
>                          for (my $b=1; $b<=$nres; $b++) {
>                             $donors{$b} = 0;
>                          }
> 
>                          my %acceptors;
>                          for (my $c=1; $c<=$nres; $c++) {
>                             $acceptors{$c} = 0;
>                          }
> 
>                          # clean up the output - up to 18 lines of
>         comments, etc.
>                          splice(@map, 0, 18);
> 
>                          # remove any "x-axis" or "y-axis" lines
>                          for (my $n=0; $n<scalar(@map); $n++) {
>                             if (($map[$n] =~ /x-axis/) || ($map[$n] =~
>         /y-axis/)) {
>                                     shift(@map);
>                                     $n--;
>                             }
>                          }
> 
>                          # There should now be $nres lines left in the file
>                          # The HB map for the last index is written first
>                (top-down in
>                       .xpm file)
>                          #   * Element 0 is index $nres, element 1 is
>         $nres-1, etc.
> 
>                          for (my $i=$nres; $i>=1; $i--) {
>                             # There will be $nframes+2 elements in @line
>         (extra two
>                       are " at
>                          beginning
>                             # and end of the line)
>                             # Establish a conversion factor and split
>         the input
>                lines
>                             my $j = $nres - $i;
>                             my @line = split('', $map[$j]);
> 
>                             # for each index, write to hash
>                             for (my $k=1; $k<=($nframes+1); $k++) {
>                                 if ($line[$k] =~ /o/) {
>                                     $hbonds{$i}++;
>                                 }
>                             }
>                          }
> 
>                          print "Processing the index file...\n";
> 
>                          # Open up the index file and work with it
>                          for (my $n=0; $n<$nres; $n++) {
>                             my @line = split(" ", $index[$n]);
>                             $donors{$n+1} = $line[0];
>                             $acceptors{$n+1} = $line[2];
>                          }
> 
> 
>                          # some arrays for donor and acceptor atom names
>                          my @donor_names;
>                          my @donor_resn;
>                          my @acceptor_names;
>                          my @acceptor_resn;
> 
>                          # Open up the structure file and work with it
>                          print "Processing coordinate file...\n";
>                          foreach $_ (@coord) {
>                             my @line = split(" ", $_);
>                             my $natom = $line[1];
>                             my $name = $line[2];
>                             my $resn = $line[3];
>                             my $resnum = $line[4];
> 
>                             if ($line[0] =~ /ATOM/) {
>                                 unless ($resn =~ /SOL/) {
>                                     for (my $z=1; $z<=$nres; $z++) {
>                                         if ($donors{$z} == $natom) {
>                                             $donor_names[$z] = $name;
>                                             $donor_resn[$z] = join('',
>         $resn,
>                $resnum);
>                                         } elsif ($acceptors{$z} == $natom) {
>                                             $acceptor_names[$z] = $name;
>                                             $acceptor_resn[$z] = join('',
>                $resn, $resnum);
>                                         }
>                                     }
>                                 }
>                             }
>                          }
> 
>                          # open a single output file for writing
>                          open(OUT, ">>summary_HBmap.dat") || die "Cannot
>         open
>                output
>                       file!\n";
>                          printf(OUT "%10s\t%10s\t%10s\t%10s\%10s\n", "#
>                   Donor", " ",
>                          "Acceptor", " ", "% Exist.");
> 
>                          for (my $o=1; $o<=$nres; $o++) {
>                             printf(OUT "%10s\t%10s\t%10s\t%10s\t%10.3f\n",
>                       $donor_resn[$o],
>                          $donor_names[$o], $acceptor_resn[$o],
>         $acceptor_names[$o],
>                          (($hbonds{$o}/$nframes)*100));
>                          }
> 
>                          close(OUT);
> 
>                          exit;
> 
> 
> 
> 
>                          Carla Jamous wrote:
> 
>                              Hi everyone,
> 
>                              I tried to analyze the H-bonds in my
>         trajectory with
>                       g-hbond and
>                              I analysed the xpm and ndx file. But now I
>         need to
>                know the
>                              percentage of existence of each hbond during my
>                       trajectory. Is
>                              there a way to do it with a command line? Or is
>                there a
>                       program
>                              (someone told me there are python programs for
>                analysis of
>                              gromacs trajectories) to extract this
>         information from
>                       the .xpm
>                              file?
> 
>                              Thank you.
> 
>                              Cheers,
>                              Carla
> 
> 
>                          --     ========================================
> 
>                          Justin A. Lemkul
>                          Ph.D. Candidate
>                          ICTAS Doctoral Scholar
>                          MILES-IGERT Trainee
>                          Department of Biochemistry
>                          Virginia Tech
>                          Blacksburg, VA
>                          jalemkul[at]vt.edu <http://vt.edu>
>         <http://vt.edu> <http://vt.edu>
>                <http://vt.edu> | (540)
> 
>                       231-9080
> 
>                        
>          http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin
> 
>                          ========================================
>                          --     gmx-users mailing list  
>          gmx-users at gromacs.org <mailto:gmx-users at gromacs.org>
>                <mailto:gmx-users at gromacs.org <mailto:gmx-users at gromacs.org>>
>                       <mailto:gmx-users at gromacs.org
>         <mailto:gmx-users at gromacs.org> <mailto:gmx-users at gromacs.org
>         <mailto:gmx-users at gromacs.org>>>
>                          <mailto:gmx-users at gromacs.org
>         <mailto:gmx-users at gromacs.org>
>                <mailto:gmx-users at gromacs.org
>         <mailto:gmx-users at gromacs.org>> <mailto:gmx-users at gromacs.org
>         <mailto:gmx-users at gromacs.org>
>                <mailto:gmx-users at gromacs.org
>         <mailto:gmx-users at gromacs.org>>>>
> 
> 
>                          http://lists.gromacs.org/mailman/listinfo/gmx-users
>                          Please search the archive at
>                        
>          http://www.gromacs.org/Support/Mailing_Lists/Search before
>                       posting!
>                          Please don't post (un)subscribe requests to the
>         list.
>                Use the www
>                          interface or send it to
>         gmx-users-request at gromacs.org <mailto:gmx-users-request at gromacs.org>
>                <mailto:gmx-users-request at gromacs.org
>         <mailto:gmx-users-request at gromacs.org>>
>                       <mailto:gmx-users-request at gromacs.org
>         <mailto:gmx-users-request at gromacs.org>
>                <mailto:gmx-users-request at gromacs.org
>         <mailto:gmx-users-request at gromacs.org>>>
>                          <mailto:gmx-users-request at gromacs.org
>         <mailto:gmx-users-request at gromacs.org>
>                <mailto:gmx-users-request at gromacs.org
>         <mailto:gmx-users-request at gromacs.org>>
>                       <mailto:gmx-users-request at gromacs.org
>         <mailto:gmx-users-request at gromacs.org>
>                <mailto:gmx-users-request at gromacs.org
>         <mailto:gmx-users-request at gromacs.org>>>>.
> 
>                          Can't post? Read
>                http://www.gromacs.org/Support/Mailing_Lists
> 
> 
> 
>                   --     ========================================
> 
>                   Justin A. Lemkul
>                   Ph.D. Candidate
>                   ICTAS Doctoral Scholar
>                   MILES-IGERT Trainee
>                   Department of Biochemistry
>                   Virginia Tech
>                   Blacksburg, VA
>                   jalemkul[at]vt.edu <http://vt.edu> <http://vt.edu>
>         <http://vt.edu> | (540)
>                231-9080
>                   http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin
> 
>                   ========================================
>                   --     gmx-users mailing list    gmx-users at gromacs.org
>         <mailto:gmx-users at gromacs.org>
>                <mailto:gmx-users at gromacs.org <mailto:gmx-users at gromacs.org>>
>                   <mailto:gmx-users at gromacs.org
>         <mailto:gmx-users at gromacs.org> <mailto:gmx-users at gromacs.org
>         <mailto:gmx-users at gromacs.org>>>
>                   http://lists.gromacs.org/mailman/listinfo/gmx-users
>                   Please search the archive at
>                   http://www.gromacs.org/Support/Mailing_Lists/Search before
>                posting!
>                   Please don't post (un)subscribe requests to the list.
>         Use the www
>                   interface or send it to gmx-users-request at gromacs.org
>         <mailto:gmx-users-request at gromacs.org>
>                <mailto:gmx-users-request at gromacs.org
>         <mailto:gmx-users-request at gromacs.org>>
>                   <mailto:gmx-users-request at gromacs.org
>         <mailto:gmx-users-request at gromacs.org>
>                <mailto:gmx-users-request at gromacs.org
>         <mailto:gmx-users-request at gromacs.org>>>.
>                   Can't post? Read
>         http://www.gromacs.org/Support/Mailing_Lists
> 
> 
> 
>            --     ========================================
> 
>            Justin A. Lemkul
>            Ph.D. Candidate
>            ICTAS Doctoral Scholar
>            MILES-IGERT Trainee
>            Department of Biochemistry
>            Virginia Tech
>            Blacksburg, VA
>            jalemkul[at]vt.edu <http://vt.edu> <http://vt.edu> | (540)
>         231-9080
>            http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin
> 
>            ========================================
>            --     gmx-users mailing list    gmx-users at gromacs.org
>         <mailto:gmx-users at gromacs.org>
>            <mailto:gmx-users at gromacs.org <mailto:gmx-users at gromacs.org>>
>            http://lists.gromacs.org/mailman/listinfo/gmx-users
>            Please search the archive at
>            http://www.gromacs.org/Support/Mailing_Lists/Search before
>         posting!
>            Please don't post (un)subscribe requests to the list. Use the www
>            interface or send it to gmx-users-request at gromacs.org
>         <mailto:gmx-users-request at gromacs.org>
>            <mailto:gmx-users-request at gromacs.org
>         <mailto:gmx-users-request at gromacs.org>>.
>            Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
> 
> 
> 
>     -- 
>     ========================================
> 
>     Justin A. Lemkul
>     Ph.D. Candidate
>     ICTAS Doctoral Scholar
>     MILES-IGERT Trainee
>     Department of Biochemistry
>     Virginia Tech
>     Blacksburg, VA
>     jalemkul[at]vt.edu <http://vt.edu> | (540) 231-9080
>     http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin
> 
>     ========================================
> 
> 

-- 
========================================

Justin A. Lemkul
Ph.D. Candidate
ICTAS Doctoral Scholar
MILES-IGERT Trainee
Department of Biochemistry
Virginia Tech
Blacksburg, VA
jalemkul[at]vt.edu | (540) 231-9080
http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin

========================================



More information about the gromacs.org_gmx-users mailing list