version 0.6.3

cran · Oct 20, 2020 · 6398b97 · 6398b97
1 parent c05981a
commit 6398b97
Show file tree

Hide file tree

Showing 21 changed files with 745 additions and 594 deletions.
diff --git a/DESCRIPTION b/DESCRIPTION
@@ -1,8 +1,8 @@
 Package: greybox
 Type: Package
 Title: Toolbox for Model Building and Forecasting
-Version: 0.6.2
-Date: 2020-09-02
+Version: 0.6.3
+Date: 2020-10-20
 Authors@R: c(person("Ivan", "Svetunkov", email = "ivan@svetunkov.ru", role = c("aut", "cre"),
                   comment="Lecturer at Centre for Marketing Analytics and Forecasting, Lancaster University, UK"),
              person("Yves R.", "Sagaert", role = c("ctb"),
@@ -21,7 +21,7 @@ Description: Implements functions and instruments for regression model building
 License: GPL (>= 2)
 Depends: R (>= 3.0.2)
 Imports: forecast, stats, graphics, utils, lamW, pracma, nloptr,
-        statmod, zoo, Matrix
+        statmod, zoo, Matrix, texreg, methods
 LinkingTo: Rcpp
 Suggests: smooth (>= 2.5.1), doMC, doParallel, foreach, testthat,
         rmarkdown, knitr
@@ -30,11 +30,11 @@ RoxygenNote: 7.1.1
 VignetteBuilder: knitr
 Encoding: UTF-8
 NeedsCompilation: yes
-Packaged: 2020-09-02 12:56:49 UTC; config
+Packaged: 2020-10-19 23:39:08 UTC; config
 Author: Ivan Svetunkov [aut, cre] (Lecturer at Centre for Marketing Analytics
     and Forecasting, Lancaster University, UK),
   Yves R. Sagaert [ctb] (Visiting Research at Centre for Marketing
     Analytics and Forecasting, Lancaster University, UK)
 Maintainer: Ivan Svetunkov <ivan@svetunkov.ru>
 Repository: CRAN
-Date/Publication: 2020-09-02 15:00:03 UTC
+Date/Publication: 2020-10-20 07:00:16 UTC
diff --git a/MD5 b/MD5
@@ -1,14 +1,14 @@
-582aeec9b523ff70680ced0e4ab00d2f *DESCRIPTION
-9dc185e419865f41105cc8d1ef347e58 *NAMESPACE
-5fe5191a7a91ecab0283f73b64625766 *NEWS
+1c368f54bacfbae8cbaa56178eb95499 *DESCRIPTION
+602cf569f03eae787802caa6319322d8 *NAMESPACE
+7b8316e0628654a2ee9a0684f684ad96 *NEWS
 b04f00196442fad912053f13321c737f *R/RcppExports.R
 7853380e821e94f7bac2096abb087cf0 *R/alaplace.R
-82f055e8f202bad5e8884dbc2ae55b3b *R/alm.R
+a760c1f2cb2d88300a197bc632e20730 *R/alm.R
 6b248b2677385d382f13544174737548 *R/association.R
 1fa4c10a4d87a433c6b56108164b21cb *R/bcnorm.R
 3d36022fbac3500f1bc7b5f7e26877d5 *R/coefbootstrap.R
 9162a49ecbca3900f5da3222c67c35ac *R/cramer.R
-cef343c41e2c9d1b6457b427dad41e32 *R/detector-dst-leap.R
+f7628d5e44724a2ed7a711addf8e06a5 *R/detector-dst-leap.R
 7fda66df665392fafebfc5fd2d0a65ab *R/determination.R
 7e5c0ce8f368d66c5305cba104a12460 *R/error-measures.R
 3d3deabe68ca8fb0f85deb754bb8091c *R/fnorm.R
@@ -21,7 +21,7 @@ cad33cb70b10fff7ba3d4540f75f8f11 *R/intermittentDemandIdentificator.R
 e7f46bc8f90bc0d86c607b76f895ff1c *R/lmCombine.R
 98b4abf47a817d0ab6e8a12e22d36053 *R/lmDynamic.R
 c98573264afdfdffc11e9071fe8c7751 *R/mcor.R
-5be084865535f8614166e41e001f8620 *R/methods.R
+b4cf6c6df0c94c74c82c9de517569720 *R/methods.R
 4eb0114d8e61620b2764a57598044d98 *R/nemenyi.R
 873bc5653f241231d5a40dc16ff5c2c0 *R/pcor.R
 8dacc087c0c87e54764518c96a9482bb *R/rmcb.R
@@ -32,35 +32,35 @@ fe30efab65007ea10c851709f6ad5c1c *R/svet.R
 4535729d68379caba1a579c1ea815e75 *R/tableplot.R
 f8f8d9a646a104835c1b4b84d907ea2f *R/temporaldummy.R
 b326bc5271a425583e7a7a8a8dbddd7a *R/tplnorm.R
-efb8ad4e5bc32000bbeab206a23e5054 *R/xregExpander.R
+44ca03914250f05ef01b28f240e06d11 *R/xregExpander.R
 bdb7ecee051a442d310f4e83e6e9865d *R/xregMultiplier.R
 964ae5b1c6a404c9bcfdb87d7512f2c3 *R/xregTransformer.R
-78765a7c239c85b2c1fcaef6eda0a758 *R/zzz.R
+52e8e62fc2962549cde8a85fb4ca5651 *R/zzz.R
 fd6e0246d580518a6088cc28ddc36615 *README.md
-00f4d8cf07fc1ff72036b846e7b96517 *build/partial.rdb
-9e9158d38e5b8c2211053037242fb57b *build/vignette.rds
+4fc7cd0454a15c88d28192a4318a144b *build/partial.rdb
+81c03b57d5b5aff8aeb802515ea022ae *build/vignette.rds
 d7701fbec190a77b73e572d81ab8402c *inst/doc/alm.R
 8e03b5f7ce2e46bcdb0c91c50b13067d *inst/doc/alm.Rmd
-df13745fb8009d2a9be3afb102274488 *inst/doc/alm.html
-9333725d54a0a462b8ccd87788916e1a *inst/doc/greybox.R
-40e93c83df7e83a8b49199f732d3fe71 *inst/doc/greybox.Rmd
-4a6fed9d98f8524dd7d8a558d728ea77 *inst/doc/greybox.html
+f60d2018be41e99faddbbc1adff88556 *inst/doc/alm.html
+75d5150f03862add1eaf3c26ef9b289a *inst/doc/greybox.R
+23410948046a62a563cf03be05a358f5 *inst/doc/greybox.Rmd
+cd7970bd6c9bbb64817e8c9ec55d8ca6 *inst/doc/greybox.html
 35f3b10303d1942a7d0726228f9fe5d4 *inst/doc/maUsingGreybox.R
 6b5fcb2644ce454ac2bf05431c2e36ef *inst/doc/maUsingGreybox.Rmd
-2b7053e7561160679c7f2f4e92cfedfe *inst/doc/maUsingGreybox.html
+a1a782795dda2fe7614c8a5bea119e4d *inst/doc/maUsingGreybox.html
 ec5a1bacfc67b5db0d6a09ec1934ec76 *inst/doc/ro.R
 af266b3c1ed360928fd6dd365a77ea4a *inst/doc/ro.Rmd
-e264d6cb6604c20ba5d4130adb0d98a0 *inst/doc/ro.html
+291c004543a1d7687e408d767ad7c393 *inst/doc/ro.html
 289e8c5fa923dd43039cfabae03445a0 *man/InformationCriteria.Rd
 e9fce46df49e844416c626c15ebc7628 *man/actuals.Rd
 2da7e8885cca6425fdac21945a64754f *man/alaplace-distribution.Rd
-2898c794e2f529ac5255df8ff293925a *man/alm.Rd
+a10e6e07a77e2f9d2acb29b3d248bc4d *man/alm.Rd
 676fa177135756674cbec993f05b30dd *man/association.Rd
 c460f16bd0eef6d2aec89adfdffad6ed *man/bcnorm-distribution.Rd
 3a3ffa1bb28ba7e0b7791f393a1dd6cf *man/coef.alm.Rd
 5dac8e2865668426cf275eb0c5467988 *man/coefbootstrap.Rd
 dd4df3dd66d4a3aafc593b5b2f3ad79b *man/cramer.Rd
-f3d27164597ee4745b9b0ca3bd2ef1e8 *man/detectdst.Rd
+99aa860c000b423afa7eadb4827f59e8 *man/detectdst.Rd
 c6b538d5f9bff5449a7fe8047bfe596c *man/determination.Rd
 348f87f260981b476b494685016f77bc *man/error-measures.Rd
 7107d89e7d1666697d37f3ffe06b2ff9 *man/errorType.Rd
@@ -91,13 +91,13 @@ c6ebc1c587ffc34c6161c9a0051dbaeb *man/polyprod.Rd
 6e75842d13a488a1f839ac7a63c658e5 *man/tableplot.Rd
 f0fae4d40f58596a7e5c993a07758b8a *man/temporaldummy.Rd
 35386b7aacb453e2d4a06e1a86e11272 *man/tplnorm-distribution.Rd
-09495fff3a561034a51fb5e478340b50 *man/xregExpander.Rd
+3356be96255ff0722e73c961bc7508eb *man/xregExpander.Rd
 ecaff1d893d0e2226b358df3b6e18cd8 *man/xregMultiplier.Rd
 d384247a8281dc3ce210c2825d232aca *man/xregTransformer.Rd
 471a2428a7493a27c2f59c352b429357 *src/RcppExports.cpp
 b0129f5614ed7877c7af10bd9dce93b0 *src/polyMult.cpp
 8e03b5f7ce2e46bcdb0c91c50b13067d *vignettes/alm.Rmd
-40e93c83df7e83a8b49199f732d3fe71 *vignettes/greybox.Rmd
+23410948046a62a563cf03be05a358f5 *vignettes/greybox.Rmd
 8401400e9158cfc020fa9477e42d1fd6 *vignettes/img/ROProcessCO.gif
 da71fe35eb5beb45dc2cc8ef16da371a *vignettes/img/ROProcessCOCI.gif
 7f2a09e79b06d202de45672aa9762d16 *vignettes/img/ROProcessNoCO.gif

diff --git a/NAMESPACE b/NAMESPACE
@@ -6,6 +6,7 @@ S3method(BICc,default)
 S3method(BICc,varest)
 S3method(actuals,alm)
 S3method(actuals,default)
+S3method(as.data.frame,summary.greybox)
 S3method(coef,greybox)
 S3method(coef,greyboxD)
 S3method(coefbootstrap,alm)
@@ -207,6 +208,7 @@ importFrom(graphics,rect)
 importFrom(graphics,text)
 importFrom(graphics,title)
 importFrom(lamW,lambertWm1)
+importFrom(methods,setMethod)
 importFrom(nloptr,nloptr)
 importFrom(pracma,hessian)
 importFrom(statmod,dinvgauss)
@@ -286,6 +288,8 @@ importFrom(stats,time)
 importFrom(stats,ts)
 importFrom(stats,var)
 importFrom(stats,vcov)
+importFrom(texreg,createTexreg)
+importFrom(texreg,extract)
 importFrom(utils,head)
 importFrom(utils,packageVersion)
 importFrom(utils,tail)

diff --git a/NEWS b/NEWS
@@ -1,3 +1,22 @@
+greybox v0.6.3 (Release data: 2020-10-20)
+==============
+
+Changes:
+* Corrections in the main greybox vignette - just fixed some typos and formatting issues.
+* Examples in alm() now rely on subset parameter.
+* xregExpander() now selects between "extrapolate" and "naive" in case of gaps="auto".
+* The as.data.frame.summary.greybox() method, which produces the matrix of coefficients. This should be handy when printing the outputs.
+* We now support texreg package, so that the ALM regression outputs can be printed in a fine quality.
+* detectdst() now complains, when it finds several missing hours
+
+Bugfixes:
+* alm() occurence!="none" would not work if the response variable was not called "y".
+* Corrected the table with levels for the predict function.
+* A fix in alm() with ar and complicated names of variables.
+* A fix in subset parameter for normal distribution.
+* The parameters are now standardised for the hessian calculation in order to get more accurate vcov.
+
+
 greybox v0.6.2 (Release data: 2020-09-02)
 ==============
 

diff --git a/R/alm.R b/R/alm.R
@@ -8,7 +8,7 @@
 #' \item \link[stats]{dnorm} - Normal distribution,
 #' \item \link[greybox]{dlaplace} - Laplace distribution,
 #' \item \link[greybox]{ds} - S-distribution,
-#' \item \link[gnorm]{dgnorm} - Generalised Normal distribution,
+#' \item dgnorm - Generalised Normal distribution,
 #' \item \link[stats]{dlogis} - Logistic Distribution,
 #' \item \link[stats]{dt} - T-distribution,
 #' \item \link[greybox]{dalaplace} - Asymmetric Laplace distribution,
@@ -190,55 +190,49 @@
 #' xreg <- cbind(rlaplace(100,10,3),rnorm(100,50,5))
 #' xreg <- cbind(100+0.5*xreg[,1]-0.75*xreg[,2]+rlaplace(100,0,3),xreg,rnorm(100,300,10))
 #' colnames(xreg) <- c("y","x1","x2","Noise")
-#' inSample <- xreg[1:80,]
-#' outSample <- xreg[-c(1:80),]
 #'
 #' # An example with Laplace distribution
-#' ourModel <- alm(y~x1+x2, inSample, distribution="dlaplace")
+#' ourModel <- alm(y~x1+x2, xreg, subset=c(1:80), distribution="dlaplace")
 #' summary(ourModel)
-#' plot(predict(ourModel,outSample))
+#' plot(predict(ourModel,xreg[-c(1:80),]))
 #'
 #' # And another one with Asymmetric Laplace distribution (quantile regression)
 #' # with optimised alpha
-#' ourModel <- alm(y~x1+x2, inSample, distribution="dalaplace")
+#' ourModel <- alm(y~x1+x2, xreg, subset=c(1:80), distribution="dalaplace")
 #' summary(ourModel)
-#' plot(predict(ourModel,outSample))
+#' plot(predict(ourModel,xreg[-c(1:80),]))
 #'
 #' # An example with AR(1) order
-#' ourModel <- alm(y~x1+x2, inSample, distribution="dnorm", ar=1)
+#' ourModel <- alm(y~x1+x2, xreg, subset=c(1:80), distribution="dnorm", ar=1)
 #' summary(ourModel)
-#' plot(predict(ourModel,outSample))
+#' plot(predict(ourModel,xreg[-c(1:80),]))
 #'
 #' ### Examples with the count data
 #' xreg[,1] <- round(exp(xreg[,1]-70),0)
-#' inSample <- xreg[1:80,]
-#' outSample <- xreg[-c(1:80),]
 #'
 #' # Negative Binomial distribution
-#' ourModel <- alm(y~x1+x2, inSample, distribution="dnbinom")
+#' ourModel <- alm(y~x1+x2, xreg, subset=c(1:80), distribution="dnbinom")
 #' summary(ourModel)
-#' predict(ourModel,outSample,interval="p",side="u")
+#' predict(ourModel,xreg[-c(1:80),],interval="p",side="u")
 #'
 #' # Poisson distribution
-#' ourModel <- alm(y~x1+x2, inSample, distribution="dpois")
+#' ourModel <- alm(y~x1+x2, xreg, subset=c(1:80), distribution="dpois")
 #' summary(ourModel)
-#' predict(ourModel,outSample,interval="p",side="u")
+#' predict(ourModel,xreg[-c(1:80),],interval="p",side="u")
 #'
 #'
 #' ### Examples with binary response variable
 #' xreg[,1] <- round(xreg[,1] / (1 + xreg[,1]),0)
-#' inSample <- xreg[1:80,]
-#' outSample <- xreg[-c(1:80),]
 #'
 #' # Logistic distribution (logit regression)
-#' ourModel <- alm(y~x1+x2, inSample, distribution="plogis")
+#' ourModel <- alm(y~x1+x2, xreg, subset=c(1:80), distribution="plogis")
 #' summary(ourModel)
-#' plot(predict(ourModel,outSample,interval="c"))
+#' plot(predict(ourModel,xreg[-c(1:80),],interval="c"))
 #'
 #' # Normal distribution (probit regression)
-#' ourModel <- alm(y~x1+x2, inSample, distribution="pnorm")
+#' ourModel <- alm(y~x1+x2, xreg, subset=c(1:80), distribution="pnorm")
 #' summary(ourModel)
-#' plot(predict(ourModel,outSample,interval="p"))
+#' plot(predict(ourModel,xreg[-c(1:80),],interval="p"))
 #'
 #' @importFrom pracma hessian
 #' @importFrom nloptr nloptr
@@ -256,6 +250,7 @@ alm <- function(formula, data, subset, na.action,
                                "plogis","pnorm"),
                 loss=c("likelihood","MSE","MAE","HAM","LASSO","RIDGE"),
                 occurrence=c("none","plogis","pnorm"),
+                # scaleFormula=NULL,
                 ar=0,# i=0,
                 parameters=NULL, fast=FALSE, ...){
 # Useful stuff for dnbinom: https://scialert.net/fulltext/?doi=ajms.2010.1.15
@@ -477,7 +472,10 @@ alm <- function(formula, data, subset, na.action,
         return(fitterReturn);
     }
 
-    CF <- function(B, distribution, loss, y, matrixXreg, recursiveModel, denominator){
+    CF <- function(B, distribution, loss, y, matrixXreg, recursiveModel, denominator, hessianCalculation=FALSE){
+        if(hessianCalculation){
+            B[] <- B * denominator;
+        }
         if(recursiveModel){
             fitterReturn <- fitterRecursive(B, distribution, y, matrixXreg);
         }
@@ -502,7 +500,7 @@ alm <- function(formula, data, subset, na.action,
                                                           scale=fitterReturn$scale, log=TRUE)-log(y[otU]),
                                    "dls" = ds(log(y[otU]), mu=fitterReturn$mu[otU], scale=fitterReturn$scale, log=TRUE)-log(y[otU]),
                                    "dlgnorm" = dgnorm(log(y[otU]), mu=fitterReturn$mu[otU], alpha=fitterReturn$scale,
-                                                     beta=fitterReturn$other, log=TRUE)-log(y[otU]),
+                                                      beta=fitterReturn$other, log=TRUE)-log(y[otU]),
                                    "dbcnorm" = dbcnorm(y[otU], mu=fitterReturn$mu[otU], sigma=fitterReturn$scale,
                                                        lambda=fitterReturn$other, log=TRUE),
                                    "dfnorm" = dfnorm(y[otU], mu=fitterReturn$mu[otU], sigma=fitterReturn$scale, log=TRUE),
@@ -1005,7 +1003,7 @@ alm <- function(formula, data, subset, na.action,
         CDF <- FALSE;
     }
 
-    if(CDF & any(y!=0 & y!=1)){
+    if(CDF && any(y!=0 & y!=1)){
         y[] <- (y!=0)*1;
     }
 
@@ -1122,6 +1120,10 @@ alm <- function(formula, data, subset, na.action,
     # The number of exogenous variables (no ARI elements)
     nVariablesExo <- nVariables;
 
+    #### The model for the scale ####
+    # if(!is.null(scaleFormula)){
+    # }
+
     #### Estimate parameters of the model ####
     if(is.null(parameters)){
         #### Add AR and I elements in the regression ####
@@ -1420,14 +1422,14 @@ alm <- function(formula, data, subset, na.action,
                                 maxtime=maxtime, xtol_abs=xtol_abs, ftol_rel=ftol_rel, ftol_abs=ftol_abs),
                       lb=BLower, ub=BUpper,
                       distribution=distribution, loss=loss, y=y, matrixXreg=matrixXreg,
-                      recursiveModel=recursiveModel, denominator=denominator);
+                      recursiveModel=recursiveModel, denominator=denominator, hessianCalculation=FALSE);
         if(recursiveModel){
             res2 <- nloptr(res$solution, CF,
                            opts=list(algorithm=algorithm, xtol_rel=xtol_rel, maxeval=maxeval, print_level=print_level,
                                 maxtime=maxtime, xtol_abs=xtol_abs, ftol_rel=ftol_rel, ftol_abs=ftol_abs),
                            lb=BLower, ub=BUpper,
                            distribution=distribution, loss=loss, y=y, matrixXreg=matrixXreg,
-                           recursiveModel=recursiveModel, denominator=denominator);
+                           recursiveModel=recursiveModel, denominator=denominator, hessianCalculation=FALSE);
             if(res2$objective<res$objective){
                 res[] <- res2;
             }
@@ -1700,9 +1702,16 @@ alm <- function(formula, data, subset, na.action,
     if(FI){
         # Only vcov is needed, no point in redoing the occurrenceModel
         occurrenceModel <- FALSE;
-        FI <- hessian(CF, B,
+        # Standartisation is needed in order to get accurate hessian
+        denominator <- apply(matrixXreg, 2, sd);
+        # No variability, substitute by 1
+        denominator[is.infinite(denominator)] <- 1;
+        # Fix the denominator for the intercept
+        denominator[denominator==0] <- B[denominator==0];
+        BNew <- B / denominator;
+        FI <- hessian(CF, BNew,
                       distribution=distribution, loss=loss, y=y, matrixXreg=matrixXreg,
-                      recursiveModel=recursiveModel, denominator=denominator);
+                      recursiveModel=recursiveModel, denominator=denominator, hessianCalculation=TRUE);
 
         if(any(is.nan(FI))){
             warning(paste0("Something went wrong and we failed to produce the covariance matrix of the parameters.\n",
@@ -1726,6 +1735,7 @@ alm <- function(formula, data, subset, na.action,
             dataWork <- cbind(y,matrixXreg);
             variablesUsed <- variablesNames;
         }
+        colnames(dataWork)[1] <- responseName;
 
         # If there are NaN values, remove the respective observations
         if(any(sapply(mf$data,is.nan))){

diff --git a/R/detector-dst-leap.R b/R/detector-dst-leap.R
@@ -2,15 +2,18 @@
 #'
 #' Functions to detect, when Daylight Saving Time and leap year start and finish
 #'
-#' The \code{detectdst} function detects, when the change for the DST starts and ends.
+#' The \code{detectdst} function detects, when the change for the DST starts and ends. It
+#' assumes that the frequency of data is not lower than hourly.
 #' The \code{detectleap} function does similar for the leap year, but flagging the 29th
 #' of February as a starting and to the 28th of February next year as the ending dates.
 #'
 #' In order for the methods to work, the object needs to be of either zoo / xts or POSIXt
 #' class and should contain valid dates.
 #'
 #' @param object Either a zoo / xts object or a vector of dates / times in POSIXt / Date
-#' class.
+#' class. Note that in order for \code{detectdst()} to work correctly, your data should
+#' not have missing observations. Otherwise it might not be possible to locate, when DST
+#' happens.
 #'
 #' @return List containing:
 #' \itemize{
@@ -45,12 +48,28 @@ detectdst.default <- function(object){
 detectdst.POSIXt <- function(object){
     # Function detects the dst days and flags them.
     # It returns a matrix with coordinates of the days that contain the change
+    DSTIssues <- FALSE;
 
     # The Spring day
     DSTTime1 <- which(diff(as.numeric(strftime(object,"%H")))==2)+1;
+    # If there are non-spring times, this might be an issue with the data
+    if(any(!(strftime(object[DSTTime1],"%m") %in% c("03","04")))){
+        DSTTime1 <- DSTTime1[strftime(object[DSTTime1],"%m") %in% c("03","04")];
+        DSTIssues[] <- TRUE;
+    }
 
     # The Autumn day
     DSTTime2 <- which(diff(as.numeric(strftime(object,"%H")))==0)+1;
+    # If there are non-spring times, this might be an issue with the data
+    if(any(!(strftime(object[DSTTime2],"%m") %in% c("10","11")))){
+        DSTTime2 <- DSTTime2[strftime(object[DSTTime2],"%m") %in% c("10","11")];
+        DSTIssues[] <- TRUE;
+    }
+
+    if(DSTIssues){
+        warning("It seems that you have missing observations. The function might return wrong dates.",
+                call.=FALSE);
+    }
 
     return(list(start=data.frame(id=DSTTime1,date=object[DSTTime1]),
                 end=data.frame(id=DSTTime2,date=object[DSTTime2])))