Allow PDBWriter (and similar) to write record_types

[xiki-tempula]

(@richardjgowers edit:) Following #1762 we now read the record type (ATOM/HETATM) of each atom. Now we need to make the PDBWriter (and PQR/PDBQT) write these records when found, otherwise defaulting to a sensible value.

Original Issue:

I have a gro file which contains a ligand and I want to write it out as a pdb file. In the new pdb file can the ligand has a Record type of HETATOM instead of ATOM?

Like

import MDAnalysis as mda
u = mda.Universe(grofile)
u.atoms.recordtype = 'HETATOM'
u.atoms.write('ligand.pdb')

[kain88-de]

Our PDB writer currently doesn't support writing HETATM records. I'm not sure we track this information in the topology either.

[xiki-tempula]

@kain88-de It seems not.
I have loaded a PDB file where the record type is HETATOM and write it into another pdb. It seems that the HETATOM information is lost.

import MDAnalysis as mda
u = mda.Universe('old.pdb')
u.atoms.write('new.pdb')

old pdb
HETATM 1 N LEU 1 32.310 13.778 14.372 1.00 0.00 SYST N 0
new.pdb
ATOM 1 N LEU 1 32.310 13.778 14.372 1.00 0.00 SYST N 0

[kain88-de]

Well yeah for the writing it is lost. Like I said we do not write HETATM records at all. This doesn't mean the topology isn't keeping track of this information internally.

[xiki-tempula]

@kain88-de Thank you for the help. I guess I will do some string manipulation at the end of my script to solve this problem.

[richardjgowers]

So what we could do is add this as a TopologyAttr and read/write it as you'd expect

…

On Thu, 11 Jan 2018, 12:42 p.m. xiki-tempula, ***@***.***> wrote: Closed #1753 <#1753>. — You are receiving this because you are subscribed to this thread. Reply to this email directly, view it on GitHub <#1753 (comment)>, or mute the thread <https://github.com/notifications/unsubscribe-auth/AI0jB3x9ifkTFy6mkoed8JKf1SB2H-dYks5tJgGXgaJpZM4RavIv> .

[mieczyslaw]

I planned to create an issue ticket with the same request so I'm glad that someone else noticed this. I have also noticed that the ticket is already 2 years old.

I can confirm that HETATM is kept in selection, but I'm not sure what's the problem with writing it out?

The simple code to reproduce the issue is as follows (the one in the first message with u.atoms.recordtype = 'HETATOM' crashes because of Cannot set arbitrary attributes to a Group):

from MDAnalysis import Universe
reference = Universe("input.pdb", "input.pdb")
hetatms = reference.select_atoms("record_type HETATM")
hetatms.atoms.write("out.pdb")

[IAlibay]

@mieczyslaw the issue currently is just a case that the PDBWriter as-is doesn't have the correct if/else construct to switch between HETATMs and ATOMs when writing the record lines. It should be relatively simple (from the looks of it maybe ~ 3 lines + tests), unfortunately it's one of those things that get added on top of very long to-do list.

I've added a help wanted tag in case this might incite others to join in and work on this.

[mieczyslaw]

@IAlibay I have just improved PDB.py code and will create a PR in a couple of minutes!

[mieczyslaw]

@IAlibay my first PR created; is this change desired also for PDBQT?

[IAlibay]

@mieczyslaw I think we'll eventually need to fix the PDBQT writer too (although I'm not super clear on where the standard is written out for PDBQT). That being said, the changes here are probably going to be more involved as I'm already seeing a ton of tests that are failing. My suggestion here would be to get #2880 closer to merging before starting on the PDBQT writer (just to avoid code backtracking in case we make decisions that should then be ported from PDB to PDBQT). Either way it's your choice though.

[mieczyslaw]

@mieczyslaw I think we'll eventually need to fix the PDBQT writer too (although I'm not super clear on where the standard is written out for PDBQT). That being said, the changes here are probably going to be more involved as I'm already seeing a ton of tests that are failing. My suggestion here would be to get #2880 closer to merging before starting on the PDBQT writer (just to avoid code backtracking in case we make decisions that should then be ported from PDB to PDBQT). Either way it's your choice though.

Good suggestion, let's focus on #2880.