Thursday, February 19, 2015

R25. Gamma Distribution in R


Like the Poisson and Exponential distribution, the Gamma distribution is also dependent on the exponential series.


Besides the calculation, we can use the function pgamma.

# ex25.R
# shape, scale parameters
k <- 7
s <- 2

get_num <- function(x,kv) {
  xp <- x/s
  coeff <- exp(-xp)
  if (kv == 0) return(coeff) 
  cv <- xp/1:kv
  coeff*prod(cv)
}

gamma.dis <- function(x, calc) {
  dis.calc <- numeric(length(x))
  for (i in 1:length(x)) {
    S <- mapply(get_num,x[i],0:(k-1))
    dis.calc[i] <- 1 - sum(S)
  }
  dis.R <- pgamma(x, shape = k, scale = s)
  if (calc) return(dis.calc) else return(dis.R)
}

main <- paste('Gamma Cumulative Distribution Function',
              'Shape = 7, Scale = 2', sep = '\n', collapse = '')
x <- seq(0,40,.1)
plot(x, gamma.dis(x,FALSE), type = 'l', col = 'red', 
     ylab = 'P', main = main)
lines(x, gamma.dis(x,TRUE), col = 'blue')

Output:

No comments:

Post a Comment