bio3d (version 2.1-3)

pca.xyz: Principal Component Analysis

Description

Performs principal components analysis (PCA) on a xyz numeric data matrix.

Usage

## S3 method for class 'xyz':
pca(xyz, subset = rep(TRUE, nrow(as.matrix(xyz))),
                 use.svd = FALSE, rm.gaps=FALSE, ...)
## S3 method for class 'pca':
print(x, nmodes=6, \dots)

Arguments

xyz

numeric matrix of Cartesian coordinates with a row per structure.

subset

an optional vector of numeric indices that selects a subset of rows (e.g. experimental structures vs molecular dynamics trajectory structures) from the full xyz matrix. Note: the full xyz is projected onto this subspa

use.svd

logical, if TRUE singular value decomposition (SVD) is called instead of eigenvalue decomposition.

rm.gaps

logical, if TRUE gap positions (with missing coordinate data in any input structure) are removed before calculation. This is equivalent to removing NA cols from xyz.

x

an object of class pca, as obtained from function pca.xyz.

nmodes

numeric, number of modes to be printed.

...

additional arguments to print.

Value

• Returns a list with the following components:
• L
  eigenvalues.
• U
  eigenvectors (i.e. the x, y, and z variable loadings).
• z
  scores of the supplied xyz on the pcs.
• au
  atom-wise loadings (i.e. xyz normalized eigenvectors).
• sdev
  the standard deviations of the pcs.
• mean
  the means that were subtracted.

References

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

See Also

pca, pca.pdbs, plot.pca, mktrj.pca, pca.tor, project.pca

Examples

#-- Read transducin alignment and structures
aln <- read.fasta(system.file("examples/transducin.fa",package="bio3d"))
pdbs <- read.fasta.pdb(aln)

# Find core
core <- core.find(pdbs, 
                  #write.pdbs = TRUE,
                  verbose=TRUE)

rm(list=c("pdbs", "core"))

#-- OR for demo purposes just read previously saved transducin data
attach(transducin)

# Previously fitted coordinates based on sub 1.0A^3 core. See core.find() function.
xyz <- pdbs$xyz
                
#-- Do PCA ignoring gap containing positions
pc.xray <- pca.xyz(xyz, rm.gaps=TRUE)

# Plot results (conformer plots & scree plot)
plot(pc.xray, col=annotation[, "color"])

## Plot atom wise loadings
plot.bio3d(pc.xray$au[,1], ylab="PC1 (A)")

## Plot loadings in relation to reference structure 1TAG
gaps.res <- gap.inspect(pdbs$ali)
pdb <- read.pdb("1tag")

ind <- grep("1TAG", pdbs$id)
res.ind <- pdbs$resno[ind, gaps.res$f.ind]
op <- par(no.readonly=TRUE)
par(mfrow = c(3, 1), cex = 0.6, mar = c(3, 4, 1, 1))
plot.bio3d(res.ind, pc.xray$au[,1], sse=pdb, ylab="PC1 (A)")
plot.bio3d(res.ind, pc.xray$au[,2], sse=pdb, ylab="PC2 (A)")
plot.bio3d(res.ind, pc.xray$au[,3], sse=pdb, ylab="PC3 (A)")
par(op)

# Write PC trajectory
a <- mktrj.pca(pc.xray, pc=1, file="pc1.pdb",
               resno = pdbs$resno[1, gaps.res$f.inds],
               resid = aa123(pdbs$ali[1, gaps.res$f.inds]) )

b <- mktrj.pca(pc.xray, pc=2, file="pc2.pdb",
               resno = pdbs$resno[1, gaps.res$f.inds],
               resid = aa123(pdbs$ali[1, gaps.res$f.inds]) )

c <- mktrj.pca(pc.xray, pc=3, file="pc3.pdb",
               resno = pdbs$resno[1, gaps.res$f.inds],
               resid = aa123(pdbs$ali[1, gaps.res$f.inds]) )

detach(transducin)