0

我正在尝试使用该atomselect命令删除.gro文件中蛋白质区域内的脂质。由于它是一个马提尼粗粒度文件，我可以将关键字resname用于残基名称和name/或type珠类型（我的伪原子）。所以没有定义默认的单字。

我用命令试了一下：

atomselect 0 "all not resname DPPE DOPE POPE POPG within 1 of resname ALA ARG ASN ASP CYS GLN GLU GLY HIS ILE LEU LYS MET PHE PRO SER THR TRP TYR VAL"

并得到以下错误：

错误）选择过早终止

ERROR) 语法错误 atomselect: 无法解析选择文本: all not resname DPPE DOPE POPE POPG 在 resname 的 1 内 ALA ARG ASN ASP CYS GLN GLU GLY HIS ILE LEU LYS MET PHE PRO SER THR TRP TYR VAL

所以显然我没有得到正确的语法。我尝试了几个不同的版本并将resname选择存储在变量中，但没有任何效果。我该如何解决这个问题？

4

1 回答 1

0

因此，我使用 bio3d 在 R 中编写了一些代码来完成任务：

library(bio3d)
library (gdata)

#read pdb input
insane <- read.pdb("nice.pdb")


#vector with residuenumbers to delete 
select_resi <- c(1000, 1136, 1026, 1252, 1449, 970, 1067, 1298, 
                 1287, 1357, 1051, 993, 1241, 1282, 1341, 1344, 
                 1048, 1154, 1205, 1274, 1465, 1322, 1418, 992)

select_resi <- sort(select_resi)

#select residues to delete
selected_resi <-  atom.select(insane, resno = select_resi)


#delete selected residues
insane$atom <- insane$atom[-selected_resi$atom,]
insane$calpha <- insane$calpha[-selected_resi$atom]
insane$xyz <- insane$xyz[-selected_resi$xyz]

#renumber residuenumbers and convert to gromacs type pdb
printstuff <- convert.pdb(insane, type = "gromacs", 
                          renumber = TRUE, first.resno = 1, first.eleno = 1)

#write pdb file
write.pdb(pdb = printstuff, file = "memb2R.pdb", xyz = printstuff$xyz, 
          type = printstuff$atom$type, resno = printstuff$atom$resno, 
          resid = printstuff$atom$resid, eleno = printstuff$atom$eleno, 
          elety = printstuff$atom$elety, end = TRUE, verbose = TRUE)

如果您.gro之后需要该格式（这就是加载 gdata 的原因），您可以继续执行此操作：

#preparing xyz to match with .gro format
xyz_vector <- insane$xyz/10 

pos_x <- c()
pos_y <- c()
pos_z <- c()

#sorting loop for x, y, z coordinates
for (i in seq(along=xyz_vector)) {
  pos_x <- c(pos_x, xyz_vector[i])
  xyz_vector <- xyz_vector[-i]
    pos_y <- c(pos_y, xyz_vector[i])
    xyz_vector <- xyz_vector[-i]
      pos_z <- c(pos_z, xyz_vector[i])
      print(i)
}

#delete redundant entries
pos_x <- pos_x[1:length(xyz_vector)]
pos_y <- pos_y[1:length(xyz_vector)]
pos_z <- pos_z[1:length(xyz_vector)]

#prepare other passing vectors for .gro vector
resi_numb <- insane$atom$resno
resi_name <- insane$atom$resid
atom_name <- insane$atom$elety
atom_numb <- insane$atom$eleno


#prepare .gro vector
gro_vec <- sprintf ("%5d%-5s%5s%5d%8.3f%8.3f%8.3f",

                      #for velocity fields add %8.4f%8.4f%8.4f to string and create 3 velocity vectors similar to pos_x etc.

                      resi_numb, resi_name, atom_name, atom_numb,
                      pos_x, pos_y, pos_z)

#transform gro_vec in matrix
gro_matrix <- as.matrix(gro_vec, byrow = TRUE)


#write output
write.fwf (gro_matrix, file = "yourfile.txt", sep = "")

记得将.txt文件重命名为.gro. 在文件的第一行写入系统名称，第二行写入原子实体，最后一行写入框向量（如果有）。不是最漂亮的代码，但我还是新手，它可以完成工作。
随时欢迎改进建议:)

于 2015-09-15T14:54:58.317 回答