# tests/testthat/test-SVD.R In bigstatsr: Statistical Tools for Filebacked Big Matrices

```################################################################################

context("SVD")

set.seed(SEED)

options(bigstatsr.downcast.warning = FALSE)

################################################################################

TOL <- 1e-3

# function for sampling scaling
sampleScale <- function() {
tmp <- sample(list(c(TRUE, FALSE),
c(TRUE, TRUE),
c(FALSE, FALSE)))[[1]]
list(center = tmp[1], scale = tmp[2])
}

# Simulating some data
N <- 73
M <- 43
x <- matrix(rnorm(N * M, mean = 100, sd = 5), N)

###############################################################################

svd1 <- big_SVD(as_FBM(x), big_scale(), ind.row = 1:40)
svd2 <- big_randomSVD(as_FBM(x), big_scale(), ind.row = 1:40, tol = 1e-5)
expect_equal(diffPCs(svd1\$u, svd2\$u), 0, tolerance = TOL)
expect_equal(diffPCs(svd1\$v, svd2\$v), 0, tolerance = TOL)

svd3 <- big_SVD(as_FBM(x), big_scale(), ind.row = 1:40 + 7, ind.col = 1:30 + 3)
svd4 <- big_randomSVD(as_FBM(x), big_scale(), ind.row = 1:40 + 7,
ind.col = 1:30 + 3, tol = 1e-8)
expect_equal(diffPCs(svd3\$u, svd4\$u), 0, tolerance = TOL)
expect_equal(diffPCs(svd3\$v, svd4\$v), 0, tolerance = TOL)

###############################################################################

test_that("equality with prcomp", {
for (t in TEST.TYPES) {
X <- `if`(t == "raw", asFBMcode(x), big_copy(x, type = t))

k <- sample(c(1, 2, 10, 20), 1)

test <- big_SVD(X, k = k)
pca <- prcomp(X[], center = FALSE, scale. = FALSE)
expect_equal(diffPCs(predict(test), pca\$x), 0, tolerance = TOL)
expect_equal(diffPCs(test\$v, pca\$rotation), 0, tolerance = TOL)

sc <- sampleScale()
test <- big_SVD(X,
fun.scaling = big_scale(center = sc\$center,
scale = sc\$scale),
k = k)
pca <- prcomp(X[], center = sc\$center, scale. = sc\$scale)
expect_equal(diffPCs(predict(test), pca\$x), 0, tolerance = TOL)
expect_equal(diffPCs(test\$v, pca\$rotation), 0, tolerance = TOL)
if (sc\$center) expect_equal(test\$center, pca\$center)
if (sc\$scale)  expect_equal(test\$scale,  pca\$scale)

expect_s3_class(p, "ggplot")
# expect_equal(p + theme_bigstatsr(1.2), MY_THEME(p, 1.2))

expect_error(predict(test, abc = 2), "Argument 'abc' not used.")
expect_error(plot(test, abc = 2), "Argument 'abc' not used.")
}
})

###############################################################################

test_that("equality with prcomp with half of the data", {
ind <- sample(N, N / 2)
ind2 <- setdiff(1:N, ind)

for (t in TEST.TYPES) {
X <- `if`(t == "raw", asFBMcode(x), big_copy(x, type = t))

k <- sample(c(1, 2, 10, 20), 1)
sc <- sampleScale()

test <- big_SVD(X,
ind.row = ind,
fun.scaling = big_scale(center = sc\$center,
scale = sc\$scale),
k = k)
pca <- prcomp(X[ind, ], center = sc\$center, scale. = sc\$scale)

expect_equal(diffPCs(predict(test), pca\$x), 0, tolerance = TOL)
expect_equal(diffPCs(test\$v, pca\$rotation), 0, tolerance = TOL)

if (sc\$center) expect_equal(test\$center, pca\$center)
if (sc\$scale)  expect_equal(test\$scale,  pca\$scale)

expect_equal(diffPCs(predict(test, X, ind.row = ind2),
predict(pca, X[ind2, ])), 0, tolerance = TOL)

expect_s3_class(p, "ggplot")
}
})

################################################################################

test_that("equality with prcomp with half of half of the data", {
ind <- sample(N, N / 2)
ind2 <- setdiff(1:N, ind)
ind.col <- sample(M, M / 2)

for (t in TEST.TYPES) {
X <- `if`(t == "raw", asFBMcode(x), big_copy(x, type = t))

k <- sample(c(1, 2, 10, 20), 1)
sc <- sampleScale()

test <- big_SVD(X,
ind.row = ind, ind.col = ind.col,
fun.scaling = big_scale(center = sc\$center,
scale = sc\$scale),
k = k)
pca <- prcomp(X[ind, ind.col], center = sc\$center, scale. = sc\$scale)

expect_equal(diffPCs(predict(test), pca\$x), 0, tolerance = TOL)
expect_equal(diffPCs(test\$v, pca\$rotation), 0, tolerance = TOL)

if (sc\$center) expect_equal(test\$center, pca\$center)
if (sc\$scale)  expect_equal(test\$scale,  pca\$scale)

expect_equal(diffPCs(predict(test, X, ind.row = ind2, ind.col = ind.col),
predict(pca, X[ind2, ind.col])), 0, tolerance = TOL)

expect_s3_class(p, "ggplot")
}
})

################################################################################

expect_s3_class(p <- plot(svd1, type = "scores"), "ggplot")
expect_length(p\$layers, 1)
expect_s3_class(p <- plot(svd1, type = "scores", scores = 3:4), "ggplot")
expect_length(p\$layers, 1)
expect_s3_class(p <- plot(svd1, type = "scores", scores = 1:8), "ggplot")
expect_length(p\$layers, 4)
expect_s3_class(p <- plot(svd1, type = "scores", scores = 1:9), "ggplot")
expect_length(p\$layers, 4)

###############################################################################
```

## Try the bigstatsr package in your browser

Any scripts or data that you put into this service are public.

bigstatsr documentation built on Oct. 14, 2022, 9:05 a.m.