# rinvgamma: sample a vector of inverse-gamma random variables
rinvgamma = function(n, shape = 1, rate = 1, scale = 1/rate) {
    return( 1 / rgamma(n,shape,rate,scale) )
}

# dinvgamma: density function of the inverse-gamma distribution
dinvgamma = function(x, shape = 1, rate = 1, scale = 1/rate) {
    return( shape^rate * gamma(rate) * exp(-shape/x) * x^(-rate-1) )
}