Mercurial > repos > marie-tremblay-metatoul > nmr_annotation
view nmr_preprocessing/ptw.R @ 2:7304ec2c9ab7 draft
Uploaded
author | marie-tremblay-metatoul |
---|---|
date | Mon, 30 Jul 2018 10:33:03 -0400 |
parents | |
children |
line wrap: on
line source
ptw <- function (ref, samp, selected.traces, init.coef = c(0, 1, 0), try = FALSE, warp.type = c("individual", "global"), optim.crit = c("WCC", "RMS"), mode = c("forward", "backward"), smooth.param = ifelse(try, 0, 1e+05), trwdth = 20, trwdth.res = trwdth, verbose = FALSE, ...) { optim.crit <- match.arg(optim.crit) warp.type <- match.arg(warp.type) mode <- match.arg(mode) if (is.vector(ref)) ref <- matrix(ref, nrow = 1) if (is.vector(samp)) samp <- matrix(samp, nrow = 1) if (nrow(ref) > 1 && nrow(ref) != nrow(samp)) stop("The number of references does not equal the number of samples") if (length(dim(ref)) > 2) stop("Reference cannot be an array") if (length(dim(samp)) > 2) stop("Sample cannot be an array") if (nrow(samp) == 1) warp.type <- "individual" r <- nrow(samp) if (!missing(selected.traces)) { samp <- samp[selected.traces, , drop = FALSE] if (nrow(ref) > 1) ref <- ref[selected.traces, , drop = FALSE] } if (is.vector(init.coef)) init.coef <- matrix(init.coef, nrow = 1) if (warp.type == "global") { if (nrow(init.coef) != 1) stop("Only one warping function is allowed with global alignment.") } else { if (nrow(init.coef) != nrow(samp)) if (nrow(init.coef) == 1) { init.coef <- matrix(init.coef, byrow = TRUE, nrow = nrow(samp), ncol = length(init.coef)) } else { stop("The number of warping functions does not match the number of samples") } } if (warp.type == "individual") { w <- matrix(0, nrow(samp), ncol(ref)) a <- matrix(0, nrow(samp), ncol(init.coef)) v <- rep(0, nrow(samp)) warped.sample <- matrix(NA, nrow = nrow(samp), ncol = ncol(samp)) for (i in 1:nrow(samp)) { if (verbose & nrow(samp) > 1) cat(ifelse(nrow(ref) == 1, paste("Warping sample", i, "with the reference \n"), paste("Warping sample", i, "with reference \n", i))) if (nrow(ref) == 1) { rfrnc <- ref } else { rfrnc <- ref[i, , drop = FALSE] } quad.res <- pmwarp(rfrnc, samp[i, , drop = FALSE], optim.crit, init.coef[i, ], try = try, mode = mode, smooth.param = smooth.param, trwdth = trwdth, trwdth.res = trwdth.res, ...) w[i, ] <- quad.res$w a[i, ] <- quad.res$a v[i] <- quad.res$v warped.sample[i, ] <- c(warp.sample(samp[i, , drop = FALSE], w[i, ], mode = mode)) } } else { if (nrow(ref) == 1) ref <- matrix(ref, nrow = nrow(samp), ncol = ncol(ref), byrow = TRUE) if (verbose) { if (nrow(ref) == 1) { cat("Simultaneous warping of samples with reference... \n") } else { cat("Simultaneous warping of samples with references... \n") } } quad.res <- pmwarp(ref, samp, optim.crit, c(init.coef), try = try, mode = mode, smooth.param = smooth.param, trwdth = trwdth, trwdth.res = trwdth.res, ...) w <- t(as.matrix(quad.res$w)) a <- t(as.matrix(quad.res$a)) v <- quad.res$v warped.sample <- t(warp.sample(samp, w, mode)) } if (verbose) cat("\nFinished.\n") result <- list(reference = ref, sample = samp, warped.sample = warped.sample, warp.coef = a, warp.fun = w, crit.value = v, optim.crit = optim.crit, mode = mode, warp.type = warp.type) class(result) <- "ptw" result }