n1 <- 10
n2 <- 10
p <- 6
y <- rnorm(n1)
a <- round(rnorm(n1*p), 2)
a[abs(a) < 0.5] <- 0
A <- matrix(a, n1,p)
A.csr <- as.matrix.csr(A)
B <- matrix(c(1.5, 0, 0, 0, -1.4, 0, 0, 0, 0, -1.4,
2, 0, -1, 0, 0, 2.1, -1.9, 1.4, 0, 0,
0,-2.3, 0, 0, -1.9, 0, 0, 0, 0, -1.4,
0, 0, 0, 0, 0, -3, 0, 1.3, 0, 1.1,
0, 0, 0, 0, 2, 0, 0, 0, -1, 0,
0, 0, -1.6,0, 0, 0, 0, 0, -1.7,0),
10L, 6L)
rcond(B) # 0.21 .. i.e., quite well conditioned
B.csr <- as.matrix.csr(B)
B.csr
## norm() : different 'type' for base R and SparseM:
(nR <- vapply(c("1", "I", "F", "M", "2"), norm, 0, x = B))
## 1 I F M 2
## 8.400000 5.300000 7.372923 3.000000 4.464801
(nSpM <- vapply(c("sup","M", "HS","F", "l1"), norm, 0, x = B.csr))
## sup M HS F l1
## 3.000000 3.000000 7.372923 7.372923 30.000000
stopifnot(all.equal(unname(nSpM[c("M", "F")]),
unname(nR [c("M", "F" )]), tolerance = 1e-14))
# matrix transposition and multiplication
BtB <- crossprod(B) # == t(B) %*% B {more efficiently}
A.csr %*% t(B.csr)
BtBs <- t(B.csr) %*% B.csr
BtBs
stopifnot(all.equal( BtB, as.matrix(BtBs), tolerance = 1e-14),
all.equal(det(BtB), print(det(BtBs)), tolerance = 1e-14))
# matrix o vector
stopifnot(all.equal(y %*% A , y %*% A.csr) ,
all.equal(A %*% 1:6, A.csr %*% 1:6)
)
# kronecker product - via kronecker() methods:
A.csr %x% matrix(1:4,2,2)
Run the code above in your browser using DataLab