2
$\begingroup$

I'm using R and its Rpdb package to work with PDB files. I use the read.pdb function but, nevertheless, I couldn't find the argument which will allow me to view the SEQRES section and the full amino-acid sequence of the protein.

For example, I want to find which amino-acid is in the 17th position of the sequence and I want to extract it from the 'SRQRES' section in the PDB file (the original length of the sequence was reduced for the sake of the example):

SEQRES   1 A  390  GLU PRO GLU TRP THR TYR PRO ARG LEU SER CYS GLN GLY          
SEQRES   2 A  390  SER THR PHE GLN LYS ALA LEU LEU ILE SER PRO HIS ARG          
SEQRES   3 A  390  PHE GLY GLU ALA ARG GLY ASN SER ALA PRO LEU ILE ILE          
SEQRES   4 A  390  ARG GLU PRO PHE ILE ALA CYS GLY PRO LYS GLU CYS LYS          
SEQRES   5 A  390  HIS PHE ALA LEU THR HIS TYR ALA ALA GLN PRO GLY GLY          
SEQRES   6 A  390  TYR TYR ASN GLY THR ARG GLU ASP ARG ASN LYS LEU ARG          
SEQRES   7 A  390  HIS LEU ILE SER VAL LYS LEU GLY LYS ILE PRO THR VAL          
SEQRES   8 A  390  GLU ASN SER ILE PHE HIS MET ALA ALA TRP SER GLY SER          
SEQRES   9 A  390  ALA CYS HIS ASP GLY ARG GLU TRP THR TYR ILE GLY VAL          
SEQRES  10 A  390  ASP GLY PRO ASP SER ASN ALA LEU ILE LYS ILE LYS TYR          

...

The original file can be found here 1a4g.pdb.

Any ideas on how to accomplish that?

$\endgroup$
2
  • $\begingroup$ Please edit your question and add a minimal, reproducible example we can use to test our answers. Ideally, give us a PDB file we can use and the code you are trying that fails. $\endgroup$
    – terdon
    Commented Jan 27, 2018 at 20:01
  • 1
    $\begingroup$ Generally libraries that read PDB files interpret only a subset of records from the PDB file. I don't know Rpdb, but if it can't read SEQRES then, well, it can't. You can either code it yourself or look for another library. $\endgroup$
    – marcin
    Commented Jan 27, 2018 at 21:15

1 Answer 1

4
$\begingroup$

I have not used Rpdb but I have used, and would recommend, bio3d.

To read the PDB file 1A4G:

library(bio3d)
pdb1a4g <- read.pdb("1a4g")

To get the SEQRES information, there are two options. Get it from the pdb object directly:

seq1a4g <- pdb1a4g$seqres

Or call pdbseq:

seq1a4g <- pdbseq(pdb1a4g)

Using the second method, the 17th position (note this is not the same as the 17th position of the original protein sequence):

seq1a4g[17]
 92 
"Q" 
$\endgroup$
4
  • $\begingroup$ Thank you. How come pdbfile$seqres and pdbseq(pdbfile) aren't the same? Which one is more reliable? $\endgroup$ Commented Feb 4, 2018 at 11:22
  • $\begingroup$ The first one is just accessing a property of the pdb object. Think of it as having a variable attached which you can read, a bit like a column of a data frame. The second is a function which has been written to access, parse and format the data. So they're both "reliable", you choose whichever is easier to work with for your use case. $\endgroup$
    – neilfws
    Commented Feb 4, 2018 at 11:53
  • $\begingroup$ Thanks, but I've noticed that both functions return a different vector (as can be seen with "1d1v.pdb" file for example). How come? $\endgroup$ Commented Feb 4, 2018 at 12:00
  • 1
    $\begingroup$ If you look at the PDB record for 1D1V, there's a section titled "MISSING RESIDUES". It tells you that residues 39-66 are missing from the structure for chain A. So pdbseq() shows you the residues present in the structure, starting from G67. $seqres returns the complete sequence that you'd see in the FASTA record, including the missing residues. $\endgroup$
    – neilfws
    Commented Feb 4, 2018 at 21:13

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.