sábado, 25 de agosto de 2012

Metro en man: medindo un PDB

Xa escribín antes sobre os PDBs, hoxe toca algo máis metido co formato PDB. Un PDB trae unha lista de átomos con un montón de datos sobre a estructura en xeral e cada átomo en particular. Para esta entrada só nos interesan as liñas de átomos, que son as que comezan con "ATOM", "HETATM" ou "ANISOU".

Parseando lixeiramente o archivo

Buscamos estas liñas:

ATOM     10 1HG1 ILE A   1       1.867  30.218  56.564  1.00  0.00

Que según a especificación do PDB,
31 - 38 Son as coordenadas ortogonales de X en Angstroms.
39 - 46 Son as coordenadas ortogonales de Y en Angstroms.
47 - 54 Son as coordenadas ortogonales de Z en Angstroms.

Así que ahí vamos a cazar esas liñas e gardalas nun array.

import numpy as np

atom = ["ATOM", "HETATM", "ANISOU"]
coords = []
n_atoms = 0
for line in pdb:
    if line.split()[0] in atom:
        coords.append(
            [float(s) for s in [line[31:38], line[39:46], line[47:54]]])

#Reformateamos as posicions como unha matriz numpy
positions = np.reshape(coords, (n_atoms, 3))

positions vale algo así como:

[[ 133.927  117.426  173.587]
 [ 133.405  117.299  174.43 ]
 [ 133.622  116.764  172.902]
 ..., 
 [ -27.042   56.598    9.173]
 [ -26.066   57.535    9.316]
 [ -27.61    56.209    8.   ]]

Calculando o tamaño do PDB

Agora simplemente temos que tomar os valores máximos e mínimos en cada eixe, e teremos as medidas do PDB:

box_size = []

for col in range(3):
    box_size.append(abs(positions[:,col].max() - positions[:,col].min()))

print box_size

box_size vale algo así como [220.609, 121.049, 178.267]. As medidas en Angstroms do PDB.