View source: R/multiscaleSVDxpts.R
jointSmoothMatrixReconstruction | R Documentation |
Joints reconstruct a n by p, q, etc matrices or predictors. The reconstruction can be regularized. # norm( x - u_x v_x^t ) + norm( y - u_y v_y^t) + .... etc
jointSmoothMatrixReconstruction(
x,
nvecs,
parameters,
iterations = 10,
subIterations = 5,
gamma = 1e-06,
sparsenessQuantile = 0.5,
positivity = FALSE,
smoothingMatrix = NA,
rowWeights = NA,
repeatedMeasures = NA,
doOrth = FALSE,
verbose = FALSE
)
x |
input list of matrices to be jointly predicted. |
nvecs |
number of basis vectors to compute |
parameters |
should be a ncomparisons by 3 matrix where the first two columns define the pair to be matched and the last column defines the weight in the objective function. |
iterations |
number of gradient descent iterations |
subIterations |
number of gradient descent iterations in sub-algorithm |
gamma |
step size for gradient descent |
sparsenessQuantile |
quantile to control sparseness - higher is sparser |
positivity |
restrict to positive or negative solution (beta) weights. choices are positive, negative or either as expressed as a string. |
smoothingMatrix |
a list containing smoothing matrices of the same length as x. |
rowWeights |
vectors of weights with size n (assumes diagonal covariance) |
repeatedMeasures |
list of repeated measurement identifiers. this will allow estimates of per identifier intercept. |
doOrth |
boolean enforce gram-schmidt orthonormality |
verbose |
boolean option |
matrix list each of size p by k is output
Avants BB
## Not run:
mat <- replicate(100, rnorm(20))
mat2 <- replicate(100, rnorm(20))
mat <- scale(mat)
mat2 <- scale(mat2)
wt <- 0.666
mat3 <- mat * wt + mat2 * (1 - wt)
params <- matrix(nrow = 2, ncol = 3)
params[1, ] <- c(1, 2, 1)
params[2, ] <- c(2, 1, 1)
x <- list((mat), (mat3))
jj <- jointSmoothMatrixReconstruction(x, 2, params,
gamma = 1e-4, sparsenessQuantile = 0.5, iterations = 10,
smoothingMatrix = list(NA, NA), verbose = TRUE
)
# latent feature
mask <- getMask(antsImageRead(getANTsRData("r16")))
spatmat <- t(imageDomainToSpatialMatrix(mask, mask))
smoomat <- knnSmoothingMatrix(spatmat, k = 27, sigma = 120.0)
lfeats <- t(replicate(100, rnorm(3)))
# map these - via matrix - to observed features
n <- sum(mask)
ofeats1 <- (lfeats + rnorm(length(lfeats), 0.0, 1.0)) %*% rbind(rnorm(n), rnorm(n), rnorm(n))
ofeats2 <- (lfeats + rnorm(length(lfeats), 0.0, 1.0)) %*% rbind(rnorm(n), rnorm(n), rnorm(n))
# only half of the matrix contains relevant data
ofeats1[, 1:round(n / 2)] <- matrix(rnorm(round(n / 2) * 100), nrow = 100)
ofeats2[, 1:round(n / 2)] <- matrix(rnorm(round(n / 2) * 100), nrow = 100)
x <- list((ofeats1), (ofeats2))
jj <- jointSmoothMatrixReconstruction(x, 2, params,
gamma = 0.0001, sparsenessQuantile = 0.75, iterations = 19,
subIterations = 11,
smoothingMatrix = list(smoomat, smoomat), verbose = TRUE
)
p1 <- ofeats2 %*% jj$v[[2]]
p2 <- ofeats1 %*% jj$v[[1]]
cor(p1, lfeats)
cor(p2, lfeats)
print(cor(rowMeans(ofeats1), lfeats))
print(cor(rowMeans(ofeats2), lfeats))
# a 2nd example with 3 modalities
imageIDs <- c("r16", "r27", "r30", "r62", "r64", "r85")
images <- list()
feature1Images <- list()
feature2Images <- list()
feature3Images <- list()
ref <- antsImageRead(getANTsRData("r16"))
for (i in 1:length(imageIDs))
{
cat("Processing image", imageIDs[i], "\n")
tar <- antsRegistration(ref, antsImageRead(getANTsRData(imageIDs[i])),
typeofTransform = "Affine"
)$warpedmov
images[[i]] <- tar
feature1Images[[i]] <- iMath(images[[i]], "Grad", 1.0)
feature2Images[[i]] <- iMath(images[[i]], "Laplacian", 1.0)
feature3Images[[i]] <- reflectImage(tar, axis = 0, tx = "Affine")$warpedmovout
}
i <- 1
mask <- getMask(antsImageRead(getANTsRData(imageIDs[i])))
mask2 <- iMath(mask, "ME", 2)
spatmat <- t(imageDomainToSpatialMatrix(mask, mask))
smoomat <- knnSmoothingMatrix(spatmat, k = 125, sigma = 100.0)
spatmat2 <- t(imageDomainToSpatialMatrix(mask2, mask2))
smoomat2 <- knnSmoothingMatrix(spatmat2, k = 125, sigma = 100.0)
params <- matrix(nrow = 6, ncol = 3)
params[1, ] <- c(1, 2, 1)
params[2, ] <- c(2, 1, 1)
params[3, ] <- c(1, 3, 1)
params[4, ] <- c(3, 1, 1)
params[5, ] <- c(2, 3, 1)
params[6, ] <- c(3, 2, 1)
mat <- imageListToMatrix(feature1Images, mask)
mat2 <- imageListToMatrix(feature2Images, mask2)
mat3 <- imageListToMatrix(feature3Images, mask)
scl <- F
x <- list(scale(mat, scale = scl), scale(mat2, scale = scl), scale(mat3, scale = scl))
slist <- list(smoomat2, smoomat, smoomat, smoomat, smoomat, smoomat2)
jj <- jointSmoothMatrixReconstruction(x, 4, params,
positivity = T,
gamma = 1e-6, sparsenessQuantile = 0.9, iterations = 10,
smoothingMatrix = slist, verbose = TRUE
)
mm <- makeImage(mask, abs(jj$v[[2]][, 1])) %>% iMath("Normalize")
plot(ref, mm, doCropping = FALSE, window.overlay = c(0.1, 1))
p1 <- mat2 %*% jj$v[[1]]
p2 <- mat %*% jj$v[[2]]
diag(cor(p1, p2))
p1 <- mat3 %*% jj$v[[5]]
p2 <- mat2 %*% jj$v[[6]]
diag(cor(p1, p2))
## End(Not run)
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.