New distibution - Asymmetric Laplace. Crazy things are comming, which…

… are relevant to issue #13 !

DESCRIPTION

@@ -1,8 +1,8 @@
 Package: greybox
 Type: Package
 Title: Toolbox for Model Building and Forecasting
-Version: 0.3.2.41007
-Date: 2018-09-30
+Version: 0.3.2.41008
+Date: 2018-10-02
 Authors@R: person("Ivan", "Svetunkov", email = "[email protected]", role = c("aut", "cre"),
                   comment="Lecturer at Centre for Marketing Analytics and Forecasting, Lancaster University, UK")
 URL: https://github.com/config-i1/greybox

NAMESPACE

@@ -40,6 +40,7 @@ S3method(print,summary.alm)
 S3method(print,summary.greybox)
 S3method(print,summary.greyboxC)
 S3method(sigma,alm)
+S3method(sigma,ets)
 S3method(sigma,greybox)
 S3method(summary,alm)
 S3method(summary,greybox)
@@ -52,6 +53,7 @@ S3method(vcov,lmGreybox)
 export(AICc)
 export(BICc)
 export(alm)
+export(dalaplace)
 export(determination)
 export(dfnorm)
 export(dlaplace)
@@ -71,13 +73,16 @@ export(nParam)
 export(pAIC)
 export(pAICc)
 export(pBIC)
+export(palaplace)
 export(pfnorm)
 export(plaplace)
 export(pointLik)
 export(ps)
+export(qalaplace)
 export(qfnorm)
 export(qlaplace)
 export(qs)
+export(ralaplace)
 export(rfnorm)
 export(rlaplace)
 export(rmc)

NEWS

@@ -1,4 +1,4 @@
-greybox v0.3.2 (Release data: 2018-09-30)
+greybox v0.3.2 (Release data: 2018-01-02)
 ==============
 
 Changes:
@@ -14,6 +14,7 @@ Changes:
 * Residuals of dnbinom and dpois are now calculated as y - fitted.
 * predict with distriubion=dnbinom in cases, when scale is not available, is now calculated based on the definition of scale via variance and mean.
 * pointLik.ets() is now calculated differently, so that sum(pointLik) is close to the logLik produced by ets() function. The problem with logLik of ets() is that it is not calculated correctly, chopping off some parts of normal distribution. Total disaster!
+* New set of distribution functions - for Asymmetric Laplace Distribution.
 
 Bugfixes:
 * predict.alm() sometimes produced NAs in the lower bound.

R/alaplace.R

@@ -0,0 +1,122 @@
+#' Asymmetric Laplace Distribution
+#'
+#' Density, cumulative distribution, quantile functions and random
+#' generation for the Asymmetric Laplace distribution with the
+#' location parameter mu, Mean Absolute Error (or Mean Absolute
+#' Deviation) equal to 2 b and asymmetry parameter alpha.
+#'
+#' When mu=0 and b=1, the Laplace distribution becomes standardized.
+#' The distribution has the following density function:
+#'
+#' f(x) = alpha (1-alpha) / b exp(-(x-mu)/b (alpha - I(x<=mu))),
+#'
+#' where I(.) is the indicator function (equal to 1 if the condition is
+#' satisfied and zero otherwise).
+#'
+#' When alpha=0.5, then the distribution becomes Symmetric Laplace, where
+#' b = 1/2 MAE.
+#'
+#' This distribution function aligns with the quantile estimates of
+#' parameters (Geraci & Bottai, 2007).
+#'
+#' Finally, both \code{palaplace} and \code{qalaplace} are returned for
+#' the lower tail of the distribution.
+#'
+#' @template author
+#' @keywords distribution
+#'
+#' @param q vector of quantiles.
+#' @param p vector of probabilities.
+#' @param n number of observations. Should be a single number.
+#' @param mu vector of location parameters (means).
+#' @param b vector of scale parameters.
+#' @param alpha value of asymmetry parameter. Varies from 0 to 1.
+#' @param log if \code{TRUE}, then probabilities are returned in
+#' logarithms.
+#'
+#' @return Depending on the function, various things are returned
+#' (usually either vector or scalar):
+#' \itemize{
+#' \item \code{dalaplace} returns the density function value for the
+#' provided parameters.
+#' \item \code{palaplace} returns the value of the cumulative function
+#' for the provided parameters.
+#' \item \code{qalaplace} returns quantiles of the distribution. Depending
+#' on what was provided in \code{p}, \code{mu} and \code{b}, this
+#' can be either a vector or a matrix, or an array.
+#' \item \code{ralaplace} returns a vector of random variables
+#' generated from the Laplace distribution. Depending on what was
+#' provided in \code{mu} and \code{b}, this can be either a vector
+#' or a matrix or an array.
+#' }
+#'
+#' @examples
+#' x <- dalaplace(c(-100:100)/10, 0, 1, 0.2)
+#' plot(x, type="l")
+#'
+#' x <- palaplace(c(-100:100)/10, 0, 1, 0.2)
+#' plot(x, type="l")
+#'
+#' qalaplace(c(0.025,0.975), 0, c(1,2), c(0.2,0.3))
+#'
+#' x <- ralaplace(1000, 0, 1, 0.2)
+#' hist(x)
+#'
+#' @references \itemize{
+#' \item Geraci Marco, Bottai Matteo (2007). Quantile regression for
+#' longitudinal data using the asymmetric Laplace distribution.
+#' Biostatistics (2007), 8, 1, pp. 140-154
+#' \url{https://doi.org/10.1093/biostatistics/kxj039}
+#' \item Yu, K., & Zhang, J. (2005). A three-parameter asymmetric
+#' laplace distribution and its extension. Communications in Statistics
+#' - Theory and Methods, 34, 1867-1879.
+#' \url{https://doi.org/10.1080/03610920500199018}
+#' }
+#'
+#' @rdname alaplace-distribution
+
+#' @rdname alaplace-distribution
+#' @export dalaplace
+#' @aliases dalaplace
+dalaplace <- function(q, mu=0, b=1, alpha=0.5, log=FALSE){
+    e <- q-mu
+    alaplaceReturn <- alpha * (1-alpha) / b * exp(-(e)/b * (alpha - ((e) <= 0)*1));
+    if(log){
+        alaplaceReturn <- log(alaplaceReturn);
+    }
+    return(alaplaceReturn);
+}
+
+#' @rdname alaplace-distribution
+#' @export palaplace
+#' @aliases palaplace
+palaplace <- function(q, mu=0, b=1, alpha=0.5){
+    Ie <- (q<=mu)*1;
+    alaplaceReturn <- 1 - Ie - (1 - Ie - alpha) * exp((Ie - alpha) / b * (q-mu));
+
+    return(alaplaceReturn);
+}
+
+#' @rdname alaplace-distribution
+#' @export qalaplace
+#' @aliases qalaplace
+qalaplace <- function(p, mu=0, b=1, alpha=0.5){
+    Ie <- (p<=alpha)*1;
+    alaplaceReturn <- mu + b / (Ie - alpha) * log((1-Ie-p)/(1-Ie-alpha));
+    if(any(p==0)){
+        alaplaceReturn[p==0] <- -Inf;
+    }
+    if(any(p==1)){
+        alaplaceReturn[p==1] <- Inf;
+    }
+
+    return(alaplaceReturn);
+}
+
+#' @rdname alaplace-distribution
+#' @export ralaplace
+#' @aliases ralaplace
+ralaplace <- function(n=1, mu=0, b=1, alpha=0.5){
+    alaplaceReturn <- qalaplace(runif(n,0,1),mu,b,alpha);
+    return(alaplaceReturn);
+}

R/alm.R

@@ -56,6 +56,8 @@
 #' @param vcovProduce whether to produce variance-covariance matrix of
 #' coefficients or not. This is done via hessian calculation, so might be
 #' computationally costly.
+#' @param ... additional parameters to pass to distribution functions
+#' (e.g. \code{alpha} value for Asymmetric Laplace distribution).
 #'
 #' @return Function returns \code{model} - the final model of the class
 #' "alm", which contains:
@@ -119,7 +121,7 @@ alm <- function(formula, data, subset, na.action,
                                "dpois","dnbinom",
                                "plogis","pnorm"),
                 occurrence=c("none","plogis","pnorm"),
-                B=NULL, vcovProduce=FALSE){
+                B=NULL, vcovProduce=FALSE, ...){
 # Useful stuff for dnbinom: https://scialert.net/fulltext/?doi=ajms.2010.1.15
 
     cl <- match.call();

R/fnorm.R

@@ -9,6 +9,9 @@
 #'
 #' f(x) = 1/sqrt(2 pi) (exp(-(x-mu)^2 / (2 sigma^2)) + exp(-(x+mu)^2 / (2 sigma^2)))
 #'
+#' Both \code{pfnorm} and \code{qfnorm} are returned for the lower
+#' tail of the distribution.
+#'
 #' @template author
 #' @keywords distribution
 #'

R/laplace.R

@@ -4,11 +4,14 @@
 #' generation for the Laplace distribution with the location parameter mu
 #' and Mean Absolute Error (or Mean Absolute Deviation) equal to b.
 #'
-#' When mu=0 and b=1, the Laplace distribution becomes standardized with.
+#' When mu=0 and b=1, the Laplace distribution becomes standardized.
 #' The distribution has the following density function:
 #'
 #' f(x) = 1/(2b) exp(-abs(x-mu) / b)
 #'
+#' Both \code{plaplace} and \code{qlaplace} are returned for the lower
+#' tail of the distribution.
+#'
 #' @template author
 #' @keywords distribution
 #'
@@ -37,10 +40,10 @@
 #' }
 #'
 #' @examples
-#' x <- dlaplace(c(-1000:1000)/10, 0, 1)
+#' x <- dlaplace(c(-100:100)/10, 0, 1)
 #' plot(x, type="l")
 #'
-#' x <- plaplace(c(-1000:1000)/10, 0, 1)
+#' x <- plaplace(c(-100:100)/10, 0, 1)
 #' plot(x, type="l")
 #'
 #' qlaplace(c(0.025,0.975), 0, c(1,2))

R/svet.R

@@ -12,14 +12,17 @@
 #'
 #' The S distribution has fat tails and large excess.
 #'
+#' Both \code{ps} and \code{qs} are returned for the lower tail of
+#' the distribution.
+#'
 #' @template author
 #' @keywords distribution
 #'
 #' @param q vector of quantiles.
 #' @param p vector of probabilities.
 #' @param n number of observations. Should be a single number.
 #' @param mu vector of location parameters (means).
-#' @param b vector of scaling parameter (which equals to ham/2).
+#' @param b vector of scaling parameter (which are equal to ham/2).
 #' @param log if \code{TRUE}, then probabilities are returned in
 #' logarithms.
 #'

README.md

@@ -15,16 +15,18 @@ The list of included functions:
 7. lmDynamic - linear regression with time varying parameters based on pAIC.
 8. determination - the function returns the vector of coefficients of determination (R^2) for the provided data. This is useful for the diagnostics of multicollinearity.
 9. qlaplace, dlaplace, rlaplace, plaplace - functions for Laplace distribution.
-10. qs, ds, rs, ps - functions for S distribution.
-11. qfonrm, dfnorm, rfnorm, pfnorm - functions for folded normal distribution.
-12. graphmaker - plots the original series, the fitted values and the forecasts.
+10. qalaplace, dalaplace, ralaplace, palaplace - functions for Asymmetric Laplace distribution.
+11. qs, ds, rs, ps - functions for S distribution.
+12. qfonrm, dfnorm, rfnorm, pfnorm - functions for folded normal distribution.
+13. graphmaker - plots the original series, the fitted values and the forecasts.
 
 Future functions:
 1. nonlinearExpander - Function produces non-linear transformations of the provided data.
 
 Methods already implemented:
 1. pointLik - point likelihood method for the time series models.
 2. pAIC - point AIC based on pointLik.
+3. pAICc, pBIC, pBICc - respective point values for the information criteria.
 3. summary - returns summary of the regression (either selected or combined).
 4. vcov - covariance matrix for combined models. This is an approximate thing. The real one is quite messy and not yet available.
 5. confint - confidence intervals for combined models.