Re: psf generation for multiple trehalose molecules

From: Ritu Arora (ritu.arora_at_mail.concordia.ca)
Date: Mon Nov 09 2020 - 10:33:55 CST

Perfect! Problem is solved. Thank you so much, Josh.

Sincerely,
Ritu
________________________________
From: Josh Vermaas <yoshi121212_at_gmail.com> on behalf of Josh Vermaas <joshua.vermaas_at_gmail.com>
Sent: Monday, November 9, 2020 10:40 AM
To: Ritu Arora <ritu.arora_at_mail.concordia.ca>; Peter Freddolino <petefred_at_umich.edu>; Namd Mailing List <namd-l_at_ks.uiuc.edu>
Subject: Re: namd-l: psf generation for multiple trehalose molecules

Hi Ritu,

Ok. Then it comes down to a matching problem. What psfgen does it is looks for the first atoms that match the specification of what is in the topology so far. So since each of your trehalose molecules is in its own segment, and you are passing the same pdb file every time (but with renamed segments), resid 1 and 2 *always* matches the first two residues that are resid 1 and 2 in the pdb. What you may want to do is pretreat your pdb so that every two residue pair has a distinct segname that matches what you'll be assigning later. The edited commands should match correctly.

-Josh

mol new output/carb_all.pdb

for { set i 0 } { $i < 177 } { incr i } {
set resid1 [expr {$i * 2}]
set resid2 [expr {$i*2 + 1}]
set sel [atomselect top "residue $resid1 $resid2"]
segment L$i {
residue 1 AGLC
residue 2 AGLC
      }
      patch 11aa L$i:1 L$i:2
      $sel set segname L$i
      $sel writepdb tmp.pdb
      regenerate angles dihedrals
      coordpdb tmp.pdb L$i
       $sel delete
 }

On 11/9/20 8:23 AM, Ritu Arora wrote:
Hi Josh,

I am sorry for providing insufficient information about my system. My system comprises of a protein, chloride ions, water, and carbohydrate molecules. The psfgen script worked fine for all the components except trehalose. I made my input using packmol. VMD autopsf builder worked fine for monosaccharides but fails to identify all the disaccharide molecules and its glycosidic linkages. So, psfgen was the final resort for my system having disaccharides.

If you wanted to do a cube, what you could do is place the centers on a 5x5x6 grid and move the selection accordingly.
I want the script to read the coordinates from my initial pdb file generated using packmol (file attached). I do not want to randomly position the carbohydrate molecules in a cube as this can cause clashes with other atoms of my system. I apologize for not being clear at first. This is the psfgen script for the whole system:

mkdir -p output

grep 'PROA' 1rex_tre.pdb > output/protein.pdb
grep 'HOH' 1rex_tre.pdb > output/water.pdb
grep 'CLA CLA' 1rex_tre.pdb > output/cla.pdb
grep 'CARB' 1rex_tre.pdb > output/carb_all.pdb

psfgen << ENDMOL

topology tre-mix.top

segment PP1 { pdb output/protein.pdb}

patch DISU PP1:6 PP1:128
patch DISU PP1:30 PP1:116
patch DISU PP1:65 PP1:81
patch DISU PP1:77 PP1:95

coordpdb output/protein.pdb PP1

pdbalias residue HOH TIP3
segment SOLV {
auto none
pdb output/water.pdb
}

pdbalias atom HOH O OH2
coordpdb output/water.pdb SOLV

segment CLA {pdb output/cla.pdb }
coordpdb output/cla.pdb CLA

mol new output/carb_all.pdb
set sel [atomselect top "all"]
for { set i 0 } { $i < 177 } { incr i } {
segment L$i {
residue 1 AGLC
residue 2 AGLC
      }
      patch 11aa L$i:1 L$i:2
      $sel set segname L$i
      $sel writepdb tmp.pdb
      regenerate angles dihedrals
      coordpdb tmp.pdb L$i
 }

guesscoord
writepsf output/1rex_tre_n.psf
writepdb output/1rex_tre_n.pdb
ENDMOL

I wonder why the script can't read the coordinates of disaccharides like other molecules or solvents of my system. I already spent a lot of time in modifying the script at my end, all in vain. Please help!

Sincerely,
Ritu

________________________________
From: Josh Vermaas <yoshi121212_at_gmail.com><mailto:yoshi121212_at_gmail.com> on behalf of Josh Vermaas <joshua.vermaas_at_gmail.com><mailto:joshua.vermaas_at_gmail.com>
Sent: Monday, November 9, 2020 9:34 AM
To: Ritu Arora <ritu.arora_at_mail.concordia.ca><mailto:ritu.arora_at_mail.concordia.ca>; Peter Freddolino <petefred_at_umich.edu><mailto:petefred_at_umich.edu>
Subject: Re: namd-l: psf generation for multiple trehalose molecules

Hi Ritu,

Your script is moving each segment by 1 Angstrom in each of x, y, and z. The new.pdb you are showing looks like what you'd expect from that result, since it forms a line where the atoms are close together. You need to ask what you *want* the result to look like, and edit the script accordingly. If you wanted to do a cube, what you could do is place the centers on a 5x5x6 grid and move the selection accordingly.

-Josh

On 11/9/20 3:22 AM, Ritu Arora wrote:
Dear Peter,

Please find attached the input and output pdb files, just in case.

I look forward to your reply. Many thanks for your kind suggestions.

Best,
Ritu
________________________________
From: Ritu Arora <ritu.arora_at_mail.concordia.ca><mailto:ritu.arora_at_mail.concordia.ca>
Sent: Monday, November 9, 2020 5:12 AM
To: Peter Freddolino <petefred_at_umich.edu><mailto:petefred_at_umich.edu>; <namd-l_at_ks.uiuc.edu><mailto:namd-l_at_ks.uiuc.edu> <namd-l_at_ks.uiuc.edu><mailto:namd-l_at_ks.uiuc.edu>
Cc: Josh Vermaas <joshua.vermaas_at_gmail.com><mailto:joshua.vermaas_at_gmail.com>; namd_at_ks.uiuc.edu<mailto:namd_at_ks.uiuc.edu> <namd_at_ks.uiuc.edu><mailto:namd_at_ks.uiuc.edu>
Subject: Re: namd-l: psf generation for multiple trehalose molecules

Hi Peter,

Omitting 'resetpsf' inside the loop helped. 'tmp.pdb' file looks fine, but all the molecules get stacked on top of each other in 'new.pdb' file.

Ritu
________________________________
From: Peter Freddolino <petefred_at_umich.edu><mailto:petefred_at_umich.edu>
Sent: Monday, November 9, 2020 1:08 AM
To: <namd-l_at_ks.uiuc.edu><mailto:namd-l_at_ks.uiuc.edu> <namd-l_at_ks.uiuc.edu><mailto:namd-l_at_ks.uiuc.edu>; Ritu Arora <ritu.arora_at_mail.concordia.ca><mailto:ritu.arora_at_mail.concordia.ca>
Cc: Josh Vermaas <joshua.vermaas_at_gmail.com><mailto:joshua.vermaas_at_gmail.com>; namd_at_ks.uiuc.edu<mailto:namd_at_ks.uiuc.edu> <namd_at_ks.uiuc.edu><mailto:namd_at_ks.uiuc.edu>
Subject: Re: namd-l: psf generation for multiple trehalose molecules

I suspect that the 'resetpsf' inside your loop isn't helping.
Best,
Peter

On Sun, Nov 8, 2020 at 10:53 PM Ritu Arora <ritu.arora_at_mail.concordia.ca<mailto:ritu.arora_at_mail.concordia.ca>> wrote:
Hi Josh,

Thank you for clarifying this. After removing alias commands and minor typos, the following script saves the coordinates and structure of only the last trehalose molecule in the output files. Something is still missing!

package require psfgen
topology top_all36_carb.rtf
mol new carb_all.pdb
set sel [atomselect top "all"]
for { set i 0 } { $i < 150 } { incr i } {
resetpsf
segment L$i {
residue 1 AGLC
residue 2 AGLC
}
patch 11aa L$i:1 L$i:2
$sel set segname L$i
$sel moveby [list 1 1 1]
$sel writepdb tmp.pdb
regenerate angles dihedrals
coordpdb tmp.pdb L$i
}
guesscoord
writepdb new.pdb
writepsf new.psf

Sincerely,
Ritu
________________________________
From: Josh Vermaas <joshua.vermaas_at_gmail.com<mailto:joshua.vermaas_at_gmail.com>>
Sent: Sunday, November 8, 2020 9:13 PM
To: NAMD list <namd-l_at_ks.uiuc.edu<mailto:namd-l_at_ks.uiuc.edu>>; Ritu Arora <ritu.arora_at_mail.concordia.ca<mailto:ritu.arora_at_mail.concordia.ca>>
Cc: namd_at_ks.uiuc.edu<mailto:namd_at_ks.uiuc.edu> <namd_at_ks.uiuc.edu<mailto:namd_at_ks.uiuc.edu>>
Subject: Re: namd-l: psf generation for multiple trehalose molecules

Hi Ritu,

Why are you aliasing at all? In patches that connect more than one residue, the first digit indicates to psfgen which residue the atom name is on. This is fundamentally why atom names can't start with a digit. I think your script otherwise looks close to what I would write.

Josh

On Sun, Nov 8, 2020, 5:26 PM Ritu Arora <ritu.arora_at_mail.concordia.ca<mailto:ritu.arora_at_mail.concordia.ca>> wrote:
Dear All,

I wish to run a simulation of 150 molecules of trehalose. I adapted a psfgen script from earlier NAMD threads, but I am not able to get the script working. The trial script for the carbohydrate section of my system follows:

package require psfgen
topology top_all36_carb.rtf
mol new carb_only.pdb
set sel [atomselect top 'all']
for {set i 0} {$i<150} {incr i} {
resetpsf
segment L$i {
residue 1 AGLC
residue 2 AGLC
}
pdbalias atom 'AGLC 1' O1 1O1
pdbalias atom 'AGLC 1' C1 1C1
pdbalias atom 'AGLC 2' C1 2C1
#
patch 11aa L$i:1 L$i:2
$sel set segname L$i
$set moveby [list 1 1 1]
$set writepdb tmp.pdb
regenerate angles dihederals
coordpdb L$i tmp.pdb
}
guesscoord
writepdb $name.pdb
writepsf $name.psf

The confusion is that there are two AGLCs in trehalose and each has either resid 1 or 2, rather than sequential increments for whole molecule of trehalose. The pdbalias command also seems to be wrong.

I also tried to split the two fragments of trehalose in two separate pdbs:

pdbalias atom 'AGLC 1' O1 1O1
pdbalias atom 'AGLC 2' C1 2C1

for {set i 0} {$i < 150} {incr i} {
segment TRE1$i {pdb output/carb1.pdb }

segment TRE2$i {pdb output/carb2.pdb }

patch 11aa TRE1$i:1 TRE2$i:2
coordpdb output/carb1.pdb TRE1$i
coordpdb output/carb2.pdb TRE2$i
}
I would be thankful for any suggestions to improvise the script.

Best,
Ritu

This archive was generated by hypermail 2.1.6 : Thu Dec 31 2020 - 23:17:14 CST