Skip to content

Commit

Permalink
ivar = NA will now skip log-ratio transform
Browse files Browse the repository at this point in the history
  • Loading branch information
tpq committed Apr 5, 2020
1 parent 614e816 commit b85112e
Show file tree
Hide file tree
Showing 16 changed files with 151 additions and 61 deletions.
4 changes: 2 additions & 2 deletions DESCRIPTION
Original file line number Diff line number Diff line change
@@ -1,6 +1,6 @@
Package: propr
Title: Calculating Proportionality Between Vectors of Compositional Data
Version: 4.2.6
Version: 4.2.7
URL: http://github.com/tpq/propr
BugReports: http://github.com/tpq/propr/issues
Authors@R: c(
Expand All @@ -22,7 +22,7 @@ Description: The bioinformatic evaluation of gene co-expression often begins wit
License: GPL-2
LazyData: true
VignetteBuilder: knitr
RoxygenNote: 6.1.1
RoxygenNote: 7.0.2
Encoding: UTF-8
Depends:
methods,
Expand Down
4 changes: 4 additions & 0 deletions NEWS.md
Original file line number Diff line number Diff line change
@@ -1,3 +1,7 @@
## propr 4.2.7
---------------------
* Allow user to skip log-ratio transformation with `ivar = NA`

## propr 4.2.6
---------------------
* Fix new R version `class` bug where matrix is also an array
Expand Down
3 changes: 2 additions & 1 deletion R/0-params.R
Original file line number Diff line number Diff line change
Expand Up @@ -12,7 +12,8 @@
#' for additive log-ratio transformation. The argument will also
#' accept feature name(s) instead of the index position(s).
#' Set to "iqlr" to use inter-quartile log-ratio transformation.
#' Ignore to use centered log-ratio transformation.
#' Ignore to use centered log-ratio transformation. If the data
#' are already transformed, skip it with \code{ivar = NA}.
#' @param select Optional. Use this to subset the final
#' proportionality matrix without altering the result.
#' Use this argument to rearrange feature order.
Expand Down
81 changes: 46 additions & 35 deletions R/1-propr.R
Original file line number Diff line number Diff line change
Expand Up @@ -128,49 +128,59 @@ propr <- function(counts, metric = c("rho", "phi", "phs", "cor", "vlr"), ivar =
select, symmetrize = FALSE, alpha, p = 100){

# Clean "count matrix"
# if(any(apply(counts, 2, function(x) all(x == 0)))){
# stop("Remove components with all zeros before proceeding.")}
if(any(counts < 0)) stop("Data may not contain negative measurements.")
if(any(is.na(counts))) stop("Remove NAs from 'counts' before proceeding.")
if("data.frame" %in% class(counts)) counts <- as.matrix(counts)
if(is.null(colnames(counts))) colnames(counts) <- as.character(1:ncol(counts))
if(is.null(rownames(counts))) rownames(counts) <- as.character(1:nrow(counts))
if(missing(alpha)){ alpha <- NA
}else{ if(!is.na(alpha)) if(alpha == 0) alpha <- NA }
ct <- counts

# Replace zeros unless alpha is provided
if(any(as.matrix(counts) == 0) & is.na(alpha)){
message("Alert: Replacing 0s with next smallest value.")
zeros <- ct == 0
ct[zeros] <- min(ct[!zeros])
}
if(any(is.na(counts))) stop("Remove NAs from 'counts' before proceeding.")

# Establish reference
use <- ivar2index(ct, ivar)
if(!is.na(ivar)){ # if a log-ratio transform is performed...

# Transform counts based on reference
if(is.na(alpha)){
if(any(counts < 0)) stop("Data may not contain negative measurements.")

# Use g(x) = Mean[log(x)] to log-ratio transform data
message("Alert: Saving log-ratio transformed counts to @logratio.")
logX <- log(ct)
logSet <- logX[, use, drop = FALSE]
ref <- rowMeans(logSet)
lr <- sweep(logX, 1, ref, "-")
if(missing(alpha)){ alpha <- NA
}else{ if(!is.na(alpha)) if(alpha == 0) alpha <- NA }
ct <- counts

}else{
# Replace zeros unless alpha is provided
if(any(as.matrix(counts) == 0) & is.na(alpha)){
message("Alert: Replacing 0s with next smallest value.")
zeros <- ct == 0
ct[zeros] <- min(ct[!zeros])
}

# Establish reference
use <- ivar2index(ct, ivar)

# Transform counts based on reference
if(is.na(alpha)){

# Since y = log(x) = [x^a-1]/a, ref = Mean[log(x)] = Mean[y]
# Calculate log(x/ref) = log(x) - log(ref) = [x^a-1]/a - [ref^a-1]/a
message("Alert: Saving alpha transformed counts to @logratio.")
aX <- (ct^alpha - 1) / alpha
aSet <- aX[, use, drop = FALSE]
ref <- rowMeans(aSet)
lr <- sweep(aX, 1, ref, "-")
# Use g(x) = Mean[log(x)] to log-ratio transform data
message("Alert: Saving log-ratio transformed counts to @logratio.")
logX <- log(ct)
logSet <- logX[, use, drop = FALSE]
ref <- rowMeans(logSet)
lr <- sweep(logX, 1, ref, "-")

}else{

# Since y = log(x) = [x^a-1]/a, ref = Mean[log(x)] = Mean[y]
# Calculate log(x/ref) = log(x) - log(ref) = [x^a-1]/a - [ref^a-1]/a
message("Alert: Saving alpha transformed counts to @logratio.")
aX <- (ct^alpha - 1) / alpha
aSet <- aX[, use, drop = FALSE]
ref <- rowMeans(aSet)
lr <- sweep(aX, 1, ref, "-")
}

}else{ # if skipping the log-ratio transform...

message("Alert: Skipping the log-ratio transformation.")
alpha <- NA
ct <- counts
lr <- ct
}

# Optionally apply select for performance gain
# Optionally reduce the number of features used in matrix algebra
if(!missing(select)){

message("Alert: Using 'select' disables permutation testing.")
Expand Down Expand Up @@ -222,7 +232,8 @@ propr <- function(counts, metric = c("rho", "phi", "phs", "cor", "vlr"), ivar =
# Set up @permutes
result@permutes <- list(NULL)
if(p > 0){
# Sample compositions so that the clr does not change

# Shuffle row-wise so that per-sample CLR does not change
message("Alert: Fixing permutations to active random seed.")
permutes <- vector("list", p)
for(ins in 1:length(permutes)) permutes[[ins]] <- t(apply(ct, 1, sample))
Expand All @@ -242,7 +253,7 @@ propr <- function(counts, metric = c("rho", "phi", "phs", "cor", "vlr"), ivar =
)

# Initialize @results -- Tally frequency of 0 counts
if(any(as.matrix(counts) == 0)){
if(any(as.matrix(counts) == 0) & !is.na(ivar)){
message("Alert: Tabulating the presence of 0 counts.")
result@results$Zeros <- ctzRcpp(as.matrix(counts)) # count 0s
}
Expand Down
1 change: 0 additions & 1 deletion R/6-propdViz.R
Original file line number Diff line number Diff line change
Expand Up @@ -49,7 +49,6 @@ pals <- function(object, k){
#' where the \code{@@matrix} slot contains \eqn{1 - \theta}.
#' Allows the user to interrogate theta using any
#' visualization built for \code{propr} objects.
#' @rdname propd
#' @export
propd2propr <- function(object, ivar){

Expand Down
3 changes: 2 additions & 1 deletion man/all.Rd

Some generated files are not rendered by default. Learn more about how customized files appear on GitHub.

11 changes: 9 additions & 2 deletions man/calculateTheta.Rd

Some generated files are not rendered by default. Learn more about how customized files appear on GitHub.

18 changes: 15 additions & 3 deletions man/getNetwork.Rd

Some generated files are not rendered by default. Learn more about how customized files appear on GitHub.

3 changes: 1 addition & 2 deletions man/getRatios.Rd

Some generated files are not rendered by default. Learn more about how customized files appear on GitHub.

3 changes: 2 additions & 1 deletion man/ivar2index.Rd

Some generated files are not rendered by default. Learn more about how customized files appear on GitHub.

9 changes: 7 additions & 2 deletions man/pra.Rd

Some generated files are not rendered by default. Learn more about how customized files appear on GitHub.

3 changes: 2 additions & 1 deletion man/propd.Rd

Some generated files are not rendered by default. Learn more about how customized files appear on GitHub.

14 changes: 11 additions & 3 deletions man/propr.Rd

Some generated files are not rendered by default. Learn more about how customized files appear on GitHub.

3 changes: 2 additions & 1 deletion man/proprPerb.Rd

Some generated files are not rendered by default. Learn more about how customized files appear on GitHub.

20 changes: 14 additions & 6 deletions man/visualize.Rd

Some generated files are not rendered by default. Learn more about how customized files appear on GitHub.

32 changes: 32 additions & 0 deletions tests/testthat/test-multiomic.R
Original file line number Diff line number Diff line change
@@ -0,0 +1,32 @@
library(propr)
data(iris)
met.rel <- iris[,1:2]
mic.rel <- iris[,3:4]

# Analyze multi-omics via back-end
clr <- function(x) sweep(log(x), 1, rowMeans(log(x)), "-")
REL <- cbind(clr(met.rel), clr(mic.rel))
pr.r <- propr:::lr2rho(as.matrix(REL))
colnames(pr.r) <- colnames(REL)
rownames(pr.r) <- colnames(REL)

# Analyze multi-omics with wrapper
ivarNA <- propr(REL, ivar = NA)

test_that("ivar = NA works as expected", {

expect_equal(
NA,
ivarNA@ivar
)

expect_equal(
pr.r,
ivarNA@matrix
)

expect_equal(
ivarNA@counts,
ivarNA@logratio
)
})

0 comments on commit b85112e

Please sign in to comment.