Use Choleski decomposition for the calculation of the inverse of matr…

…ices in alm. Related to #14

DESCRIPTION

@@ -1,7 +1,7 @@
 Package: greybox
 Type: Package
 Title: Toolbox for Model Building and Forecasting
-Version: 0.3.1.41006
+Version: 0.3.1.41007
 Date: 2018-08-09
 Authors@R: person("Ivan", "Svetunkov", email = "[email protected]", role = c("aut", "cre"),
                   comment="Lecturer at Centre for Marketing Analytics and Forecasting, Lancaster University, UK")

NEWS

@@ -11,6 +11,7 @@ Changes:
 * alm() can be constructed with the provided vector of parameters (needed for vcov method).
 * We now use well-known analytical solutions for the cases of distribution="norm" of alm() and other functions.
 * Code of lmCombine and lmDynamic is slightly simplified.
+* We now use Choleski decomposition for the calculation of the inverse of matrices in alm.
 
 Bugfixes:
 * Fixed a bug with the style="line" in rmc(), where the grouping would be wrong in cases, when one method significantly differs from the others.

R/alm.R

@@ -220,10 +220,12 @@ alm <- function(formula, data, subset=NULL,  na.action,
 
     if(is.null(A)){
         if(distribution=="lnorm"){
-            A <- as.vector(solve(t(matrixXreg) %*% matrixXreg, tol=1e-10) %*% t(matrixXreg) %*% log(y));
+            # A <- as.vector(solve(t(matrixXreg) %*% matrixXreg, t(matrixXreg) %*% log(y)));
+            A <- as.vector(chol2inv(chol(t(matrixXreg) %*% matrixXreg)) %*% t(matrixXreg) %*% log(y));
         }
         else{
-            A <- as.vector(solve(t(matrixXreg) %*% matrixXreg, tol=1e-10) %*% t(matrixXreg) %*% y);
+            # A <- as.vector(solve(t(matrixXreg) %*% matrixXreg, t(matrixXreg) %*% y));
+            A <- as.vector(chol2inv(chol(t(matrixXreg) %*% matrixXreg)) %*% t(matrixXreg) %*% y);
         }
 
         # Although this is not needed in case of distribution="norm", we do that in a way, for the code consistency purposes
@@ -259,7 +261,8 @@ alm <- function(formula, data, subset=NULL,  na.action,
                 vcovMatrix <- diag(nVariables);
             }
             else{
-                vcovMatrix <- solve(vcovMatrix);
+                # vcovMatrix <- solve(vcovMatrix, diag(nVariables));
+                vcovMatrix <- chol2inv(chol(vcovMatrix));
             }
 
             # Sometimes the diagonal elements in the covariance matrix are negative because likelihood is not fully maximised...