6

我是python编码的初学者。我正在研究结构坐标。我有具有 xyz 坐标信息的 pdb 结构（最后三个列） 

ATOM      1  N   SER A   1      27.130   7.770  34.390    
ATOM      2  1H  SER A   1      27.990   7.760  34.930     
ATOM      3  2H  SER A   1      27.160   6.960  33.790    
ATOM      4  3H  SER A   1      27.170   8.580  33.790    
ATOM      5  CA  SER A   1      25.940   7.780  35.250    
ATOM      6  CB  SER A   1      25.980   9.090  36.020    
ATOM      7  OG  SER A   1      26.740  10.100  35.320    
ATOM      8  HG  SER A   1      26.750  10.940  35.860    
ATOM      9  C   SER A   1      24.640   7.790  34.460    
ATOM     10  O   SER A   1      24.530   8.510  33.500    
ATOM     11  N   CYS A   2      23.590   7.070  34.760    
ATOM     12  H   CYS A   2      23.590   6.550  35.610    
ATOM     13  CA  CYS A   2      22.420   7.010  33.900    
ATOM     14  CB  CYS A   2      21.620   5.760  34.270    
ATOM     15  SG  CYS A   2      22.480   4.210  33.970    
ATOM     16  C   CYS A   2      21.590   8.220  34.040    
ATOM     17  O   CYS A   2      21.370   8.690  35.160
  • 我的结构中有 1000 个原子。
  • 我有两个疑问。

如何从 xyz 坐标计算结构的质心。
从质心我想画一个半径为 20 厘米的球体。 

I try this


from __future__ import division
import math as mean
import numpy as nx
from string import*


infile = open('file.pdb', 'r')           #open my file
text1 = infile.read().split('\n')
infile.close()

text = []
for i in text1:
if i != '':
    text.append(i)

for j in text:
x1 = eval(replace(j[30:38], ' ', ''))         #extract x-coordinate
y1 = eval(replace(j[38:46], ' ', ''))         #extract y-coordinate
z1 = eval(replace(j[46:54], ' ', ''))         #extract z-coordinate

idcord = []
idcord.append(x1); idcord.append(y1); idcord.append(z1)

centroid = nx.mean(idcord)
print centroid

它给出了每个原子的质心（xyz）我需要一个中心点如何？？？？？？

4

2 回答 2

17

首先，读取文件的更简单方法是使用 numpy 的genfromtxt函数。您不需要导入字符串，也不需要遍历所有行并附加文本或计算字符。

from __future__ import division
import numpy as nx

data = nx.genfromtxt('file.pdb')

然后，最后三列可以访问为：

data[:, -3:]

其中第一个:表示“所有行”，-3:表示从倒数第三列到最后一列。

因此，您可以对它们进行平均：

nx.mean(data[:,-3:], axis=0)

其中axis=0参数告诉nx.mean沿第一个 ( 0th) 轴取平均值。它看起来像这样：

In : data[:,-3:]
Out: 
array([[ 27.13,   7.77,  34.39],
       [ 27.99,   7.76,  34.93],
       [ 27.16,   6.96,  33.79],
       [ 27.17,   8.58,  33.79],
       [ 25.94,   7.78,  35.25],
       [ 25.98,   9.09,  36.02],
       [ 26.74,  10.1 ,  35.32],
       [ 26.75,  10.94,  35.86],
       [ 24.64,   7.79,  34.46],
       [ 24.53,   8.51,  33.5 ],
       [ 23.59,   7.07,  34.76],
       [ 23.59,   6.55,  35.61],
       [ 22.42,   7.01,  33.9 ],
       [ 21.62,   5.76,  34.27],
       [ 22.48,   4.21,  33.97],
       [ 21.59,   8.22,  34.04],
       [ 21.37,   8.69,  35.16]])

In : np.mean(data[:,-3:], axis=0)
Out: array([ 24.74647059,   7.81117647,  34.64823529])

其他一些事情：

1）删除这一行：import math as mean，它导入整个math模块并重命名它mean。您的意图是从模块from math import mean中导入mean功能。math但是在您的代码中，您最终还是使用了( ) 模块中的math函数，因此您从未使用过该版本。numpynxmath

2）您的循环没有缩进，这意味着您要么错误地粘贴到 StackOverflow 中，要么您的循环不正确地缩进。可能这就是您的代码实际的样子：

for j in text:
    x1 = eval(replace(j[30:38], ' ', ''))         #extract x-coordinate
    y1 = eval(replace(j[38:46], ' ', ''))         #extract y-coordinate
    z1 = eval(replace(j[46:54], ' ', ''))         #extract z-coordinate

    idcord = []
    idcord.append(x1); idcord.append(y1); idcord.append(z1)

    centroid = nx.mean(idcord)
    print centroid

但问题是每次idcord循环通过时都会设置为一个空列表，并为每个粒子计算一个新的质心。如果您像上面那样一次导入数据文件，您甚至根本不需要循环。实际上，您的整个代码可以是：

from __future__ import division
import numpy as nx

data = nx.genfromtxt('file.pdb')
nx.mean(data[:,-3:], axis=0)
于 2013-09-10T14:08:12.717 回答
1

尝试这个

import numpy as nx
X = nx.rand(10,3)   # generate some number
centroid = nx.mean(X)
print centroid
于 2013-09-10T09:30:18.777 回答