Read a Protein Data Bank (PDB) coordinate file.

read.pdb(file, maxlines = -1, multi = FALSE, rm.insert = FALSE,
         rm.alt = TRUE, ATOM.only = FALSE, hex = FALSE, verbose = TRUE)

read.pdb2(file, maxlines = -1, multi = FALSE, rm.insert = FALSE,
         rm.alt = TRUE, ATOM.only = FALSE, verbose = TRUE)

# S3 method for pdb
print(x, printseq=TRUE, ...)

# S3 method for pdb
summary(object, printseq=FALSE, ...)

Arguments

file

a single element character vector containing the name of the PDB file to be read, or the four letter PDB identifier for online file access.

maxlines

the maximum number of lines to read before giving up with large files. By default if will read up to the end of input on the connection.

multi

logical, if TRUE multiple ATOM records are read for all models in multi-model files and their coordinates returned.

rm.insert

logical, if TRUE PDB insert records are ignored.

rm.alt

logical, if TRUE PDB alternate records are ignored.

ATOM.only

logical, if TRUE only ATOM/HETATM records are stored. Useful for speed enhancements with large files where secondary structure, biological unit and other remark records are not required.

hex

logical, if TRUE enable parsing of hexadecimal atom numbers (> 99.999) and residue numbers (> 9.999) (e.g. from VMD). Note that numbering is assumed to be consecutive (with no missing numbers) and the hexadecimals should start at atom number 100.000 and residue number 10.000 and proceed to the end of file.

verbose

print details of the reading process.

x

a PDB structure object obtained from read.pdb.

object

a PDB structure object obtained from read.pdb.

printseq

logical, if TRUE the PDB ATOM sequence will be printed to the screen. See also pdbseq.

...

additional arguments to ‘print’.

Details

read.pdb is a re-implementation (using Rcpp) of the slower but more tested R implementation of the same function (called read.pdb2 since bio3d-v2.3).

maxlines may be set so as to restrict the reading to a portion of input files. Note that the preferred means of reading large multi-model files is via binary DCD or NetCDF format trajectory files (see the read.dcd and read.ncdf functions).

Value

Returns a list of class "pdb" with the following components:

atom

a data.frame containing all atomic coordinate ATOM and HETATM data, with a row per ATOM/HETATM and a column per record type. See below for details of the record type naming convention (useful for accessing columns).

helix

‘start’, ‘end’ and ‘length’ of H type sse, where start and end are residue numbers “resno”.

sheet

‘start’, ‘end’ and ‘length’ of E type sse, where start and end are residue numbers “resno”.

seqres

sequence from SEQRES field.

xyz

a numeric matrix of class "xyz" containing the ATOM and HETATM coordinate data.

calpha

logical vector with length equal to nrow(atom) with TRUE values indicating a C-alpha “elety”.

remark

a list object containing information taken from 'REMARK' records of a "pdb". It can be used for building biological units (See biounit).

call

the matched call.

References

Grant, B.J. et al. (2006) Bioinformatics 22, 2695--2696.

For a description of PDB format (version3.3) see:
http://www.wwpdb.org/documentation/format33/v3.3.html.

Author

Barry Grant

Note

For both atom and het list components the column names can be used as a convenient means of data access, namely: Atom serial number “eleno” , Atom type “elety”, Alternate location indicator “alt”, Residue name “resid”, Chain identifier “chain”, Residue sequence number “resno”, Code for insertion of residues “insert”, Orthogonal coordinates “x”, Orthogonal coordinates “y”, Orthogonal coordinates “z”, Occupancy “o”, and Temperature factor “b”. See examples for further details.

See also

atom.select, write.pdb, trim.pdb, cat.pdb, read.prmtop, as.pdb, read.dcd, read.ncdf, read.fasta.pdb, read.fasta, biounit

Examples

## Read a PDB file from the RCSB online database #pdb <- read.pdb("4q21") ## Read a PDB file from those included with the package pdb <- read.pdb( system.file("examples/1hel.pdb", package="bio3d") ) ## Print a brief composition summary pdb
#> #> Call: read.pdb(file = system.file("examples/1hel.pdb", package = "bio3d")) #> #> Total Models#: 1 #> Total Atoms#: 1001, XYZs#: 3003 Chains#: 1 (values: A) #> #> Protein Atoms#: 1001 (residues/Calpha atoms#: 129) #> Nucleic acid Atoms#: 0 (residues/phosphate atoms#: 0) #> #> Non-protein/nucleic Atoms#: 0 (residues: 0) #> Non-protein/nucleic resid values: [ none ] #> #> Protein sequence: #> KVFGRCELAAAMKRHGLDNYRGYSLGNWVCAAKFESNFNTQATNRNTDGSTDYGILQINS #> RWWCNDGRTPGSRNLCNIPCSALLSSDITASVNCAKKIVSDGNGMNAWVAWRNRCKGTDV #> QAWIRGCRL #> #> + attr: atom, xyz, seqres, helix, sheet, #> calpha, remark, call
## Examine the storage format (or internal *str*ucture) str(pdb)
#> List of 8 #> $ atom :'data.frame': 1001 obs. of 16 variables: #> ..$ type : chr [1:1001] "ATOM" "ATOM" "ATOM" "ATOM" ... #> ..$ eleno : int [1:1001] 1 2 3 4 5 6 7 8 9 10 ... #> ..$ elety : chr [1:1001] "N" "CA" "C" "O" ... #> ..$ alt : chr [1:1001] NA NA NA NA ... #> ..$ resid : chr [1:1001] "LYS" "LYS" "LYS" "LYS" ... #> ..$ chain : chr [1:1001] "A" "A" "A" "A" ... #> ..$ resno : int [1:1001] 1 1 1 1 1 1 1 1 1 2 ... #> ..$ insert: chr [1:1001] NA NA NA NA ... #> ..$ x : num [1:1001] 3.294 2.388 2.438 2.406 0.949 ... #> ..$ y : num [1:1001] 10.2 10.5 12 12.9 10.1 ... #> ..$ z : num [1:1001] 10.27 9.17 8.89 9.81 9.56 ... #> ..$ o : num [1:1001] 1 1 1 1 1 1 1 1 1 1 ... #> ..$ b : num [1:1001] 11.18 9.68 14 14 13.29 ... #> ..$ segid : chr [1:1001] NA NA NA NA ... #> ..$ elesy : chr [1:1001] "N" "C" "C" "O" ... #> ..$ charge: chr [1:1001] NA NA NA NA ... #> $ xyz : 'xyz' num [1, 1:3003] 3.29 10.16 10.27 2.39 10.53 ... #> $ seqres: Named chr [1:129] "LYS" "VAL" "PHE" "GLY" ... #> ..- attr(*, "names")= chr [1:129] "A" "A" "A" "A" ... #> $ helix :List of 4 #> ..$ start: Named num [1:7] 4 25 79 88 104 109 120 #> .. ..- attr(*, "names")= chr [1:7] "" "" "" "" ... #> ..$ end : Named num [1:7] 16 34 84 100 108 114 124 #> .. ..- attr(*, "names")= chr [1:7] "" "" "" "" ... #> ..$ chain: chr [1:7] "A" "A" "A" "A" ... #> ..$ type : chr [1:7] "1" "1" "5" "1" ... #> $ sheet :List of 4 #> ..$ start: Named num [1:3] 42 49 56 #> .. ..- attr(*, "names")= chr [1:3] "" "" "" #> ..$ end : Named num [1:3] 46 54 60 #> .. ..- attr(*, "names")= chr [1:3] "" "" "" #> ..$ chain: chr [1:3] "A" "A" "A" #> ..$ sense: chr [1:3] "0" "-1" "-1" #> $ calpha: logi [1:1001] FALSE TRUE FALSE FALSE FALSE FALSE ... #> $ remark:List of 1 #> ..$ biomat:List of 4 #> .. ..$ num : int 1 #> .. ..$ chain :List of 1 #> .. .. ..$ : chr "A" #> .. ..$ mat :List of 1 #> .. .. ..$ :List of 1 #> .. .. .. ..$ A: num [1:3, 1:4] 1 0 0 0 1 0 0 0 1 0 ... #> .. ..$ method: chr "AUTHOR" #> $ call : language read.pdb(file = system.file("examples/1hel.pdb", package = "bio3d")) #> - attr(*, "class")= chr [1:2] "pdb" "sse"
## Print data for the first four atom pdb$atom[1:4,]
#> type eleno elety alt resid chain resno insert x y z o b #> 1 ATOM 1 N <NA> LYS A 1 <NA> 3.294 10.164 10.266 1 11.18 #> 2 ATOM 2 CA <NA> LYS A 1 <NA> 2.388 10.533 9.168 1 9.68 #> 3 ATOM 3 C <NA> LYS A 1 <NA> 2.438 12.049 8.889 1 14.00 #> 4 ATOM 4 O <NA> LYS A 1 <NA> 2.406 12.898 9.815 1 14.00 #> segid elesy charge #> 1 <NA> N <NA> #> 2 <NA> C <NA> #> 3 <NA> C <NA> #> 4 <NA> O <NA>
## Print some coordinate data head(pdb$atom[, c("x","y","z")])
#> x y z #> 1 3.294 10.164 10.266 #> 2 2.388 10.533 9.168 #> 3 2.438 12.049 8.889 #> 4 2.406 12.898 9.815 #> 5 0.949 10.101 9.559 #> 6 -0.050 10.621 8.573
## Or coordinates as a numeric vector #head(pdb$xyz) ## Print C-alpha coordinates (can also use 'atom.select' function) head(pdb$atom[pdb$calpha, c("resid","elety","x","y","z")])
#> resid elety x y z #> 2 LYS CA 2.388 10.533 9.168 #> 11 VAL CA 2.524 13.840 7.282 #> 18 PHE CA -1.021 15.173 7.169 #> 29 GLY CA -2.645 17.273 4.511 #> 33 ARG CA -4.106 20.560 5.674 #> 44 CYS CA -7.767 19.572 5.189
inds <- atom.select(pdb, elety="CA") head( pdb$atom[inds$atom, ] )
#> type eleno elety alt resid chain resno insert x y z o b #> 2 ATOM 2 CA <NA> LYS A 1 <NA> 2.388 10.533 9.168 1 9.68 #> 11 ATOM 11 CA <NA> VAL A 2 <NA> 2.524 13.840 7.282 1 10.02 #> 18 ATOM 18 CA <NA> PHE A 3 <NA> -1.021 15.173 7.169 1 15.34 #> 29 ATOM 29 CA <NA> GLY A 4 <NA> -2.645 17.273 4.511 1 13.97 #> 33 ATOM 33 CA <NA> ARG A 5 <NA> -4.106 20.560 5.674 1 11.63 #> 44 ATOM 44 CA <NA> CYS A 6 <NA> -7.767 19.572 5.189 1 12.93 #> segid elesy charge #> 2 <NA> C <NA> #> 11 <NA> C <NA> #> 18 <NA> C <NA> #> 29 <NA> C <NA> #> 33 <NA> C <NA> #> 44 <NA> C <NA>
## The atom.select() function returns 'indices' (row numbers) ## that can be used for accessing subsets of PDB objects, e.g. inds <- atom.select(pdb,"ligand") pdb$atom[inds$atom,]
#> [1] type eleno elety alt resid chain resno insert x y #> [11] z o b segid elesy charge #> <0 rows> (or 0-length row.names)
pdb$xyz[inds$xyz]
#> numeric(0)
## See the help page for atom.select() function for more details. if (FALSE) { ## Print SSE data for helix and sheet, ## see also dssp() and stride() functions print.sse(pdb) pdb$helix pdb$sheet$start ## Print SEQRES data pdb$seqres ## SEQRES as one letter code aa321(pdb$seqres) ## Where is the P-loop motif in the ATOM sequence inds.seq <- motif.find("G....GKT", pdbseq(pdb)) pdbseq(pdb)[inds.seq] ## Where is it in the structure inds.pdb <- atom.select(pdb,resno=inds.seq, elety="CA") pdb$atom[inds.pdb$atom,] pdb$xyz[inds.pdb$xyz] ## View in interactive 3D mode #view(pdb) }