# This file is part of the R package mkin # mkin is free software: you can redistribute it and/or modify it under the # terms of the GNU General Public License as published by the Free Software # Foundation, either version 3 of the License, or (at your option) any later # version. # This program is distributed in the hope that it will be useful, but WITHOUT # ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS # FOR A PARTICULAR PURPOSE. See the GNU General Public License for more # details. # You should have received a copy of the GNU General Public License along with # this program. If not, see <http://www.gnu.org/licenses/> #' Function to perform isometric log-ratio transformation #' #' This implementation is a special case of the class of isometric log-ratio #' transformations. #' #' @aliases ilr invilr #' @param x A numeric vector. Naturally, the forward transformation is only #' sensible for vectors with all elements being greater than zero. #' @return The result of the forward or backward transformation. The returned #' components always sum to 1 for the case of the inverse log-ratio #' transformation. #' @author René Lehmann and Johannes Ranke #' @seealso Another implementation can be found in R package #' \code{robCompositions}. #' @references Peter Filzmoser, Karel Hron (2008) Outlier Detection for #' Compositional Data Using Robust Methods. Math Geosci 40 233-248 #' @keywords manip #' @examples #' #' # Order matters #' ilr(c(0.1, 1, 10)) #' ilr(c(10, 1, 0.1)) #' # Equal entries give ilr transformations with zeros as elements #' ilr(c(3, 3, 3)) #' # Almost equal entries give small numbers #' ilr(c(0.3, 0.4, 0.3)) #' # Only the ratio between the numbers counts, not their sum #' invilr(ilr(c(0.7, 0.29, 0.01))) #' invilr(ilr(2.1 * c(0.7, 0.29, 0.01))) #' # Inverse transformation of larger numbers gives unequal elements #' invilr(-10) #' invilr(c(-10, 0)) #' # The sum of the elements of the inverse ilr is 1 #' sum(invilr(c(-10, 0))) #' # This is why we do not need all elements of the inverse transformation to go back: #' a <- c(0.1, 0.3, 0.5) #' b <- invilr(a) #' length(b) # Four elements #' ilr(c(b[1:3], 1 - sum(b[1:3]))) # Gives c(0.1, 0.3, 0.5) #' #' @export ilr <- function(x) { z <- vector() for (i in 1:(length(x) - 1)) { z[i] <- sqrt(i/(i+1)) * log((prod(x[1:i]))^(1/i) / x[i+1]) } return(z) } #' @rdname ilr #' @export invilr<-function(x) { D <- length(x) + 1 z <- c(x, 0) y <- rep(0, D) s <- sqrt(1:D*2:(D+1)) q <- z/s y[1] <- sum(q[1:D]) for (i in 2:D) { y[i] <- sum(q[i:D]) - sqrt((i-1)/i) * z[i-1] } z <- vector() for (i in 1:D) { z[i] <- exp(y[i])/sum(exp(y)) } # Work around a numerical problem with NaN values returned by the above # Only works if there is only one NaN value: replace it with 1 # if the sum of the other components is < 1e-10 if (sum(is.na(z)) == 1 && sum(z[!is.na(z)]) < 1e-10) z = ifelse(is.na(z), 1, z) return(z) }