chjackson / flexsurv

The flexsurv R package for flexible parametric survival and multi-state modelling
http://chjackson.github.io/flexsurv/
53 stars 28 forks source link

Extend the ability of hgengamma with simple change to how hgengamma is calculated #94

Closed swihart closed 2 years ago

swihart commented 3 years ago

If we calculate the log scale density and survival for generalized gamma, then exponentiate, we can avoid NaN for a wider range of mu values.

library(flexsurv)
#hgengamma <- function (x, mu = 0, sigma = 1, Q) 
#{
#  dgengamma(x = x, mu = mu, sigma = sigma, Q = Q)/pgengamma(q = x, 
#                                                            mu = mu, sigma = sigma, Q = Q, #lower.tail = FALSE)
#}

hgengamma2 <- function (x, mu = 0, sigma = 1, Q) 
{
  logdens <- dgengamma(x = x, mu = mu, sigma = sigma, Q = Q, log=TRUE)
  logsurv <- pgengamma(q = x, mu = mu, sigma = sigma, Q = Q, lower.tail = FALSE, log.p=TRUE)

  loghaz <- logdens - logsurv
  haz <- exp(loghaz)
  haz
}

## we get same answer for moderate mu
hgengamma( 10, mu= -1, sigma=1.5, Q=1) 
hgengamma2(10, mu= -1, sigma=1.5, Q=1) 

## we are able to calculate answers for more extreme mu
hgengamma( 10, mu= -14, sigma=1.5, Q=1) ## gives NaN
hgengamma2(10, mu= -14, sigma=1.5, Q=1) ## gives Answer
chjackson commented 3 years ago

Thanks - good suggestion. Just implemented it for this and several related distributions in https://github.com/chjackson/flexsurv-dev/commit/9d1160351759c9feafbbe3bef62a0c654f7c1c29 .