Why does RMSIP value show greater than 1?
I am running all aaenma
and using gap.inspect
option to compare modes of proteins with varying length. When I calculate RMSIP value, why am I getting value greater than 1. What could be the reason? Is it because of vector length mismatch.
Comments (7)
-
-
reporter I am using the following code.
library(bio3d) a <- pdbaln(c("4YUG.pdb","4YUG_s1_bp.pdb","4YUG_1_tf.pdb")) b <- read.all(a) gaps.pos <- gap.inspect(b$all) modes <- aanma.pdbs(pdbs = b, outmodes='noh') modes1 <- matrix(NA, nrow=ncol(b$all), ncol=ncol(modes[[1]]$U)) modes2 <- matrix(NA, nrow=ncol(b$all), ncol=ncol(modes[[2]]$U)) modes3 <- matrix(NA, nrow=ncol(b$all), ncol=ncol(modes[[3]]$U)) modes1[!is.gap(b$all[1, ]), ] <- modes[[1]]$U modes2[!is.gap(b$all[2, ]), ] <- modes[[2]]$U modes3[!is.gap(b$all[3, ]), ] <- modes[[3]]$U rm1 <- rmsip(modes1[gaps.pos$f.inds, ], modes2[gaps.pos$f.inds, ], subset=10) rm2 <- rmsip(modes1[gaps.pos$f.inds, ], modes3[gaps.pos$f.inds, ], subset=10) print(rm2$rmsip)
I have used
gaps.pos$f.inds
to avoid dimension mismatch error in rmsip calculation for 4YUG and 4YUG_s1_bp protein pair. Whereas there is no gaps in the alignment between 4YUG and 4YUG_1_tf pairs. When I calculate rmsip for 4YUG and 4YUG_1_tf pair, it gives value greater than 1. Since I am new to this work, it is totally perplexing. -
reporter - attached 4YUG_1_tf.pdb
- attached 4YUG.pdb
- attached 4YUG_s1_bp.pdb
-
The modes are no longer orthogonal when subsetting like this. You need to calculate the effective hessian and diagonalize it. This is done internally in aanma.pdbs() (and nma.pdbs()) when
outmodes = "calpha"
, but not whenoutmodes = "noh"
as aligning heavy atoms is not implemented. -
reporter Thank you very much for the clarification. I can now understand my mistake.
-
Thanks for explaining, Lars. What about normalized the modes? For example
normalize.vector(modes1[gaps.pos$f.inds, ])
, and then using normalized modes to dormsip()
should give a value below 1. -
- changed status to resolved
- Log in to comment
Hi,
Can you tell what exactly the commands you used?