diff --git a/R/dist_spec.R b/R/dist_spec.R
index 37d368362..a3e1a8c09 100644
--- a/R/dist_spec.R
+++ b/R/dist_spec.R
@@ -116,9 +116,12 @@ discrete_pmf <- function(distribution =
 #' )
 #' dist1 + dist1
 #'
-#' # An uncertain gamma distribution with mean 3 and sd 2
+#' # An uncertain gamma distribution with shape and rate normally distributed
+#' # as Normal(3, 0.5) and Normal(2, 0.5) respectively
 #' dist2 <- Gamma(
-#'   mean = Normal(3, 0.5), sd = Normal(2, 0.5), max = 20
+#'   shape = Normal(3, 0.5),
+#'   rate = Normal(2, 0.5),
+#'   max = 20
 #' )
 #' dist1 + dist2
 `+.dist_spec` <- function(e1, e2) {
@@ -197,9 +200,12 @@ discrete_pmf <- function(distribution =
 #' )
 #' dist1 + dist1
 #'
-#' # An uncertain gamma distribution with mean 3 and sd 2
+#' # An uncertain gamma distribution with shape and rate normally distributed
+#' # as Normal(3, 0.5) and Normal(2, 0.5) respectively
 #' dist2 <- Gamma(
-#'   mean = Normal(3, 0.5), sd = Normal(2, 0.5), max = 20
+#'   shape = Normal(3, 0.5),
+#'   rate = Normal(2, 0.5),
+#'   max = 20
 #' )
 #' c(dist1, dist2)
 c.dist_spec <- function(...) {
@@ -257,9 +263,12 @@ c.dist_spec <- function(...) {
 #' dist1 <- LogNormal(mean = 5, sd = 1, max = 20)
 #' mean(dist1)
 #'
-#' # An uncertain gamma distribution with mean 3 and sd 2
+#' # An uncertain gamma distribution with shape and rate normally distributed
+#' # as Normal(3, 0.5) and Normal(2, 0.5) respectively
 #' dist2 <- Gamma(
-#'  mean = Normal(3, 0.5), sd = Normal(2, 0.5), max = 20
+#'   shape = Normal(3, 0.5),
+#'   rate = Normal(2, 0.5),
+#'   max = 20
 #' )
 #' mean(dist2)
 #'
@@ -393,8 +402,13 @@ sd.default <- function(x, ...) {
 #' dist1 <- Gamma(mean = 5, sd = 1, max = 20)
 #' max(dist1)
 #'
-#' # An uncertain lognormal distribution with mean 3 and sd 2
-#' dist2 <- LogNormal(mean = Normal(3, 0.5), sd = Normal(2, 0.5), max = 20)
+#' # An uncertain lognormal distribution with meanlog and sdlog normally
+#' # distributed as Normal(3, 0.5) and Normal(2, 0.5) respectively
+#' dist2 <- LogNormal(
+#'   meanlog = Normal(3, 0.5),
+#'   sdlog = Normal(2, 0.5),
+#'   max = 20
+#' )
 #' max(dist2)
 #'
 #' # The max the sum of two distributions
@@ -447,8 +461,13 @@ discretise <- function(x, ...) {
 #' # A fixed gamma distribution with mean 5 and sd 1.
 #' dist1 <- Gamma(mean = 5, sd = 1, max = 20)
 #'
-#' # An uncertain lognormal distribution with mean 3 and sd 2
-#' dist2 <- LogNormal(mean = Normal(3, 0.5), sd = Normal(2, 0.5), max = 20)
+#' # An uncertain lognormal distribution with meanlog and sdlog normally
+#' # distributed as Normal(3, 0.5) and Normal(2, 0.5) respectively
+#' dist2 <- LogNormal(
+#'   meanlog = Normal(3, 0.5),
+#'   sdlog = Normal(2, 0.5),
+#'   max = 20
+#' )
 #'
 #' # The maxf the sum of two distributions
 #' discretise(dist1 + dist2, strict = FALSE)
@@ -531,8 +550,13 @@ collapse <- function(x, ...) {
 #' # A fixed gamma distribution with mean 5 and sd 1.
 #' dist1 <- Gamma(mean = 5, sd = 1, max = 20)
 #'
-#' # An uncertain lognormal distribution with mean 3 and sd 2
-#' dist2 <- LogNormal(mean = 3, sd = 2, max = 20)
+#' # An uncertain lognormal distribution with meanlog and sdlog normally
+#' # distributed as Normal(3, 0.5) and Normal(2, 0.5) respectively
+#' dist2 <- LogNormal(
+#'   meanlog = Normal(3, 0.5),
+#'   sdlog = Normal(2, 0.5),
+#'   max = 20
+#' )
 #'
 #' # The maxf the sum of two distributions
 #' collapse(discretise(dist1 + dist2))
@@ -588,9 +612,10 @@ collapse.multi_dist_spec <- function(x, ...) {
 #' dist1 <- LogNormal(mean = 1.5, sd = 0.5, max = 20)
 #' print(dist1)
 #'
-#' # An uncertain gamma distribution with mean 3 and sd 2
+#' # An uncertain gamma distribution with shape and rate normally distributed
+#' # as Normal(3, 0.5) and Normal(2, 0.5) respectively
 #' dist2 <- Gamma(
-#'   mean = Normal(3, 0.5), sd = Normal(2, 0.5), max = 20
+#'   shape = Normal(3, 0.5), rate = Normal(2, 0.5), max = 20
 #' )
 #' print(dist2)
 print.dist_spec <- function(x, ...) {
@@ -680,9 +705,12 @@ print.dist_spec <- function(x, ...) {
 #' # Plot discretised distribution with 0.01 day discretisation window
 #' plot(dist1, res = 0.01, cumulative = FALSE)
 #'
-#' # An uncertain gamma distribution with mean 3 and sd 2
+#' # An uncertain gamma distribution with shape and rate normally distributed
+#' # as Normal(3, 0.5) and Normal(2, 0.5) respectively
 #' dist2 <- Gamma(
-#'   mean = Normal(3, 0.5), sd = Normal(2, 0.5), max = 20
+#'   shape = Normal(3, 0.5),
+#'   rate = Normal(2, 0.5),
+#'   max = 20
 #' )
 #' plot(dist2)
 #'
@@ -781,9 +809,12 @@ plot.dist_spec <- function(x, samples = 50L, res = 1, cumulative = TRUE, ...) {
 #' @examples
 #' dist1 <- LogNormal(mean = 1.6, sd = 0.5, max = 20)
 #'
-#' # An uncertain gamma distribution with mean 3 and sd 2
+#' # An uncertain gamma distribution with shape and rate normally distributed
+#' # as Normal(3, 0.5) and Normal(2, 0.5) respectively
 #' dist2 <- Gamma(
-#'   mean = Normal(3, 0.5), sd = Normal(2, 0.5), max = 20
+#'   shape = Normal(3, 0.5),
+#'   rate = Normal(2, 0.5),
+#'   max = 20
 #' )
 #'
 #' # Multiple distributions
@@ -828,9 +859,12 @@ fix_parameters <- function(x, ...) {
 #' @importFrom rlang arg_match
 #' @method fix_parameters dist_spec
 #' @examples
-#' # An uncertain gamma distribution with mean 3 and sd 2
-#' dist <- LogNormal(
-#'   meanlog = Normal(3, 0.5), sdlog = Normal(2, 0.5), max = 20
+#' # An uncertain gamma distribution with shape and rate normally distributed
+#' # as Normal(3, 0.5) and Normal(2, 0.5) respectively
+#' dist <- Gamma(
+#'   shape = Normal(3, 0.5),
+#'   rate = Normal(2, 0.5),
+#'   max = 20
 #' )
 #'
 #' fix_parameters(dist)
@@ -891,8 +925,13 @@ is_constrained <- function(x, ...) {
 #' # A fixed gamma distribution with mean 5 and sd 1.
 #' dist1 <- Gamma(mean = 5, sd = 1, max = 20)
 #'
-#' # An uncertain lognormal distribution with mean 3 and sd 2
-#' dist2 <- LogNormal(mean = Normal(3, 0.5), sd = Normal(2, 0.5), max = 20)
+#' # An uncertain lognormal distribution with meanlog and sdlog normally
+#' # distributed as Normal(3, 0.5) and Normal(2, 0.5) respectively
+#' dist2 <- LogNormal(
+#'   meanlog = Normal(3, 0.5),
+#'   sdlog = Normal(2, 0.5),
+#'   max = 20
+#' )
 #'
 #' # both distributions are constrained and therefore so is the sum
 #' is_constrained(dist1 + dist2)
@@ -952,7 +991,8 @@ is_constrained.multi_dist_spec <- function(x, ...) {
 #' @examples
 #' LogNormal(mean = 4, sd = 1)
 #' LogNormal(mean = 4, sd = 1, max = 10)
-#' LogNormal(mean = Normal(4, 1), sd = 1, max = 10)
+#' # If specifying uncertain parameters, use the natural parameters
+#' LogNormal(meanlog = Normal(4, 1), sdlog = 1, max = 10)
 LogNormal <- function(meanlog, sdlog, mean, sd, ...) {
   params <- as.list(environment())
   return(new_dist_spec(params, "lognormal", ...))
diff --git a/R/opts.R b/R/opts.R
index 04d9e5577..613a6452f 100644
--- a/R/opts.R
+++ b/R/opts.R
@@ -42,8 +42,8 @@
 #' # An uncertain gamma distributed generation time
 #' generation_time_opts(
 #'   Gamma(
-#'     mean = Normal(mean = 3, sd = 1),
-#'     sd = Normal(mean = 2, sd = 0.5),
+#'     shape = Normal(mean = 3, sd = 1),
+#'     rate = Normal(mean = 2, sd = 0.5),
 #'     max = 14
 #'   )
 #' )
@@ -186,7 +186,11 @@ secondary_opts <- function(type = c("incidence", "prevalence"), ...) {
 #' delay_opts()
 #'
 #' # A single delay that has uncertainty
-#' delay <- LogNormal(mean = Normal(1, 0.2), sd = Normal(0.5, 0.1), max = 14)
+#' delay <- LogNormal(
+#'   meanlog = Normal(1, 0.2),
+#'   sdlog = Normal(0.5, 0.1),
+#'   max = 14
+#' )
 #' delay_opts(delay)
 #'
 #' # A single delay without uncertainty
diff --git a/man/Distributions.Rd b/man/Distributions.Rd
index 583b0193d..6db8e8911 100644
--- a/man/Distributions.Rd
+++ b/man/Distributions.Rd
@@ -77,7 +77,8 @@ probability mass function is given directly as a numeric vector.
 \examples{
 LogNormal(mean = 4, sd = 1)
 LogNormal(mean = 4, sd = 1, max = 10)
-LogNormal(mean = Normal(4, 1), sd = 1, max = 10)
+# If specifying uncertain parameters, use the natural parameters
+LogNormal(meanlog = Normal(4, 1), sdlog = 1, max = 10)
 Gamma(mean = 4, sd = 1)
 Gamma(shape = 16, rate = 4)
 Gamma(shape = Normal(16, 2), rate = Normal(4, 1))
diff --git a/man/c.dist_spec.Rd b/man/c.dist_spec.Rd
index 45cffbf55..ae603e9e2 100644
--- a/man/c.dist_spec.Rd
+++ b/man/c.dist_spec.Rd
@@ -27,9 +27,12 @@ dist1 <- LogNormal(
 )
 dist1 + dist1
 
-# An uncertain gamma distribution with mean 3 and sd 2
+# An uncertain gamma distribution with shape and rate normally distributed
+# as Normal(3, 0.5) and Normal(2, 0.5) respectively
 dist2 <- Gamma(
-  mean = Normal(3, 0.5), sd = Normal(2, 0.5), max = 20
+  shape = Normal(3, 0.5),
+  rate = Normal(2, 0.5),
+  max = 20
 )
 c(dist1, dist2)
 }
diff --git a/man/collapse.Rd b/man/collapse.Rd
index 97ed0d1c2..e3d67be00 100644
--- a/man/collapse.Rd
+++ b/man/collapse.Rd
@@ -25,8 +25,13 @@ in the <dist_spec>.
 # A fixed gamma distribution with mean 5 and sd 1.
 dist1 <- Gamma(mean = 5, sd = 1, max = 20)
 
-# An uncertain lognormal distribution with mean 3 and sd 2
-dist2 <- LogNormal(mean = 3, sd = 2, max = 20)
+# An uncertain lognormal distribution with meanlog and sdlog normally
+# distributed as Normal(3, 0.5) and Normal(2, 0.5) respectively
+dist2 <- LogNormal(
+  meanlog = Normal(3, 0.5),
+  sdlog = Normal(2, 0.5),
+  max = 20
+)
 
 # The maxf the sum of two distributions
 collapse(discretise(dist1 + dist2))
diff --git a/man/delay_opts.Rd b/man/delay_opts.Rd
index b28a47b72..dc383a541 100644
--- a/man/delay_opts.Rd
+++ b/man/delay_opts.Rd
@@ -44,7 +44,11 @@ functions.
 delay_opts()
 
 # A single delay that has uncertainty
-delay <- LogNormal(mean = Normal(1, 0.2), sd = Normal(0.5, 0.1), max = 14)
+delay <- LogNormal(
+  meanlog = Normal(1, 0.2),
+  sdlog = Normal(0.5, 0.1),
+  max = 14
+)
 delay_opts(delay)
 
 # A single delay without uncertainty
diff --git a/man/discretise.Rd b/man/discretise.Rd
index e6c2cda1a..a76646d4a 100644
--- a/man/discretise.Rd
+++ b/man/discretise.Rd
@@ -42,8 +42,13 @@ from the difference of two censored events.
 # A fixed gamma distribution with mean 5 and sd 1.
 dist1 <- Gamma(mean = 5, sd = 1, max = 20)
 
-# An uncertain lognormal distribution with mean 3 and sd 2
-dist2 <- LogNormal(mean = Normal(3, 0.5), sd = Normal(2, 0.5), max = 20)
+# An uncertain lognormal distribution with meanlog and sdlog normally
+# distributed as Normal(3, 0.5) and Normal(2, 0.5) respectively
+dist2 <- LogNormal(
+  meanlog = Normal(3, 0.5),
+  sdlog = Normal(2, 0.5),
+  max = 20
+)
 
 # The maxf the sum of two distributions
 discretise(dist1 + dist2, strict = FALSE)
diff --git a/man/extract_single_dist.Rd b/man/extract_single_dist.Rd
index e3407cf6b..e6761e96c 100644
--- a/man/extract_single_dist.Rd
+++ b/man/extract_single_dist.Rd
@@ -20,9 +20,12 @@ A single \code{dist_spec} object
 \examples{
 dist1 <- LogNormal(mean = 1.6, sd = 0.5, max = 20)
 
-# An uncertain gamma distribution with mean 3 and sd 2
+# An uncertain gamma distribution with shape and rate normally distributed
+# as Normal(3, 0.5) and Normal(2, 0.5) respectively
 dist2 <- Gamma(
-  mean = Normal(3, 0.5), sd = Normal(2, 0.5), max = 20
+  shape = Normal(3, 0.5),
+  rate = Normal(2, 0.5),
+  max = 20
 )
 
 # Multiple distributions
diff --git a/man/fix_parameters.Rd b/man/fix_parameters.Rd
index 04a8a4136..a43df5bbd 100644
--- a/man/fix_parameters.Rd
+++ b/man/fix_parameters.Rd
@@ -25,9 +25,12 @@ If the given \verb{<dist_spec>} has any uncertainty, it is removed and the
 corresponding distribution converted into a fixed one.
 }
 \examples{
-# An uncertain gamma distribution with mean 3 and sd 2
-dist <- LogNormal(
-  meanlog = Normal(3, 0.5), sdlog = Normal(2, 0.5), max = 20
+# An uncertain gamma distribution with shape and rate normally distributed
+# as Normal(3, 0.5) and Normal(2, 0.5) respectively
+dist <- Gamma(
+  shape = Normal(3, 0.5),
+  rate = Normal(2, 0.5),
+  max = 20
 )
 
 fix_parameters(dist)
diff --git a/man/generation_time_opts.Rd b/man/generation_time_opts.Rd
index 1200101c7..8306d7a61 100644
--- a/man/generation_time_opts.Rd
+++ b/man/generation_time_opts.Rd
@@ -79,8 +79,8 @@ generation_time_opts(Gamma(mean = 3, sd = 2, max = 14))
 # An uncertain gamma distributed generation time
 generation_time_opts(
   Gamma(
-    mean = Normal(mean = 3, sd = 1),
-    sd = Normal(mean = 2, sd = 0.5),
+    shape = Normal(mean = 3, sd = 1),
+    rate = Normal(mean = 2, sd = 0.5),
     max = 14
   )
 )
diff --git a/man/is_constrained.Rd b/man/is_constrained.Rd
index 7a097ab6d..39803ee62 100644
--- a/man/is_constrained.Rd
+++ b/man/is_constrained.Rd
@@ -23,8 +23,13 @@ Logical; TRUE if \code{x} is constrained
 # A fixed gamma distribution with mean 5 and sd 1.
 dist1 <- Gamma(mean = 5, sd = 1, max = 20)
 
-# An uncertain lognormal distribution with mean 3 and sd 2
-dist2 <- LogNormal(mean = Normal(3, 0.5), sd = Normal(2, 0.5), max = 20)
+# An uncertain lognormal distribution with meanlog and sdlog normally
+# distributed as Normal(3, 0.5) and Normal(2, 0.5) respectively
+dist2 <- LogNormal(
+  meanlog = Normal(3, 0.5),
+  sdlog = Normal(2, 0.5),
+  max = 20
+)
 
 # both distributions are constrained and therefore so is the sum
 is_constrained(dist1 + dist2)
diff --git a/man/max.dist_spec.Rd b/man/max.dist_spec.Rd
index c81f811ef..516a0cc3a 100644
--- a/man/max.dist_spec.Rd
+++ b/man/max.dist_spec.Rd
@@ -25,8 +25,13 @@ in parameters)
 dist1 <- Gamma(mean = 5, sd = 1, max = 20)
 max(dist1)
 
-# An uncertain lognormal distribution with mean 3 and sd 2
-dist2 <- LogNormal(mean = Normal(3, 0.5), sd = Normal(2, 0.5), max = 20)
+# An uncertain lognormal distribution with meanlog and sdlog normally
+# distributed as Normal(3, 0.5) and Normal(2, 0.5) respectively
+dist2 <- LogNormal(
+  meanlog = Normal(3, 0.5),
+  sdlog = Normal(2, 0.5),
+  max = 20
+)
 max(dist2)
 
 # The max the sum of two distributions
diff --git a/man/mean.dist_spec.Rd b/man/mean.dist_spec.Rd
index 5448ca1b8..8be98f2d0 100644
--- a/man/mean.dist_spec.Rd
+++ b/man/mean.dist_spec.Rd
@@ -25,9 +25,12 @@ distributions combined in the passed <dist_spec>.
 dist1 <- LogNormal(mean = 5, sd = 1, max = 20)
 mean(dist1)
 
-# An uncertain gamma distribution with mean 3 and sd 2
+# An uncertain gamma distribution with shape and rate normally distributed
+# as Normal(3, 0.5) and Normal(2, 0.5) respectively
 dist2 <- Gamma(
- mean = Normal(3, 0.5), sd = Normal(2, 0.5), max = 20
+  shape = Normal(3, 0.5),
+  rate = Normal(2, 0.5),
+  max = 20
 )
 mean(dist2)
 
diff --git a/man/plot.dist_spec.Rd b/man/plot.dist_spec.Rd
index 27698649a..7e1a8d19e 100644
--- a/man/plot.dist_spec.Rd
+++ b/man/plot.dist_spec.Rd
@@ -33,9 +33,12 @@ plot(dist1)
 # Plot discretised distribution with 0.01 day discretisation window
 plot(dist1, res = 0.01, cumulative = FALSE)
 
-# An uncertain gamma distribution with mean 3 and sd 2
+# An uncertain gamma distribution with shape and rate normally distributed
+# as Normal(3, 0.5) and Normal(2, 0.5) respectively
 dist2 <- Gamma(
-  mean = Normal(3, 0.5), sd = Normal(2, 0.5), max = 20
+  shape = Normal(3, 0.5),
+  rate = Normal(2, 0.5),
+  max = 20
 )
 plot(dist2)
 
diff --git a/man/plus-.dist_spec.Rd b/man/plus-.dist_spec.Rd
index 16e967faa..091b7e536 100644
--- a/man/plus-.dist_spec.Rd
+++ b/man/plus-.dist_spec.Rd
@@ -26,9 +26,12 @@ dist1 <- LogNormal(
 )
 dist1 + dist1
 
-# An uncertain gamma distribution with mean 3 and sd 2
+# An uncertain gamma distribution with shape and rate normally distributed
+# as Normal(3, 0.5) and Normal(2, 0.5) respectively
 dist2 <- Gamma(
-  mean = Normal(3, 0.5), sd = Normal(2, 0.5), max = 20
+  shape = Normal(3, 0.5),
+  rate = Normal(2, 0.5),
+  max = 20
 )
 dist1 + dist2
 }
diff --git a/man/print.dist_spec.Rd b/man/print.dist_spec.Rd
index c3ea43cfb..b6925af9f 100644
--- a/man/print.dist_spec.Rd
+++ b/man/print.dist_spec.Rd
@@ -24,9 +24,10 @@ functions of fixed delay distributions combined in the passed <dist_spec>.
 dist1 <- LogNormal(mean = 1.5, sd = 0.5, max = 20)
 print(dist1)
 
-# An uncertain gamma distribution with mean 3 and sd 2
+# An uncertain gamma distribution with shape and rate normally distributed
+# as Normal(3, 0.5) and Normal(2, 0.5) respectively
 dist2 <- Gamma(
-  mean = Normal(3, 0.5), sd = Normal(2, 0.5), max = 20
+  shape = Normal(3, 0.5), rate = Normal(2, 0.5), max = 20
 )
 print(dist2)
 }