TreeAge Pro 2016 implementation of the Gompertz distribution restricts both rate (lambda) and shape (gamma) parameters to be strictly positive. Some statistical survival packages (such as Flexsurv in R) have relaxed the restriction for the shape parameter and it is possible to obtain negative values of shape from survival studies.
When sampling from Gompertz distribution with a negative shape there is a probability of receiving infinite values. The probability is dependent on the actual values of shape and rate parameters. Specifically the probability of sampling infinity is equal to exp(rate / shape).
Here is an example from R flexsurv:
> y <- rgompertz(100000,-0.02369,0.08501)
> max(y)
[1] Inf
> mean(y)
[1] Inf
> min(y)
[1] 7.819454e-06
We can remove the infinite values from the y set as follows:
> x <- y[is.finite(y)]
> length(x)
[1] 97270
> max(x)
[1] 420.1374
> mean(x)
[1] 15.39763
> sqrt(var(x))
[1] 20.89373
Note that the maximum value of the positive samples is 420.1374. If the distribution is supposed to represent age at death, that would obviously be a problematic value. Regression analysis with limited time window may generate such parameter estimates. The distribution fitted into that time window may have a very good fit, but the tail of the fitted distribution extrapolates into the future with potentially unreasonable values. It is important to examine the behavior of the distributions to make sure that only reasonable values are used.
It is possible to use custom function within TreeAge Pro to generate samples that correspond to Gompertz distribution with negative shape parameter and make sure that only finite (but still potentially large values) will be returned. Refer to the CustomGompertz.trex model which implements custom Gompertz sampling. A uniform distribution with lower bound parameter equal to exp(rate / shape) is used for negative shape. For positive shape values both built in and custom Gompertz will return very similar results. When shape is negative the samples will be similar to the R example above (after the infinities where removed). Please note that it is possible to obtain all infinities with some extreme parameter values:
> y <- rgompertz(100000,-20.00,0.0000005)
> min(y)
[1] Inf
Using such values of shape and rate within custom Gompertz distribution will result in an error.
Additionally, if truncated values are desired, such as in the example above truncating the samples to maximum of e.g. 110.0. It is possible to do that with a min(CustomGompertz; 110) statement as shown in
Comments
Please sign in to leave a comment.