Nothing
context("Testing aanma.pdbs()")
test_that("aanma based eNMA works", {
## commented lines:
## update after changing from ff=aaenm to ff=aaenm2
## new mass weighting approach (~may 2016)
## test updated 24-08-16
skip_on_cran()
skip_on_travis()
"mysign" <- function(a,b) {
if(all(sign(a)==sign(b)))
return(1)
else
return(-1)
}
pdbdir <- tempdir()
invisible( capture.output( capture.output(pdbfiles <- get.pdb(c("1TND_A", "1TAG", "1AS0", "1AS2"), path=pdbdir, split=TRUE), type="message") ))
invisible( capture.output(aln <- pdbaln(pdbfiles)) )
invisible( capture.output(pdbs <- read.all(aln)) )
gaps <- gap.inspect(pdbs$xyz)
## Calc modes (use default forcefield aanma2)
invisible(capture.output(modes <- aanma.pdbs(pdbs, pfc.fun=load.enmff("aaenm2"), ncore=1)))
## check dimensions
expect_that(dim(modes$U), equals(c(939, 933, 4)))
expect_that(dim(modes$L), equals(c(4, 933)))
expect_that(dim(modes$fluctuations), equals(c(4, 313)))
expect_that(dim(modes$rmsip), equals(c(4, 4)))
## structure 1- mode1:
## U1 <- c(0.04766153, -0.007993432, -0.06513697, 0.04313399, -0.008026397, -0.03755748)
U1 <- c(0.04785697, -0.007613113, -0.063552590, 0.043594060, -0.007665827, -0.037405791)
nowU1 <- head(modes$U.subspace[,1,1], n=6)
expect_that(nowU1 * mysign(U1, nowU1), equals(U1, tolerance=1e-4))
## structure 1- mode2:
## U2 <- c(0.01151413, -0.01098311, -0.05809988, 0.006112359, -0.0164621, -0.04298916)
U2 <- c(-0.01051829, 0.01169047, 0.05700469, -0.00545407, 0.01683833, 0.04307286)
nowU2 <- head(modes$U.subspace[,2,1], n=6)
expect_that(nowU2 * mysign(U2, nowU2), equals(U2, tolerance=1e-4))
## structure 4- mode3:
##U3 <- c(-0.1188249, -0.05850153, 0.02082409, -0.09123394, -0.02270938, 0.00987617)
U3 <- c(0.113108253, 0.057547242, -0.013566069, 0.088356656, 0.025218140, -0.006994051)
nowU3 <- head(modes$U.subspace[,3,4], n=6)
expect_that(nowU3 * mysign(U3, nowU3), equals(U3, tolerance=1e-4))
## structure 4-mode1 - tail:
## U1 <- c(0.009886232, -0.0006358995, -0.05189479, 0.009095127, 0.0005291993, -0.06936025)
U1 <- c(1.008168e-02, -7.579627e-05, -5.113142e-02, 8.872045e-03, 1.126475e-03, -6.340656e-02)
nowU1 <- tail(modes$U.subspace[,1,4], n=6)
expect_that(nowU1 * mysign(U1, nowU1), equals(U1, tolerance=1e-4))
## Fluctuations:
## f1 <- c(0.4259852, 0.3173786, 0.2339003, 0.2126361, 0.1965271, 0.1865991)
## f2 <- c(0.5467513, 0.4003687, 0.2712613, 0.2266669, 0.2006127, 0.1855429)
## f4 <- c(0.8295543, 0.4094441, 0.281545, 0.2527504, 0.2091488, 0.1931585)
f1 <- c(0.2937191, 0.2410464, 0.1696626, 0.1553645, 0.1328955, 0.1237624)
f2 <- c(0.4076411, 0.3114183, 0.2020671, 0.1662545, 0.1409823, 0.1274508)
f4 <- c(0.6278258, 0.3151669, 0.2104897, 0.1819106, 0.1426998, 0.1333115)
expect_that(modes$fluctuations[1,1:6], equals(f1, tolerance=1e-3))
expect_that(modes$fluctuations[2,1:6], equals(f2, tolerance=1e-3))
expect_that(modes$fluctuations[4,1:6], equals(f4, tolerance=1e-3))
## Orthognal
expect_that(as.numeric(modes$U.subspace[,1,1] %*% modes$U.subspace[,1,1]),
equals(1, tolerance=1e-6))
expect_that(as.numeric(modes$U.subspace[,1,1] %*% modes$U.subspace[,2,1]),
equals(0, tolerance=1e-6))
## RMSIP
## rmsips <- c(1.0000, 0.8988234, 0.9627764, 0.915605)
rmsips <- c(1.0000, 0.8927, 0.9600, 0.9158)
expect_that(as.vector(modes$rmsip[1,]), equals(rmsips, tolerance=1e-3))
## Multicore (same arguments as above!)
invisible(capture.output(mmc <- aanma.pdbs(pdbs, pfc.fun=load.enmff("aaenm2"), ncore=NULL)))
expect_that(mmc$fluctuations, equals(modes$fluctuations, tolerance=1e-6))
expect_that(mmc$U.subspace, equals(modes$U.subspace, tolerance=1e-6))
## Calc modes with rm.gaps=FALSE
invisible(capture.output(modes <- aanma.pdbs(pdbs, pfc.fun=load.enmff("aaenm2"), rm.gaps=FALSE, ncore=NULL)))
## structure 1-mode1 - tail:
## U1 <- c(-0.05024617, 0.009081938, 0.01764626, -0.09325497, 0.01754876, 0.01996309)
U1 <- c(0.051308400, -0.008862397, -0.017957255, 0.085623903, -0.015406574, -0.018774554)
nowU1 <- tail(modes$U.subspace[,1,1], n=6)
expect_that(nowU1 * mysign(U1, nowU1), equals(U1, tolerance=1e-4))
## structure 2-mode1 - tail:
## U1 <- c(0.008412425, -0.005413275, -0.06655021, NA, NA, NA)
U1 <- c(0.007581776, -0.005342612, -0.060598339, NA, NA, NA)
nowU1 <- modes$U.subspace[940:945,1,2]
U1[is.na(U1)] <- 0
nowU1[is.na(nowU1)] <- 0
expect_that(nowU1 * mysign(nowU1, U1), equals(U1, tolerance=1e-4))
## fluctuations
na.expected <- c(3, 4, 1258, 1259, 1262, 1263, 1266, 1267, 1268, 1270, 1271, 1272, 1274, 1275, 1276,
1278, 1279, 1280, 1282, 1283, 1284, 1286, 1287, 1288, 1290, 1291, 1292)
expect_that(which(is.na(modes$fluctuations)), equals(na.expected))
## f1 <- c(0.9114199, 0.4010512, 0.2977837, 0.2202269)
## f4 <- c(0.4401748, 0.5743443, 0.6644086, rep(NA, 7))
f1 <- c(0.6954173, 0.2722280, 0.2242001, 0.1581307)
f4 <- c(0.3369227, 0.4506579, 0.5365029, rep(NA, 7))
expect_that(modes$fluctuations[1,1:4], equals(f1, tolerance=1e-3))
expect_that(tail(modes$fluctuations[4,], n=10), equals(f4, tolerance=1e-3))
## Calc modes with mass=FALSE and temp=NULL -- use aaenm (not default ff)
invisible(capture.output(modes <- aanma.pdbs(pdbs, pfc.fun=load.enmff("aaenm"), mass=FALSE, temp=NULL, ncore=NULL)))
## structure 1- mode1:
U1 <- c(0.04185427, -0.007185472, -0.05865138, 0.04704491, -0.009310229, -0.0432674)
nowU1 <- head(modes$U.subspace[,1,1], n=6)
expect_that(nowU1 * mysign(U1, nowU1), equals(U1, tolerance=1e-4))
## structure 1- mode2:
U2 <- c(-0.006541777, 0.01006403, 0.04607164, -0.00293836, 0.01800689, 0.043626)
nowU2 <- head(modes$U.subspace[,2,1], n=6)
expect_that(nowU2 * mysign(U2, nowU2), equals(U2, tolerance=1e-4))
## structure 4- mode3:
U3 <- c(-0.09327228, -0.04605804, 0.01065858, -0.08065877, -0.02255212, 0.005277435)
nowU3 <- head(modes$U.subspace[,3,4], n=6)
expect_that(nowU3 * mysign(U3, nowU3), equals(U3, tolerance=1e-4))
})
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.