--- a
+++ b/R/dirichlet.R
@@ -0,0 +1,63 @@
+#' Generate a matrix of `n` rows of sets of probabilities 
+#' generated from the Dirichlet distribution, each row 
+#' summing to 1. Uses the alternative \eqn{\mu} \eqn{\phi} 
+#' parameterisation for the Dirichlet distribution, 
+#' representing means and precision.
+#' 
+#' @param n Number of sets of probabilities (defaults to 1)
+#' @param mu Vector of mean values for each probability in the set
+#' (defaults to c(0.001, 0.029. 0.7)). Must be greater than 0 and 
+#' finite, and contain at least two values.
+#' @param phi Parameter representing precision, where precision is 
+#' 1/variance. Must be positive and finite. Defaults to 10.
+#'
+#' @examples rdirichlet_alt(n = 3, mu = c(0.001, 0.029, 0.7), phi = 10)
+#' 
+#' @export 
+#' 
+#' @import checkmate
+#' @importFrom stats rgamma
+#' 
+rdirichlet_alt <- function(
+  n = 1, 
+  mu = c(0.3, 0.3, 0.3), 
+  phi = 10
+) {
+
+  ### Check and convert inputs
+
+  # True if n is "close to an integer"
+  checkmate::assert_integerish(
+    n, lower = 1, , upper = 10^7, len = 1, any.missing = FALSE
+  )
+
+  # mu should be a numeric vector in the range [0, 1)
+  checkmate::assert_vector(
+    mu, min.len = 2, strict = TRUE, any.missing = FALSE
+  )
+  checkmate::assert_numeric(
+    mu, lower = 10^-7, upper = 1 - 10^-7
+  )
+
+  # phi should be a positive number
+  checkmate::assert_number(phi, lower = 10^-7, finite = TRUE)
+
+  # Make n actually an integer
+  n <- as.integer(n)
+
+  # Shape parameter (alpha) from mu and phi; alpha_0 = phi
+  alpha <- phi * mu / mu[1]
+  no_probs <- length(mu)
+
+  # Dirichlet is a set of normalised independent gamma(alpha, 1)
+  draws_mx <- matrix(
+    stats::rgamma(n * no_probs, alpha), 
+    ncol = no_probs, 
+    byrow = TRUE
+  )
+  # Each set should sum to 1
+  draws_mx <- draws_mx / rowSums(draws_mx)
+
+  return(draws_mx)
+}
+