Additional Technical Details

Bias and CI coverage and interval width

AICc based model weights

Aside from considerations of bias, coverage and confidence interval width, it is also prudent to examine how the weights of the different distributions changed with sample size, for data generated using the six different default distributions, to more fully investigate sample size issues associated with the use of the mixture distributions. This was examined using the simulation study across the larger combined ssddata and EnviroTox datasets (353 datasets) to ensure a wide range of potential representations of each of the six default distributions was considered.

Below is a plot of the mean AICc weights as a function of sample size (N) obtained for data simulated using the best fit distribution to 353 datasets. Results are shown separately for the six different source distributions, with the upper plot (A) showing the AICc weight of the source (data generating) distributions, when fit using the default set of six distributions using ssdtools; and the lower plot (B) showing the AICc weight of the lognormal-lognormal mixture distribution for each of the source (data generating) distributions.

We found that the AICc weights for the five unimodal distributions were relatively similar (~0.2) for very low N. This is because for small N it is difficult to discern differences between the distributions in the candidate list. The weights increase to above 0.5 as N increases (i.e. their converge to the true underlying generating distribution at high N, upper plot). For very low sample sizes (N=5, 6 or 7) the source (data-generating) distribution is not preferentially weighted by the AICc (upper plot, A).

For the lnorm_lnorm mixture, AICc weights can be >0.5 at relatively N (>8, lower plot, B). We also looked specifically at the AICc weight of the lognormal-lognormal mixture as one of six distributions in the default candidate set, across simulation based on all six source generating distributions. This was done to examine the potential for erroneously highly weighting the lognormal-lognormal mixture distribution by chance, when data are generated instead using one of the five unimodal distributions. For all the unimodal source distributions the lognormal-lognormal mixture never has high AICc weight, even at very high sample sizes (N=256, lower plot, B). The AICc weights are particularly low for the lnorm_lnorm mixture at low sample sizes for all the uni-modal source distributions (a desirable property) (lower plot, B).

Conclusions

Overall, the results suggest that (relative) bias as a function of N behaves as expected. The N=6 recommendation appears to be well-supported as coverage is particularly low for N=5, but acceptable for N=6. The lognormal-lognormal mixture AICc weights suggest that this distribution will only be preferred (i.e. have a high AICc weight) when (i) there is clear evidence of bimodality in the source data; and (ii) large N.

Based on these results, our recommendation is that only a single minimum sample size of N=6 be adopted (for both unimodal and bimodal), since our results suggest any gains in increasing this to 7 are minimal.

Burr III distribution

The probability density function, \({f_X}(x;b,c,k)\) and cumulative distribution function, \({F_X}(x;b,c,k)\) for the Burr III distribution (also known as the Dagum distribution) as used in ssdtools are:

Inverse Pareto distribution

Let \(X \sim Burr(b,c,k)\) have the pdf given in the box above. It is well known that the distribution of \(Y = \frac{1}{X}\) is the inverse Burr distribution (also known as the SinghMaddala distribution) for which:\[\begin{array}{*{20}{c}} {{f_Y}(y;b,c,k) = \frac{{c{\kern 1pt} {\kern 1pt} k{{\left( {\frac{y}{b}} \right)}^c}}}{{y{\kern 1pt} {{\left[ {1 + {{\left( {\frac{y}{b}} \right)}^c}} \right]}^{k + 1}}}}}&{b,c,k,y > 0} \end{array}\]

\[\begin{array}{*{20}{c}} {{F_Y}(y;b,c,k) = 1 - \frac{1}{{{{\left[ {1 + {{\left( {\frac{y}{b}} \right)}^c}} \right]}^k}}}}&{b,c,k,y > 0} \end{array}\]

We now consider the limiting distribution when \(c \to \infty\) and \(k \to 0\) in such a way that the product \(ck\) remains constant, i.e. \(ck = \lambda\).

Now, \[\begin{array}{l} \mathop {\mathop {\lim }\limits_{(c,k) \to (\infty ,0)} }\limits_{ck = \lambda } \left\{ {{F_Y}(y;b,c,k)} \right\} = 1 - \mathop {\mathop {\lim }\limits_{(c,k) \to (\infty ,0)} }\limits_{ck = \lambda } \frac{1}{{{{\left[ {1 + {{\left( {\frac{y}{b}} \right)}^c}} \right]}^k}}}\\ \\ and\\ \\ \mathop {\mathop {\lim }\limits_{(c,k) \to (\infty ,0)} }\limits_{ck = \lambda } {\left[ {1 + {{\left( {\frac{y}{b}} \right)}^c}} \right]^k} = \mathop {\mathop {\lim }\limits_{(c,k) \to (\infty ,0)} }\limits_{ck = \lambda } \left\{ {{{\left( {\frac{y}{b}} \right)}^{ck}}{{\left[ {1 + {{\left( {\frac{b}{y}} \right)}^c}} \right]}^k}} \right\}\\ \\ and\\ \\ \mathop {\mathop {\lim }\limits_{(c,k) \to (\infty ,0)} }\limits_{ck = \lambda } \left\{ {{{\left( {\frac{y}{b}} \right)}^{ck}}{{\left[ {1 + {{\left( {\frac{b}{y}} \right)}^c}} \right]}^k}} \right\} = \mathop {\mathop {\lim }\limits_{(c,k) \to (\infty ,0)} }\limits_{ck = \lambda } \left\{ {{{\left( {\frac{y}{b}} \right)}^{ck}}} \right\}\mathop {\mathop {\lim }\limits_{(c,k) \to (\infty ,0)} }\limits_{ck = \lambda } \left\{ {{{\left[ {1 + {{\left( {\frac{b}{y}} \right)}^c}} \right]}^k}} \right\}\\ = \mathop {\mathop {\lim }\limits_{(c,k) \to (\infty ,0)} }\limits_{ck = \lambda } \left\{ {{{\left( {\frac{y}{b}} \right)}^{ck}}} \right\}\; \cdot \,1\\ = {\left( {\frac{y}{b}} \right)^\lambda } \end{array}\]

Therefore, \[\begin{array}{*{20}{c}} {\mathop {\mathop {\lim }\limits_{(c,k) \to (\infty ,0)} }\limits_{ck = \lambda } \left\{ {{F_Y}(y;b,c,k)} \right\} = 1 - {{\left( {\frac{b}{y}} \right)}^\lambda }}&{y \ge b} \end{array}\]

which we recognise as the (American) Pareto distribution. So, if the limiting distribution of \(Y = \frac{1}{X}\) is a Pareto distribution, then the limiting distribution of \(X = \frac{1}{Y}\) is the (American) inverse Pareto distribution:

\[\begin{array}{l} {f_X}\left( {x;\alpha ,\beta } \right) = \lambda {b^\lambda }{x^{\lambda - 1}};{\rm{ }}0 \le x \le {\textstyle{1 \over b}};{\rm{ }}\lambda {\rm{,}}b > 0\\ {F_X}\left( {x;\alpha ,\beta } \right) = {\left( {xb} \right)^\lambda };{\rm{ }}0 \le x \le {\textstyle{1 \over b}};{\rm{ }}\lambda {\rm{,}}b > 0 \end{array}\]

For completeness, the MLEs of this distribution have closed-form expressions and are given by: \[\begin{array}{l} \hat \lambda = {\left[ {\ln \left( {\frac{{{g_X}}}{{\hat b}}} \right)} \right]^{ - 1}}\\ \hat b = \frac{1}{{\max \left\{ {{X_i}} \right\}}}{\rm{ }} \end{array}\]

and \({\rm{ }}{g_X}\)is the geometric mean of the data.

Inverse Weibull distribution

Let \(X \sim Burr(b,c,k)\) have the pdf given in the box above. We make the transformation \[Y = \frac{{b{\kern 1pt} {k^{\tfrac{1}{c}}}{\kern 1pt} \theta }}{X}\] where \(\theta\) is a parameter (constant). The distribution of \(Y\) is also a Burr distribution and has cdf \[{G_Y}\left( y \right) = 1 - \frac{1}{{{{\left[ {1 + {{\left( {\frac{y}{{{k^{\tfrac{1}{c}}}{\kern 1pt} \theta }}} \right)}^c}} \right]}^k}}}\].

We are interested in the limiting behaviour of this Burr distribution as \(k \to \infty\).
Now,\[\mathop {\lim }\limits_{k \to \infty } {G_Y}\left( y \right) = 1 - \mathop {\lim }\limits_{k \to \infty } {\left[ {1 + {{\left( {\frac{y}{{{k^{\tfrac{1}{c}}}{\kern 1pt} \theta }}} \right)}^c}} \right]^{ - k}}\]

\[{ = 1 - \mathop {\lim }\limits_{k \to \infty } {{\left[ {1 + \frac{{{{\left( {\frac{y}{\theta }} \right)}^c}}}{{k{\kern 1pt} }}} \right]}^{ - k}}}\]

\[\begin{matrix} =1-\exp \left[ -{{\left( \frac{y}{\theta } \right)}^{c}} \right] \\ \left\{ \text{using the fact that }\underset{n\to \infty }{\mathop{\lim }}\,{{\left( 1+{}^{z}\!\!\diagup\!\!{}_{n}\; \right)}^{-n}}={{e}^{-z}} \right\} \\ \end{matrix}\]

We recognise the last expression as the cdf of a Weibull distribution with parameters \(c\) and \(\theta\).