On Thu, 2014-09-18 at 09:52 +0000, Eric Lang wrote:
> Hi Jason,
>
> Awesome, thanks a lot!
> I am glad it helped to fix this bug.
>
> With the new version of Parmed, I can now convert the PSF and PDB
> files that include my protein and the ligands, which is really great.
>
> However, I have a new problem, but at least this time I have
> identified the cause:
>
> I am working with a big solvated system and thus I have a very large
> number of atoms: more than 99,999...
> Therefore the atom numbering switches to hexadecimal to conserve the PDB column sizes:
> ATOM 99998 H1 TIP3W3580 -31.552 26.681 43.524 1.00 0.00 WT4 H
> ATOM 99999 H2 TIP3W3580 -31.867 26.800 44.896 1.00 0.00 WT4 H
> ATOM 186a0 OH2 TIP3W3583 -14.658 19.993 51.182 1.00 0.00 WT4 O
> ATOM 186a1 H1 TIP3W3583 -15.178 20.681 50.585 1.00 0.00 WT4 H
> ATOM 186a2 H2 TIP3W3583 -14.335 19.304 50.572 1.00 0.00 WT4 H
>
> When I try to use Parmed on the PSF/PDB of my whole solvated system, I have the following error:
>
> parmed_commands.cmdloop()
> File "/usr/lib/python2.7/cmd.py", line 142, in cmdloop
> stop = self.onecmd(line)
> File "/usr/lib/python2.7/cmd.py", line 221, in onecmd
> return func(arg)
> File "<string>", line 1, in <lambda>
> File "/usr/local/lib/python2.7/dist-packages/ParmedTools/parmed_cmd.py", line 141, in _normaldo
> action.execute()
> File "/usr/local/lib/python2.7/dist-packages/ParmedTools/ParmedActions.py", line 3516, in execute
> crd = ChemicalSystem.load_from_pdb(self.crdfile)
> File "/usr/local/lib/python2.7/dist-packages/chemistry/system.py", line 206, in load_from_pdb
> return cls.load_from_open_pdb(open(pdb, 'r'))
> File "/usr/local/lib/python2.7/dist-packages/chemistry/system.py", line 262, in load_from_open_pdb
> altloc=altloc),
> File "/usr/local/lib/python2.7/dist-packages/chemistry/system.py", line 72, in __init__
> super(Atom, self).__init__(*args, **kwargs)
> File "/usr/local/lib/python2.7/dist-packages/chemistry/system.py", line 39, in __init__
> self._types[self._ordered_args[i]](arg))
> ValueError: invalid literal for int() with base 10: '186a0'
>
> So it seems that Parmed do not recognised the hexadecimal numbering.
Ack. VMD is the only program I'm aware of that uses hexadecimal
numbering rather than just, for instance, recycling back to 1. I've
noticed that when psfgen gets its hands on a system, it can do some
wonky things with numbering. If you don't use psfgen and just save the
PDB, then the residue number does cycle using hexadecimal. My test case
had a system that was large enough to cycle residues, but not large
enough to cycle atom indexes, which is why it didn't handle your case.
This should be fixed now. Go ahead and try the updated version.
>
> I had a quick look at system.py and I saw that you have a test to see if the residue number is hexadecimal, but I don't think I saw something for the atom number. Do you think this can be easily fixed? If not is there any other way to get round of this problem?
>
> By the way regarding the hexadecimal test for the residue number, I had a look at my PDB files, and I saw that VMD (or PSFGEN or SOLVATE, I am not sure which one), do not switch to hexadecimal numbering after 9999 for the residue number, instead it restarts the numbering from 1 but changes the "segname" number (see for example below where the segname change from WAT1 to WAT2).
>
> ATOM 51835 OH2 TIP3W9999 -9.810 2.513 -32.630 1.00 0.00 WAT1 O
> ATOM 51836 H1 TIP3W9999 -9.976 2.772 -33.536 1.00 0.00 WAT1 H
> ATOM 51837 H2 TIP3W9999 -9.078 3.060 -32.345 1.00 0.00 WAT1 H
> ATOM 51838 OH2 TIP3W 1 -13.494 -46.453 -25.565 1.00 0.00 WAT2 O
> ATOM 51839 H1 TIP3W 1 -13.825 -46.812 -24.741 1.00 0.00 WAT2 H
> ATOM 51840 H2 TIP3W 1 -13.746 -45.531 -25.553 1.00 0.00 WAT2 H
>
> Later in the PDB file, VMD doesn't wait to reach 9999 to start a new segname and do not necessary restart at 1, see for example:
>
> ATOM 1958b OH2 TIP3W9260 0.636 11.509 56.128 1.00 0.00 WT4 O
> ATOM 1958c H1 TIP3W9260 1.069 10.681 55.969 1.00 0.00 WT4 H
> ATOM 1958d H2 TIP3W9260 0.543 11.679 57.135 1.00 0.00 WT4 H
> ATOM 1958e OH2 TIP3W 2 13.925 6.927 -52.530 1.00 0.00 WT5 O
> ATOM 1958f H1 TIP3W 2 14.107 6.026 -52.541 1.00 0.00 WT5 H
> ATOM 19590 H2 TIP3W 2 13.089 6.988 -52.946 1.00 0.00 WT5 H
> ATOM 19591 OH2 TIP3W 3 22.542 -1.590 -54.002 1.00 0.00 WT5 O
> ATOM 19592 H1 TIP3W 3 22.312 -1.568 -53.025 1.00 0.00 WT5 H
> ATOM 19593 H2 TIP3W 3 21.824 -1.187 -54.404 1.00 0.00 WT5 H
>
> So I was just wondering if this matters or not for the conversion of
> the PSF-PDB files into INPCRD-PRMTOP files and if so, if Parmed takes
> this into account?
>
> Last question: I cannot open in VMD the prmtop file generated by
> Parmed. I have noticed that with Chamber one can add the -vmd flag to
> generate "vmd-friendly" prmtop files. Does Parmed do the same? If not
> how can I visualize it?
Use the PSF file ;). I didn't bother implementing the -vmd flag because
VMD has no trouble using the PSF file instead (and the atom ordering
between the PSF and prmtop is guaranteed to be identical -- at least for
ParmEd).
I don't like the idea of having a VMD-compatible chamber topology since
it will actually run with sander or pmemd and give reasonable answers...
it just won't implement the parts of the CHARMM force field that are not
in the AMBER force field (e.g., impropers, urey-bradleys, 1-4 vdW
modifications, and CMAPs). And this is not a mistake that only
inexperienced users can make (I can easily envision myself making this
mistake), so I do what I can to minimize its likelihood.
You may be able to change my mind if you can show me a use-case for the
vmd-compatible prmtop that is not satisfied equally well with the PSF
file :).
Thanks for your help getting these issues worked out!
Jason
--
Jason M. Swails
BioMaPS,
Rutgers University
Postdoctoral Researcher
_______________________________________________
AMBER mailing list
AMBER.ambermd.org
http://lists.ambermd.org/mailman/listinfo/amber
Received on Thu Sep 18 2014 - 06:00:04 PDT
