From 57dcff2209ff3400d00bacbd5d875583d8cf8f36 Mon Sep 17 00:00:00 2001 From: Waldir Leoncio Date: Wed, 15 Jan 2025 11:07:44 +0100 Subject: [PATCH] Replaced `randga()` with `rgamma()` (closes #9) --- NAMESPACE | 1 + R/randdir.R | 4 ++-- R/randga.R | 61 -------------------------------------------------- man/randdir.Rd | 3 --- man/randga.Rd | 25 --------------------- 5 files changed, 3 insertions(+), 91 deletions(-) delete mode 100644 R/randga.R delete mode 100644 man/randga.Rd diff --git a/NAMESPACE b/NAMESPACE index 3094ed5..0b79163 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -36,6 +36,7 @@ importFrom(matlab2r,times) importFrom(matlab2r,uiputfile) importFrom(matlab2r,zeros) importFrom(methods,is) +importFrom(stats,rgamma) importFrom(stats,runif) importFrom(stats,sd) importFrom(utils,read.delim) diff --git a/R/randdir.R b/R/randdir.R index 7ee626d..066db2e 100644 --- a/R/randdir.R +++ b/R/randdir.R @@ -3,11 +3,11 @@ #' @examples rBAPS:::randdir(matrix(c(10, 30, 60), 3), 3) #' @param counts shape parameter #' @param nc number of rows on output -#' @seealso randga +#' @importFrom stats rgamma randdir <- function(counts, nc) { svar <- zeros(nc, 1) for (i in 1:nc) { - svar[i, 1] <- randga(counts[i, 1], 1) + svar[i, 1] <- rgamma(1L, shape = counts[i, 1], rate = 1) } svar <- svar / sum(svar) return(svar) diff --git a/R/randga.R b/R/randga.R deleted file mode 100644 index 24fb960..0000000 --- a/R/randga.R +++ /dev/null @@ -1,61 +0,0 @@ -#' @title Generates random number from a Gamma distribution -#' @description Generates one random number from shape parameter a and rate -#' parameter b -#' @param a shape -#' @param b rate -#' @return One realization of Gamma(a, b) -#' @details The generated random variable has mean a / b. It will be -#' positively-skewed for small values, but converges to a symmetric distribution -#' for very large numbers of a and b. -randga <- function(a, b) { - flag <- 0 - if (a > 1) { - c1 <- a - 1 - c2 <- (a - (1 / (6 * a))) / c1 - c3 <- 2 / c1 - c4 <- c3 + 2 - c5 <- 1 / sqrt(a) - U1 <- 1 - while (flag == 0) { - if (a <= 2.5) { - U1 <- rand() - U2 <- rand() - } else { - while (!(U1 > 0 & U1 < 1)) { - U1 <- rand() - U2 <- rand() - U1 <- U2 + c5 * (1 - 1.86 * U1) - } - } - W <- c2 * U2 / U1 - if (c3 * U1 + W + (1 / W) <= c4) { - flag <- 1 - g <- c1 * W / b - } else if (c3 * log(U1) - log(W) + W < 1) { - flag <- 1 - g <- c1 * W / b - } else { - U1 <- -1 - } - } - } else if (a == 1) { - g <- sum(-(1 / b) * log(rand(a, 1))) - } else { - while (flag == 0) { - U <- rand(2, 1) - if (U[1] > exp(1) / (a + exp(1))) { - g <- -log(((a + exp(1)) * (1 - U[1])) / (a * exp(1))) - if (U[2] <= g^(a - 1)) { - flag <- 1 - } - } else { - g <- ((a + exp(1)) * U[1] / ((exp(1))^(1 / a))) - if (U[2] <= exp(-g)) { - flag <- 1 - } - } - } - g <- g / b - } - return(g) -} diff --git a/man/randdir.Rd b/man/randdir.Rd index 59c69b4..812a160 100644 --- a/man/randdir.Rd +++ b/man/randdir.Rd @@ -20,6 +20,3 @@ Generates random numbers \examples{ rBAPS:::randdir(matrix(c(10, 30, 60), 3), 3) } -\seealso{ -randga -} diff --git a/man/randga.Rd b/man/randga.Rd deleted file mode 100644 index 7045249..0000000 --- a/man/randga.Rd +++ /dev/null @@ -1,25 +0,0 @@ -% Generated by roxygen2: do not edit by hand -% Please edit documentation in R/randga.R -\name{randga} -\alias{randga} -\title{Generates random number from a Gamma distribution} -\usage{ -randga(a, b) -} -\arguments{ -\item{a}{shape} - -\item{b}{rate} -} -\value{ -One realization of Gamma(a, b) -} -\description{ -Generates one random number from shape parameter a and rate -parameter b -} -\details{ -The generated random variable has mean a / b. It will be -positively-skewed for small values, but converges to a symmetric distribution -for very large numbers of a and b. -}