triangular Distribution: Implementation for R

triangular.R

### https://de.wikipedia.org/wiki/Dreiecksverteilung
dtriangular <- function(x,a,c,b){
  
  if(a>=b){stop("a must be smaller than b")}
  if(a>=c){stop("a must be smaller than c")}
  if(c>=b){stop("c must be smaller than b")}
  if(length(x) >1){
    d <-  NULL
    for(i in 1:length(x)){
      d<- c(d, dtriangular(x[i],a,c,b))
    }
    return(d)
  }else{
    if(x >= a && x < c){ 
      return ((2*(x-a)/((c-a)*(b-a))))
    }else if( x == c){
      return (2 / (b-a))
    }else if( x > c && x<= b){
      return ((2*(b-x)/((b-a)*(b-c))))
    }else{
      return (0)
    }
  }
}
ptriangular <- function(x,a,c,b){
  
  if(a>=b){stop("a must be smaller than b")}
  if(a>=c){stop("a must be smaller than c")}
  if(c>=b){stop("c must be smaller than b")}
  if(length(x) >1){
    p <-  NULL
    for(i in 1:length(x)){
      p<- c(p, ptriangular(x[i],a,c,b))
    }
    return(p)
  }else{
    if(x< a){
      return (0)
    }else if(x >= a && x < c){ 
      return ((x-a)*(x-a) / (b-a) / (c-a))
    }else if( x == c){
      return ((c-a) / (b-a))
    }else if( x > c && x<= b){
      return (1-(b-x)*(b-x) / (b-a) / (b-c))
    }else if(x>b){
      return (1)
    }else{
      return (NA)
    }
  }
}

qtriangular <- function(y,a,c,b){
  
  if(a>=b){stop("a must be smaller than b")}
  if(a>=c){stop("a must be smaller than c")}
  if(c>=b){stop("c must be smaller than b")}
  if(length(y) >1){
    q <-  NULL
    for(i in 1:length(y)){
      q<- c(q, qtriangular(y[i],a,c,b))
    }
    return(q)
  }else{
    if(y>=0 && y<= (c-a)/(b-a)){ 
      return(a+sqrt(y*(b-a)*(c-a)))
    }else if(y>=(c-a)/(b-a)  && y<=1){ 
      return(b-sqrt((b-a)*(b-c))*sqrt(1-y)) 
    }else{
      return (NA)
    } 
  }
}
rtriangular <- function(n, a,c,b){
  if(a>=b){stop("a must be smaller than b")}
  if(a>=c){stop("a must be smaller than c")}
  if(c>=b){stop("c must be smaller than b")}
   
 
  return( qtriangular(runif(n),a,c,b))
}