# NOT RUN {
##-- Basic sequence alignemnt
seqs <- get.seq(c("4q21_A", "1ftn_A"))
aln <- seqaln(seqs)
##-- add a sequence to the (profile) alignment
seq <- get.seq("1tnd_A")
aln <- seqaln(seq, profile=aln)
##-- Read a folder/directory of PDB files
#pdb.path <- "my_dir_of_pdbs"
#files <- list.files(path=pdb.path ,
# pattern=".pdb",
# full.names=TRUE)
##-- Use online files
files <- get.pdb(c("4q21","1ftn"), URLonly=TRUE)
##-- Extract and store sequences
raw <- NULL
for(i in 1:length(files)) {
pdb <- read.pdb(files[i])
raw <- seqbind(raw, pdbseq(pdb) )
}
##-- Align these sequences
aln <- seqaln(raw, id=files, outfile="seqaln.fa")
##-- Read Aligned PDBs storing coordinate data
pdbs <- read.fasta.pdb(aln)
## Sequence identity
seqidentity(aln)
## Note that all the above can be done with the pdbaln() function:
#pdbs <- pdbaln(files)
##-- For identical sequences with masking use a custom matrix
aa <- seqbind(c("X","C","X","X","A","G","K"),
c("C","-","A","X","G","X","X","K"))
aln <- seqaln(aln=aln, id=c("a","b"), outfile="temp.fas", protein=TRUE,
extra.args= paste("-matrix",
system.file("matrices/custom.mat", package="bio3d"),
"-gapopen -3.0 ",
"-gapextend -0.5",
"-center 0.0") )
# }
Run the code above in your browser using DataLab