diff --git a/DESCRIPTION b/DESCRIPTION
index 8f7a632..81d099a 100644
--- a/DESCRIPTION
+++ b/DESCRIPTION
@@ -33,6 +33,7 @@ Imports:
rstantools (>= 2.4.0),
dplyr,
magrittr,
+ matrixStats,
purrr,
stringr,
R6 (>= 2.4.1)
diff --git a/NAMESPACE b/NAMESPACE
index 340e413..2329fbf 100644
--- a/NAMESPACE
+++ b/NAMESPACE
@@ -8,6 +8,7 @@ export(check_numeric)
export(distribution.binomial)
export(distribution.exponential)
export(distribution.gamma)
+export(distribution.lognormal)
export(distribution.mixture)
export(distribution.negative_binomial)
export(distribution.normal)
@@ -36,6 +37,7 @@ importFrom(data.table,set)
importFrom(data.table,setnames)
importFrom(dplyr,if_else)
importFrom(magrittr,"%>%")
+importFrom(matrixStats,colMaxs)
importFrom(purrr,map_chr)
importFrom(purrr,map_dbl)
importFrom(rstan,sampling)
diff --git a/R/distribution_continuous.R b/R/distribution_continuous.R
index c888ae8..810c67f 100644
--- a/R/distribution_continuous.R
+++ b/R/distribution_continuous.R
@@ -701,4 +701,137 @@ distribution.continuous.normal.class <- R6.class(
distribution.normal <- function( mean, sd ){
distribution.continuous.normal.class$new( mean = mean,
sd = sd )
+}
+
+################################################################################/
+# distribution.continuous.lognormal
+################################################################################/
+#' Class: `distribution.continuous.lognormal.class`
+#' @description Derived class for a lognormally-distributed random variable.
+#'
+#' @param meanlog the mean of log(X)
+#' @param sdlog the standard deviation of log(X)
+#' @param x vector of quantiles.
+#' @param q vector of quantiles.
+#' @param p vector of probabilities.
+#' @param n number of observations. If `length( n ) > 1`, the length is
+#' taken to be the number required.
+#' @param log logical; if TRUE, probabilities p are given as `log(p)`.
+#' @param log.p logical; if TRUE, probabilities p are given as `log(p)`.
+#' @param lower.tail logical; if TRUE (default), probabilities are
+#' \eqn{P[ X \leq x ]}, otherwise, \eqn{P[X>x]}.
+#'
+#' @field interfaces The list of available class interfaces
+#' @field mean the mean of the distribution
+#' @field sd the standard deviation of the distribution
+#' @field var the variance of the distribution
+distribution.continuous.lognormal.class <- R6.class(
+ classname = "distribution.continuous.lognormal.class",
+ inherit = distribution.continuous.class,
+ private = list(
+ .name = "lognormal",
+ .param_names = c( "meanlog", "sdlog" ),
+ .check_params = function( params ){
+ # Check that params contains all elements of private$.param_names
+ super$.check_params( params )
+
+ if ( !is.numeric( params$meanlog ) )
+ stop( "`params$meanlog` must be a numeric value.")
+ if ( params$sdlog < 0 || !is.numeric( params$sdlog ) )
+ stop( "`params$sdlog` must be a non-negative numeric value.")
+
+ return( NULL )
+ }
+ ),
+ public = list(
+ ############################################################################/
+ # initialize
+ ############################################################################/
+ #' @description Create a new object of class
+ #' `distribution.continuous.normal.class`
+ initialize = function( meanlog, sdlog ){
+ super$initialize( support = c( 0, Inf ) )
+ self$params <- list( meanlog = meanlog,
+ sdlog = sdlog )
+ },
+ ############################################################################/
+ # density
+ ############################################################################/
+ #' @description Density function for a lognormal random variable with mean log(X)
+ #' `$params$meanlog` and standard deviation log(X) `$params$sdlog`.
+ d = function( x, log = FALSE ){
+ stats::dlnorm( x, meanlog = private$.params$meanlog, sdlog = private$.params$sdlog,
+ log = log )
+ },
+ ############################################################################/
+ # distribution function
+ ############################################################################/
+ #' @description Cumulative density function for a lognormal random variable
+ #' with mean log(X) `$params$meanlog` and standard deviation log(X) `$params$sdlog`.
+ p = function( q, lower.tail = TRUE, log.p = FALSE ){
+ stats::plnorm( q, meanlog = private$.params$meanlog, sdlog = private$.params$sdlog,
+ lower.tail = lower.tail, log.p = log.p )
+ },
+ ############################################################################/
+ # quantile function
+ ############################################################################/
+ #' @description Quantile function for a lognormal random variable with mean
+ #' log(X) `$params$meanlog` and standard deviation log(X) `$params$sdlog`.
+ q = function( p, lower.tail = TRUE, log.p = FALSE ){
+ stats::qlnorm( p, meanlog = private$.params$meanlog, sdlog = private$.params$sdlog,
+ lower.tail = lower.tail, log.p = log.p )
+ },
+ ############################################################################/
+ # random deviates
+ ############################################################################/
+ #' @description Generates random deviates for a lognormal random variable with
+ #' mean log(X) `$params$meanlog` and standard deviation logx(X) `$params$sdlog` .
+ r = function( n ){
+ stats::rlnorm( n, meanlog = private$.params$meanlog, sdlog = private$.params$sdlog )
+ }
+ ),
+ active = list(
+ ############################################################################/
+ # mean
+ ############################################################################/
+ mean = function( val ){
+ if( !missing( val ) )
+ stop( "cannot set `$mean`" )
+ return( exp( private$.params$meanlog + 0.5 * private$.params$sdlog^2 ) )
+ },
+ ############################################################################/
+ # standard deviation
+ ############################################################################/
+ sd = function( val ){
+ if( !missing( val ) )
+ stop( "cannot set `$sd`" )
+ return( sqrt( self$var ) )
+ },
+ ############################################################################/
+ # variance
+ ############################################################################/
+ var = function( val ){
+ if( !missing( val ) )
+ stop( "cannot set `$var`" )
+ sigma2 <- private$.params$sdlog^2
+ mu <- private$.params$meanlog
+ return( ( exp( sigma2) - 1 ) * exp( 2 * mu + sigma2 ) )
+ }
+ )
+)
+
+#' distribution.lognormal
+#'
+#' Constructor function for an object of class `distribution.continuous.lognormal.class`
+#'
+#' @param meanlog The mean of log(X)
+#' @param sdlog The standard deviation log(X)
+#'
+#' @returns An object of class [[distribution.continuous.lognormal.class]]
+#'
+#' @seealso [Mastiff-Distributions]
+#' @export
+distribution.lognormal <- function( meanlog, sdlog ){
+ distribution.continuous.lognormal.class$new( meanlog = meanlog,
+ sdlog = sdlog )
}
\ No newline at end of file
diff --git a/R/distribution_mixture.R b/R/distribution_mixture.R
index 72d4c6b..98b40ed 100644
--- a/R/distribution_mixture.R
+++ b/R/distribution_mixture.R
@@ -7,11 +7,22 @@
#' Class: `distribution.mixture.class`
#' @description Class to describe the mixture of distributions
#'
+#' @param distributions list of distributions in the mixture
+#' @param weights vector of weights of the distributions in the mixture
+#' @param x vector of quantiles.
#' @param q vector of quantiles.
#' @param p vector of probabilities.
+#' @param n number of observations. If `length( n ) > 1`, the length is
+#' taken to be the number required.
+#' @param log logical; if TRUE, probabilities p are given as `log(p)`.
#' @param log.p logical; if TRUE, probabilities p are given as `log(p)`.
#' @param lower.tail logical; if TRUE (default), probabilities are \eqn{P[ X \leq x ]},
#' otherwise, \eqn{P[X>x]}.
+#'
+#' @field n_distributions the number of distributions in the mixture
+#' @field distributions the distributions in the mixture
+#' @field weights the weights of the distributions in the mixture
+#' @field interfaces the list of available class interfaces
#'
#' @include R6_class.R
#' @include distribution_R6_class.R
@@ -21,7 +32,16 @@ distribution.mixture.class <- R6.class(
private = list(
.distributions = NULL,
.n_distributions = NULL,
- .weights = NULL
+ .weights = NULL,
+ .logExpAdd = function( ds, w ) {
+ eps = 1e-16
+ nw = ncol( ds )
+ ns = nrow( ds )
+ ds <- t( ds ) + matrix( log( w + eps ), nrow = nw, ncol = ns )
+ max <- matrixStats::colMaxs( ds )
+ ds <- ds - matrix( rep( max, each = nw ), nrow = nw, ncol = ns )
+ return( log( colSums( exp( ds ) ) ) + max )
+ }
),
public = list(
############################################################################/
@@ -46,9 +66,14 @@ distribution.mixture.class <- R6.class(
#' @description Density function for a random variable of the mixture
d = function( x, log = FALSE ){
ds <- matrix( unlist( lapply( private$.distributions,
- function( d ) d$d( x ) ) ),
+ function( d ) d$d( x, log = log ) ) ),
ncol = self$n_distributions )
- return( ( ds %*% self$weights )[, 1 ] )
+
+ if( log ) {
+ return( private$.logExpAdd( ds, self$weights ) )
+ } else {
+ return( ( ds %*% self$weights )[, 1 ] )
+ }
},
##############################################################################/
# cumulative distribution function
@@ -58,7 +83,12 @@ distribution.mixture.class <- R6.class(
ps <- matrix( unlist( lapply( private$.distributions,
function( d ) d$p( q, lower.tail = lower.tail, log.p = log.p ) ) ),
ncol = self$n_distributions )
- return( ( ps %*% self$weights )[, 1 ] )
+
+ if( log.p ) {
+ return( private$.logExpAdd( ps, self$weights ) )
+ } else {
+ return( ( ps %*% self$weights )[, 1 ] )
+ }
},
##############################################################################/
# quantile function
@@ -66,10 +96,6 @@ distribution.mixture.class <- R6.class(
#' @description Evaluates the quantile function of the mixture
q = function( p, lower.tail = TRUE, log.p = FALSE ){
super$q( p, lower.tail = lower.tail, log.p = log.p )
- # qs <- matrix( unlist( lapply( private$.distributions,
- # function( d ) d$q( p, lower.tail = lower.tail, log.p = log.p ) ) ),
- # ncol = self$n_distributions )
- # return( ( qs %*% self$weights )[, 1 ] )
},
############################################################################/
# random deviates
@@ -90,20 +116,15 @@ distribution.mixture.class <- R6.class(
}
),
active = list(
- support = function( val ){
- private$.staticReturn( val, "support" )
- },
############################################################################/
# n_distributions
############################################################################/
- #' @description Number of distributions in the mixture
n_distributions = function( val ){
private$.staticReturn( val, "n_distributions" )
},
############################################################################/
# distributions
############################################################################/
- #' @description The distributions of the mixture
distributions = function( val ){
if( missing( val ) )
return( private$.staticReturn( val, "distributions" ) )
@@ -121,7 +142,6 @@ distribution.mixture.class <- R6.class(
############################################################################/
# weights
############################################################################/
- #' @description The weights of each distribution in the mixture
weights = function( new = NA ){
if( length( new) == 1 )
if( is.na( new ) ) {
diff --git a/man/distribution.continuous.lognormal.class.Rd b/man/distribution.continuous.lognormal.class.Rd
new file mode 100644
index 0000000..0288d30
--- /dev/null
+++ b/man/distribution.continuous.lognormal.class.Rd
@@ -0,0 +1,158 @@
+% Generated by roxygen2: do not edit by hand
+% Please edit documentation in R/distribution_continuous.R
+\name{distribution.continuous.lognormal.class}
+\alias{distribution.continuous.lognormal.class}
+\title{Class: \code{distribution.continuous.lognormal.class}}
+\description{
+Derived class for a lognormally-distributed random variable.
+}
+\section{Super classes}{
+\code{mastiff::R6.class.class} -> \code{\link[mastiff:distribution.abstract.class]{mastiff::distribution.abstract.class}} -> \code{\link[mastiff:distribution.continuous.class]{mastiff::distribution.continuous.class}} -> \code{distribution.continuous.lognormal.class}
+}
+\section{Active bindings}{
+\if{html}{\out{}}
+\describe{
+\item{\code{interfaces}}{The list of available class interfaces}
+
+\item{\code{mean}}{the mean of the distribution}
+
+\item{\code{sd}}{the standard deviation of the distribution}
+
+\item{\code{var}}{the variance of the distribution}
+}
+\if{html}{\out{
}}
+}
+\section{Methods}{
+\subsection{Public methods}{
+\itemize{
+\item \href{#method-distribution.continuous.lognormal.class-new}{\code{distribution.continuous.lognormal.class$new()}}
+\item \href{#method-distribution.continuous.lognormal.class-d}{\code{distribution.continuous.lognormal.class$d()}}
+\item \href{#method-distribution.continuous.lognormal.class-p}{\code{distribution.continuous.lognormal.class$p()}}
+\item \href{#method-distribution.continuous.lognormal.class-q}{\code{distribution.continuous.lognormal.class$q()}}
+\item \href{#method-distribution.continuous.lognormal.class-r}{\code{distribution.continuous.lognormal.class$r()}}
+\item \href{#method-distribution.continuous.lognormal.class-clone}{\code{distribution.continuous.lognormal.class$clone()}}
+}
+}
+\if{html}{\out{
}}
+\if{html}{\out{}}
+\if{latex}{\out{\hypertarget{method-distribution.continuous.lognormal.class-new}{}}}
+\subsection{Method \code{new()}}{
+Create a new object of class
+\code{distribution.continuous.normal.class}
+\subsection{Usage}{
+\if{html}{\out{}}\preformatted{distribution.continuous.lognormal.class$new(meanlog, sdlog)}\if{html}{\out{
}}
+}
+
+\subsection{Arguments}{
+\if{html}{\out{}}
+\describe{
+\item{\code{meanlog}}{the mean of log(X)}
+
+\item{\code{sdlog}}{the standard deviation of log(X)}
+}
+\if{html}{\out{
}}
+}
+}
+\if{html}{\out{
}}
+\if{html}{\out{}}
+\if{latex}{\out{\hypertarget{method-distribution.continuous.lognormal.class-d}{}}}
+\subsection{Method \code{d()}}{
+Density function for a lognormal random variable with mean log(X)
+\verb{$params$meanlog} and standard deviation log(X) \verb{$params$sdlog}.
+\subsection{Usage}{
+\if{html}{\out{}}\preformatted{distribution.continuous.lognormal.class$d(x, log = FALSE)}\if{html}{\out{
}}
+}
+
+\subsection{Arguments}{
+\if{html}{\out{}}
+\describe{
+\item{\code{x}}{vector of quantiles.}
+
+\item{\code{log}}{logical; if TRUE, probabilities p are given as \code{log(p)}.}
+}
+\if{html}{\out{
}}
+}
+}
+\if{html}{\out{
}}
+\if{html}{\out{}}
+\if{latex}{\out{\hypertarget{method-distribution.continuous.lognormal.class-p}{}}}
+\subsection{Method \code{p()}}{
+Cumulative density function for a lognormal random variable
+with mean log(X) \verb{$params$meanlog} and standard deviation log(X) \verb{$params$sdlog}.
+\subsection{Usage}{
+\if{html}{\out{}}\preformatted{distribution.continuous.lognormal.class$p(q, lower.tail = TRUE, log.p = FALSE)}\if{html}{\out{
}}
+}
+
+\subsection{Arguments}{
+\if{html}{\out{}}
+\describe{
+\item{\code{q}}{vector of quantiles.}
+
+\item{\code{lower.tail}}{logical; if TRUE (default), probabilities are
+\eqn{P[ X \leq x ]}, otherwise, \eqn{P[X>x]}.}
+
+\item{\code{log.p}}{logical; if TRUE, probabilities p are given as \code{log(p)}.}
+}
+\if{html}{\out{
}}
+}
+}
+\if{html}{\out{
}}
+\if{html}{\out{}}
+\if{latex}{\out{\hypertarget{method-distribution.continuous.lognormal.class-q}{}}}
+\subsection{Method \code{q()}}{
+Quantile function for a lognormal random variable with mean
+log(X) \verb{$params$meanlog} and standard deviation log(X) \verb{$params$sdlog}.
+\subsection{Usage}{
+\if{html}{\out{}}\preformatted{distribution.continuous.lognormal.class$q(p, lower.tail = TRUE, log.p = FALSE)}\if{html}{\out{
}}
+}
+
+\subsection{Arguments}{
+\if{html}{\out{}}
+\describe{
+\item{\code{p}}{vector of probabilities.}
+
+\item{\code{lower.tail}}{logical; if TRUE (default), probabilities are
+\eqn{P[ X \leq x ]}, otherwise, \eqn{P[X>x]}.}
+
+\item{\code{log.p}}{logical; if TRUE, probabilities p are given as \code{log(p)}.}
+}
+\if{html}{\out{
}}
+}
+}
+\if{html}{\out{
}}
+\if{html}{\out{}}
+\if{latex}{\out{\hypertarget{method-distribution.continuous.lognormal.class-r}{}}}
+\subsection{Method \code{r()}}{
+Generates random deviates for a lognormal random variable with
+mean log(X) \verb{$params$meanlog} and standard deviation logx(X) \verb{$params$sdlog} .
+\subsection{Usage}{
+\if{html}{\out{}}\preformatted{distribution.continuous.lognormal.class$r(n)}\if{html}{\out{
}}
+}
+
+\subsection{Arguments}{
+\if{html}{\out{}}
+\describe{
+\item{\code{n}}{number of observations. If \code{length( n ) > 1}, the length is
+taken to be the number required.}
+}
+\if{html}{\out{
}}
+}
+}
+\if{html}{\out{
}}
+\if{html}{\out{}}
+\if{latex}{\out{\hypertarget{method-distribution.continuous.lognormal.class-clone}{}}}
+\subsection{Method \code{clone()}}{
+The objects of this class are cloneable with this method.
+\subsection{Usage}{
+\if{html}{\out{}}\preformatted{distribution.continuous.lognormal.class$clone(deep = FALSE)}\if{html}{\out{
}}
+}
+
+\subsection{Arguments}{
+\if{html}{\out{}}
+\describe{
+\item{\code{deep}}{Whether to make a deep clone.}
+}
+\if{html}{\out{
}}
+}
+}
+}
diff --git a/man/distribution.lognormal.Rd b/man/distribution.lognormal.Rd
new file mode 100644
index 0000000..6514767
--- /dev/null
+++ b/man/distribution.lognormal.Rd
@@ -0,0 +1,22 @@
+% Generated by roxygen2: do not edit by hand
+% Please edit documentation in R/distribution_continuous.R
+\name{distribution.lognormal}
+\alias{distribution.lognormal}
+\title{distribution.lognormal}
+\usage{
+distribution.lognormal(meanlog, sdlog)
+}
+\arguments{
+\item{meanlog}{The mean of log(X)}
+
+\item{sdlog}{The standard deviation log(X)}
+}
+\value{
+An object of class [\link{distribution.continuous.lognormal.class}]
+}
+\description{
+Constructor function for an object of class \code{distribution.continuous.lognormal.class}
+}
+\seealso{
+\link{Mastiff-Distributions}
+}
diff --git a/man/distribution.mixture.class.Rd b/man/distribution.mixture.class.Rd
index f703920..b6bc0dd 100644
--- a/man/distribution.mixture.class.Rd
+++ b/man/distribution.mixture.class.Rd
@@ -9,6 +9,19 @@ Class to describe the mixture of distributions
\section{Super classes}{
\code{mastiff::R6.class.class} -> \code{\link[mastiff:distribution.abstract.class]{mastiff::distribution.abstract.class}} -> \code{\link[mastiff:distribution.continuous.class]{mastiff::distribution.continuous.class}} -> \code{distribution.mixture.class}
}
+\section{Active bindings}{
+\if{html}{\out{}}
+\describe{
+\item{\code{n_distributions}}{the number of distributions in the mixture}
+
+\item{\code{distributions}}{the distributions in the mixture}
+
+\item{\code{weights}}{the weights of the distributions in the mixture}
+
+\item{\code{interfaces}}{the list of available class interfaces}
+}
+\if{html}{\out{
}}
+}
\section{Methods}{
\subsection{Public methods}{
\itemize{
@@ -29,6 +42,15 @@ Create a new object of class \code{distribution.mixture.class}
\if{html}{\out{}}\preformatted{distribution.mixture.class$new(distributions, weights)}\if{html}{\out{
}}
}
+\subsection{Arguments}{
+\if{html}{\out{}}
+\describe{
+\item{\code{distributions}}{list of distributions in the mixture}
+
+\item{\code{weights}}{vector of weights of the distributions in the mixture}
+}
+\if{html}{\out{
}}
+}
}
\if{html}{\out{
}}
\if{html}{\out{}}
@@ -39,6 +61,15 @@ Density function for a random variable of the mixture
\if{html}{\out{}}\preformatted{distribution.mixture.class$d(x, log = FALSE)}\if{html}{\out{
}}
}
+\subsection{Arguments}{
+\if{html}{\out{}}
+\describe{
+\item{\code{x}}{vector of quantiles.}
+
+\item{\code{log}}{logical; if TRUE, probabilities p are given as \code{log(p)}.}
+}
+\if{html}{\out{
}}
+}
}
\if{html}{\out{
}}
\if{html}{\out{}}
@@ -93,6 +124,14 @@ Generates random samples of the mixture
\if{html}{\out{}}\preformatted{distribution.mixture.class$r(n)}\if{html}{\out{
}}
}
+\subsection{Arguments}{
+\if{html}{\out{}}
+\describe{
+\item{\code{n}}{number of observations. If \code{length( n ) > 1}, the length is
+taken to be the number required.}
+}
+\if{html}{\out{
}}
+}
}
\if{html}{\out{
}}
\if{html}{\out{}}
diff --git a/tests/testthat/test-distribution_continuous.R b/tests/testthat/test-distribution_continuous.R
index 26754fd..e548b2f 100644
--- a/tests/testthat/test-distribution_continuous.R
+++ b/tests/testthat/test-distribution_continuous.R
@@ -246,3 +246,42 @@ test_that( "distribution.normal constructs a valid class", {
foo = 1 ) )
})
})
+
+test_that( "distribution.lognormal constructs a valid class", {
+ withr::with_seed( 123, {
+ n <- 1e5
+ mu <- 0.3
+ sigma <- 1.2
+ tol <- 3 / sqrt( n )
+
+ X <- distribution.lognormal( meanlog = mu, sdlog = sigma )
+
+ # Test that density is correct for initial rate parameter
+ x <- 1.2
+ expect_equal( X$d( x = x ), exp( -(log(x)-mu)^2 / 2 / sigma^2) / sqrt( 2 * pi ) /sigma / x )
+ expect_equal( mean( X$r( n ) ), X$mean, tolerance = tol * X$mean )
+ expect_equal( sd( X$r( n ) ), X$sd, tolerance = tol * X$sd )
+
+ # Test that $params can be updated via named list
+ expect_no_error( X$params <- list( meanlog = 1.5, sdlog = 3 ) )
+ expect_equal( X$params$meanlog, 1.5 )
+ expect_equal( X$params$sdlog, 3 )
+
+ # Test that elements of $params can be updated by name
+ expect_no_error( X$params$meanlog <- 2 )
+ expect_no_error( X$params$sdlog <- 10 )
+ expect_equal( X$params$meanlog, 2 )
+ expect_equal( X$params$sdlog, 10 )
+
+ # Test that invalid values of $params fail (via private$.check_params())
+ expect_error( X$params$mean <- 'a' )
+ expect_error( X$params$sd <- -1 )
+ expect_error( X$params$sd <- 'a' )
+
+ # Test that incorrectly named list $params is rejected
+ expect_error( X$params <- list( foo = 1 ) )
+ expect_error( X$params <- list( 1 ) )
+ expect_error( X$params <- list( meanlog = 1, foo = 1 ) )
+ })
+})
+
diff --git a/tests/testthat/test-distribution_mixutre.R b/tests/testthat/test-distribution_mixutre.R
index 618063a..19ab307 100644
--- a/tests/testthat/test-distribution_mixutre.R
+++ b/tests/testthat/test-distribution_mixutre.R
@@ -50,11 +50,21 @@ test_that( "Test mixture distribution", {
expect_equal( mix$d( c(1,2) ), mix$distributions[[idx]]$d( c(1,2) ) )
expect_equal( mix$p( c(0.25,0.75) ), mix$distributions[[idx]]$p( c(0.25,0.75) ) )
expect_equal( mix$q( c(0.25,0.75) ), mix$distributions[[idx]]$q( c(0.25,0.75) ) )
+
+ # add log version
+ expect_equal( mix$d( c(1,2), log = T ), mix$distributions[[idx]]$d( c(1,2), log = T ) )
+ expect_equal( mix$p( c(1,2), log = T ), mix$distributions[[idx]]$p( c(1,2), log = T ) )
}
- # check the distribution of random draws from mixture agree with CDF
+ # check log version for meaningful mixtures
weights <- seq( 1:mix$n_distributions )
mix$weights <- weights / sum( weights )
+ dd <- unlist( lapply( 1:mix$n_distributions, function( idx ) mix$distributions[[idx]]$d( 1.1 ) ) )
+ expect_equal( mix$d( 1.1, log = T ), log( sum( dd * mix$weights ) ) )
+ pd <- unlist( lapply( 1:mix$n_distributions, function( idx ) mix$distributions[[idx]]$p( 1.1 ) ) )
+ expect_equal( mix$p( 1.1, log = T ), log( sum( pd * mix$weights ) ) )
+
+ # check the distribution of random draws from mixture agree with CDF
n_samples <- 1e4
xs <- mix$r( n_samples )