Convenience wrapper for eigenanatomy decomposition.

Decomposes a matrix into sparse eigenevectors to maximize explained variance. Note: we do not scale the matrices internally. We leave scaling choices to the user.

sparseDecom(
  inmatrix = NA,
  inmask = NULL,
  sparseness = 0.1,
  nvecs = 10,
  its = 5,
  cthresh = 50,
  statdir = NA,
  z = 0,
  smooth = 0,
  initializationList = list(),
  mycoption = 0,
  robust = 0,
  ell1 = 1,
  getSmall = 0,
  verbose = FALSE,
  powerit = 0,
  priorWeight = 0,
  maxBased = FALSE
)

Arguments

inmatrix

n by p input images , subjects or time points by row , spatial variable lies along columns

inmask

optional antsImage mask

sparseness

lower values equal more sparse

nvecs

number of vectors

its

number of iterations

cthresh

cluster threshold

statdir

place on disk to save results

z

u penalty, experimental

smooth

smoothness eg 0.5

initializationList

see initializeEigenanatomy

mycoption

0, 1 or 2 all produce different output 0 is combination of 1 (spatial orthogonality) and 2 (subject space orthogonality)

robust

rank transform input data - good for data checking

ell1

the ell1 grad descent param

getSmall

try to get smallest evecs (bool)

verbose

activates verbose output

powerit

alternative power iteration implementation, faster

priorWeight

scalar weight typically in range zero to two

maxBased

boolean that chooses max-based thresholding

Value

outputs a decomposition of a population or time series matrix

Author

Avants BB

Examples

mat <- replicate(100, rnorm(20))
mydecom <- sparseDecom(mat)
mat <- scale(mat)
mydecom2 <- sparseDecom(mat)
# params that lead to algorithm similar to NMF
mydecom3 <- sparseDecom(mat, z = 1, sparseness = 1)

if (FALSE) { # \dontrun{
# for prediction
if (usePkg("randomForest") & usePkg("spls") & usePkg("BGLR")) {
  data(lymphoma) # from spls
  training <- sample(rep(c(TRUE, FALSE), 31))
  sp <- 0.02
  myz <- 0
  ldd <- sparseDecom(lymphoma$x[training, ],
    nvecs = 5, sparseness = (sp),
    mycoption = 1, z = myz
  ) # NMF style
  traindf <- data.frame(
    lclass = as.factor(lymphoma$y[training]),
    eig = lymphoma$x[training, ] %*% as.matrix(ldd$eigenanatomyimages)
  )
  testdf <- data.frame(
    lclass = as.factor(lymphoma$y[!training]),
    eig = lymphoma$x[!training, ] %*% as.matrix(ldd$eigenanatomyimages)
  )
  myrf <- randomForest(lclass ~ ., data = traindf)
  predlymp <- predict(myrf, newdata = testdf)
  print(paste(
    "N-errors:", sum(abs(testdf$lclass != predlymp)),
    " non-zero ", sum(abs(ldd$eigenanatomyimages) > 0)
  ))
  # compare to http://arxiv.org/pdf/0707.0701v2.pdf
  # now SNPs
  data(mice)
  snps <- quantifySNPs(mice.X, shiftit = TRUE)
  numericalpheno <- as.matrix(mice.pheno[, c(4, 5, 13, 15)])
  nfolds <- 6
  train <- sample(rep(c(1:nfolds), 1800 / nfolds))
  train <- (train < 4)
  lrmat <- lowrankRowMatrix(as.matrix(snps[train, ]), 50)
  lrmat <- scale(lrmat)
  snpd <- sparseDecom(lrmat - min(lrmat), nvecs = 20, sparseness = (0.001), z = -1)
  projmat <- as.matrix(snpd$eig)
  snpse <- as.matrix(snps[train, ]) %*% projmat
  traindf <- data.frame(bmi = numericalpheno[train, 3], snpse = snpse)
  snpse <- as.matrix(snps[!train, ]) %*% projmat
  testdf <- data.frame(bmi = numericalpheno[!train, 3], snpse = snpse)
  myrf <- randomForest(bmi ~ ., data = traindf)
  preddf <- predict(myrf, newdata = testdf)
  cor.test(preddf, testdf$bmi)
  plot(preddf, testdf$bmi)
} # check for packages
# prior-based example
set.seed(123)
ref <- antsImageRead(getANTsRData("r16"))
ref <- iMath(ref, "Normalize")
mi <- antsImageRead(getANTsRData("r27"))
mi2 <- antsImageRead(getANTsRData("r30"))
mi3 <- antsImageRead(getANTsRData("r62"))
mi4 <- antsImageRead(getANTsRData("r64"))
mi5 <- antsImageRead(getANTsRData("r85"))
refmask <- getMask(ref)
refmask <- iMath(refmask, "ME", 2) # just to speed things up
ilist <- list(mi, mi2, mi3, mi4, mi5)
for (i in 1:length(ilist))
{
  ilist[[i]] <- iMath(ilist[[i]], "Normalize")
  mytx <- antsRegistration(
    fixed = ref, moving = ilist[[i]],
    typeofTransform = c("Affine")
  )
  mywarpedimage <- antsApplyTransforms(
    fixed = ref, moving = ilist[[i]],
    transformlist = mytx$fwdtransforms
  )
  ilist[[i]] <- mywarpedimage
}
mat <- imageListToMatrix(ilist, refmask)
kmseg <- kmeansSegmentation(ref, 3, refmask)
initlist <- list()
for (k in 1:3) {
  initlist[[k]] <-
    thresholdImage(kmseg$probabilityimages[[k]], 0.1, Inf) *
      kmseg$probabilityimages[[k]]
}
eanat <- sparseDecom(mat,
  inmask = refmask, ell1 = 0.1,
  sparseness = 0.0, smooth = 0.5, verbose = 1,
  initializationList = initlist, cthresh = 25,
  nvecs = 3, priorWeight = 0.5
)
ee <- matrixToImages(eanat$eigenanatomyimages, refmask)
eseg <- eigSeg(refmask, ee)
priormat <- imageListToMatrix(initlist, refmask)
cor(t(eanat$eigenanatomyimages), t(priormat))
plot(ref, eseg)
} # }