Exponentiated Weibull Distribution I'm trying to understand Exponentiated Weibull (EW) function ExpoWeibull from package reliaR , however the function only use a single shape parameter and a scale parameter. From what I've read about the EW distribution, a second shape parameter is used. Is there an existing implementation of EW for R, that accepts two shape parameters? If not, how does one go about implementing the distribution function?  
 A: I disagree with Brad. The function computed in the reliaR package is correct. It's just the same formula like on Wikipedia. You'll just have to substitute (x/lambda)=x1, k=alpha1 and Alpha=Theta. Just take the cumulative distribution, substitute the 3 variables and derive it! Don't Forget to Substitute dx too! Or have a look at http://www.academia.edu/6178638/Estimation_for_the_Parameters_of_the_Exponentiated_Weibull_Distribution_Based_on_Progressive_Hybrid_Censored_Samples. Page 1714 There you can see: both cumulative function are equal
A: As far as I can tell, the function in that package is wrong. I had to code the distribution myself. If lambda is the scale parameter and x and alpha are the shape parameters (going by the notation in Wikipedia for simplicity):
dexpweib=function(x, alpha, k, lambda, log=FALSE){
if ((!is.numeric(x)) || (!is.numeric(alpha))|| (!is.numeric(k)) || (!is.numeric(lambda))){ 
    stop("non-numeric argument to mathematical function")}
if((length(alpha)!= 1) || (length(k)!=1) || (length(lambda)!=1)){
    stop("Non-x parameters must be atomic")}
if ((min(x) <= 0) || (alpha <= 0) || (k <= 0) || (lambda <= 0)){ 
    stop("Invalid arguments. Must be > 0 ")}
u <- exp(-(x/lambda)^k)
pdf <- ((alpha*k)/(lambda^k))*(x^(k-1))*u*((1-u)^(alpha-1))
if (log) 
    pdf <- log(pdf)
return(pdf)

}
pexpweib=function(q, alpha, k, lambda, log.p=FALSE){
if ((!is.numeric(q)) || (!is.numeric(alpha))|| (!is.numeric(k)) || (!is.numeric(lambda))){ 
    stop("non-numeric argument to mathematical function")}
if((length(alpha)!= 1) || (length(k)!=1) || (length(lambda)!=1)){
    stop("Non-q parameters must be atomic")}
if ((min(q) <= 0) || (alpha <= 0) || (k <= 0) || (lambda <= 0)){ 
    stop("Invalid arguments. Must be > 0 ")}
u <- exp(-(q/lambda)^k)
cdf <- (1-u)^alpha
if (log.p) 
    cdf <- log(cdf)
return(cdf)

}
