Automated build.

.setup/cpp/cTMed-direct-std.cpp

@@ -18,10 +18,14 @@ double DirectStd(const arma::mat& phi, const arma::mat& sigma,
   arma::mat total = arma::expmat(delta_t * phi);
   arma::mat phi_hashtag = arma::kron(phi, I) + arma::kron(I, phi);
   arma::vec sigma_vec = arma::vectorise(sigma);
-  arma::vec psi_vec = arma::inv(phi_hashtag) *
-                      (arma::expmat(phi_hashtag * delta_t) - J) * sigma_vec;
+  // arma::vec psi_vec = arma::inv(phi_hashtag) * (arma::expmat(phi_hashtag *
+  // delta_t) - J) * sigma_vec;
+  arma::vec psi_vec = arma::solve(
+      phi_hashtag, (arma::expmat(phi_hashtag * delta_t) - J) * sigma_vec);
+  // arma::mat total_cov = arma::reshape(arma::inv(J - arma::kron(total, total))
+  // * psi_vec, phi.n_rows, phi.n_cols);
   arma::mat total_cov =
-      arma::reshape(arma::inv(J - arma::kron(total, total)) * psi_vec,
+      arma::reshape(arma::solve(J - arma::kron(total, total), psi_vec),
                     phi.n_rows, phi.n_cols);
   arma::mat sd_row = arma::diagmat(arma::sqrt(total_cov.diag()));
   arma::mat sd_col_inv = arma::diagmat(1.0 / arma::sqrt(total_cov.diag()));

.setup/cpp/cTMed-total-cov.cpp → .setup/cpp/cTMed-exp-cov.cpp

@@ -1,20 +1,20 @@
 // -----------------------------------------------------------------------------
-// edit .setup/cpp/cTMed-total-cov.cpp
+// edit .setup/cpp/cTMed-exp-cov.cpp
 // Ivan Jacob Agaloos Pesigan
 // -----------------------------------------------------------------------------
 
 #include <RcppArmadillo.h>
 // [[Rcpp::depends(RcppArmadillo)]]
-// [[Rcpp::export(.TotalCov)]]
-arma::mat TotalCov(const arma::mat& phi, const arma::mat& sigma,
-                   const double& delta_t) {
+// [[Rcpp::export(.ExpCov)]]
+arma::mat ExpCov(const arma::mat& phi, const arma::mat& sigma,
+                 const double& delta_t) {
   arma::mat I = arma::eye(phi.n_rows, phi.n_cols);
   arma::mat J = arma::eye(phi.n_rows * phi.n_cols, phi.n_rows * phi.n_cols);
-  arma::mat total = arma::expmat(delta_t * phi);
+  arma::mat beta = arma::expmat(delta_t * phi);
   arma::mat phi_hashtag = arma::kron(phi, I) + arma::kron(I, phi);
   arma::vec sigma_vec = arma::vectorise(sigma);
   arma::vec psi_vec = arma::inv(phi_hashtag) *
                       (arma::expmat(phi_hashtag * delta_t) - J) * sigma_vec;
-  return arma::reshape(arma::inv(J - arma::kron(total, total)) * psi_vec,
+  return arma::reshape(arma::inv(J - arma::kron(beta, beta)) * psi_vec,
                        phi.n_rows, phi.n_cols);
 }

.setup/cpp/cTMed-exp-mean.cpp

@@ -0,0 +1,20 @@
+// -----------------------------------------------------------------------------
+// edit .setup/cpp/cTMed-exp-mean.cpp
+// Ivan Jacob Agaloos Pesigan
+// -----------------------------------------------------------------------------
+
+#include <RcppArmadillo.h>
+// [[Rcpp::depends(RcppArmadillo)]]
+// [[Rcpp::export(.ExpMean)]]
+Rcpp::NumericVector ExpMean(const arma::mat& phi, const arma::vec& iota,
+                            const double& delta_t) {
+  if (arma::all(iota == 0)) {
+    return Rcpp::NumericVector(phi.n_rows, 0.0);
+  }
+  arma::mat I = arma::eye(phi.n_rows, phi.n_cols);
+  arma::mat beta = arma::expmat(delta_t * phi);
+  arma::vec alpha = arma::solve(phi, (beta - I) * iota);
+  arma::vec mu = arma::solve((I - beta), alpha);
+  Rcpp::NumericVector output(mu.begin(), mu.end());
+  return output;
+}

.setup/cpp/cTMed-indirect-std.cpp

@@ -18,10 +18,14 @@ double IndirectStd(const arma::mat& phi, const arma::mat& sigma,
   arma::mat total = arma::expmat(delta_t * phi);
   arma::mat phi_hashtag = arma::kron(phi, I) + arma::kron(I, phi);
   arma::vec sigma_vec = arma::vectorise(sigma);
-  arma::vec psi_vec = arma::inv(phi_hashtag) *
-                      (arma::expmat(phi_hashtag * delta_t) - J) * sigma_vec;
+  // arma::vec psi_vec = arma::inv(phi_hashtag) * (arma::expmat(phi_hashtag *
+  // delta_t) - J) * sigma_vec;
+  arma::vec psi_vec = arma::solve(
+      phi_hashtag, (arma::expmat(phi_hashtag * delta_t) - J) * sigma_vec);
+  // arma::mat total_cov = arma::reshape(arma::inv(J - arma::kron(total, total))
+  // * psi_vec, phi.n_rows, phi.n_cols);
   arma::mat total_cov =
-      arma::reshape(arma::inv(J - arma::kron(total, total)) * psi_vec,
+      arma::reshape(arma::solve(J - arma::kron(total, total), psi_vec),
                     phi.n_rows, phi.n_cols);
   arma::mat sd_row = arma::diagmat(arma::sqrt(total_cov.diag()));
   arma::mat sd_col_inv = arma::diagmat(1.0 / arma::sqrt(total_cov.diag()));

.setup/cpp/cTMed-mc-phi-sigma-i.cpp

@@ -6,52 +6,50 @@
 #include <RcppArmadillo.h>
 // [[Rcpp::depends(RcppArmadillo)]]
 // [[Rcpp::export(.MCPhiSigmaI)]]
-Rcpp::List MCPhiSigmaI(const arma::mat& phi, const arma::mat& vcov_phi_vec,
-                       const arma::mat& sigma, const arma::mat& vcov_sigma_vech,
+Rcpp::List MCPhiSigmaI(const arma::vec& theta, const arma::mat& vcov_theta,
                        bool test_phi = true) {
   Rcpp::List output(2);
-  // phi
-  arma::mat phi_i(phi.n_rows, phi.n_cols, arma::fill::none);
-  arma::vec phi_vec = arma::vectorise(phi);
-  arma::vec phi_vec_i(phi.n_rows * phi.n_cols, arma::fill::none);
+  arma::uword n = theta.n_elem;
+  arma::uword p = (-1 + std::sqrt(1 + 24 * n)) / 6;
+  arma::uword q = (p * (p + 1)) / 2;
+  arma::uword index = 0;
+  arma::vec v_i(n, arma::fill::none);
+  arma::mat phi_i(p, p, arma::fill::none);
+  arma::vec phi_vec_i(p * p, arma::fill::none);
+  arma::mat sigma_i(p, p, arma::fill::none);
+  arma::vec sigma_vech_i(q, arma::fill::none);
+  arma::vec eigval;
+  arma::mat eigvec;
   bool run = true;
   while (run) {
     // generate data
-    phi_vec_i = arma::mvnrnd(phi_vec, vcov_phi_vec);
-    phi_i = arma::reshape(phi_vec_i, phi.n_rows, phi.n_cols);
+    v_i = arma::mvnrnd(theta, vcov_theta);
+    phi_vec_i = v_i(arma::span(0, (p * p) - 1));
+    sigma_vech_i = v_i(arma::span(p * p, n - 1));
     // test phi
+    phi_i = arma::reshape(phi_vec_i, p, p);
     if (test_phi) {
       if (TestPhi(phi_i)) {
         run = false;
       }
     } else {
       run = false;
     }
-  }
-  // sigma
-  arma::uword q = (sigma.n_rows * (sigma.n_rows + 1)) / 2;
-  arma::vec sigma_vech(q, arma::fill::none);
-  arma::uword index = 0;
-  for (arma::uword j = 0; j < sigma.n_cols; ++j) {
-    for (arma::uword i = j; i < sigma.n_rows; ++i) {
-      sigma_vech(index++) = sigma(i, j);
-    }
-  }
-  arma::vec sigma_vech_i = arma::mvnrnd(sigma_vech, vcov_sigma_vech);
-  arma::mat sigma_i(sigma.n_rows, sigma.n_rows, arma::fill::zeros);
-  index = 0;
-  for (arma::uword i = 0; i < sigma.n_rows; ++i) {
-    for (arma::uword j = i; j < sigma.n_cols; ++j) {
-      sigma_i(i, j) = sigma_vech_i(index);
-      sigma_i(j, i) = sigma_vech_i(index);
-      index++;
+    if (run == false) {
+      // test sigma
+      index = 0;
+      for (arma::uword i = 0; i < p; ++i) {
+        for (arma::uword j = i; j < p; ++j) {
+          sigma_i(i, j) = sigma_vech_i(index);
+          sigma_i(j, i) = sigma_vech_i(index);
+          index++;
+        }
+      }
+      arma::eig_sym(eigval, eigvec, sigma_i);
+      eigval.transform([](double val) { return std::max(val, 1e-8); });
+      sigma_i = eigvec * arma::diagmat(eigval) * eigvec.t();
     }
   }
-  arma::vec eigval;
-  arma::mat eigvec;
-  arma::eig_sym(eigval, eigvec, sigma_i);
-  eigval.transform([](double val) { return std::max(val, 1e-8); });
-  sigma_i = eigvec * arma::diagmat(eigval) * eigvec.t();
   output[0] = phi_i;
   output[1] = sigma_i;
   return output;

.setup/cpp/cTMed-med-std-s.cpp

@@ -31,10 +31,15 @@ arma::mat MedStds(const arma::mat& phi, const arma::mat& sigma,
   arma::mat direct_std(phi.n_rows, phi.n_cols, arma::fill::none);
   for (arma::uword t = 0; t < delta_t.n_rows; t++) {
     total = arma::expmat(delta_t[t] * phi);
-    psi_vec = arma::inv(phi_hashtag) *
-              (arma::expmat(phi_hashtag * delta_t[t]) - J) * sigma_vec;
-    total_cov = arma::reshape(arma::inv(J - arma::kron(total, total)) * psi_vec,
-                              phi.n_rows, phi.n_cols);
+    // psi_vec = arma::inv(phi_hashtag) * (arma::expmat(phi_hashtag *
+    // delta_t[t]) - J) * sigma_vec;
+    psi_vec = arma::solve(
+        phi_hashtag, (arma::expmat(phi_hashtag * delta_t[t]) - J) * sigma_vec);
+    // total_cov = arma::reshape(arma::inv(J - arma::kron(total, total)) *
+    // psi_vec, phi.n_rows, phi.n_cols);
+    total_cov =
+        arma::reshape(arma::solve(J - arma::kron(total, total), psi_vec),
+                      phi.n_rows, phi.n_cols);
     sd_row = arma::diagmat(arma::sqrt(total_cov.diag()));
     sd_col_inv = arma::diagmat(1.0 / arma::sqrt(total_cov.diag()));
     total_std = sd_row * total * sd_col_inv;

.setup/cpp/cTMed-med-std-vec.cpp

@@ -9,7 +9,7 @@
 Rcpp::NumericVector MedStdVec(const arma::vec& v, const double& delta_t,
                               const arma::uword& from, const arma::uword& to,
                               const arma::vec& med) {
-  arma::uword q = v.size();
+  arma::uword q = v.n_elem;
   arma::uword p = (-1 + std::sqrt(1 + 24 * q)) / 6;
   arma::mat phi = arma::mat(v.subvec(0, p * p - 1)).reshape(p, p);
   arma::vec sigma_vech = v.subvec(p * p, q - 1);
@@ -33,10 +33,15 @@ Rcpp::NumericVector MedStdVec(const arma::vec& v, const double& delta_t,
   arma::mat total = arma::expmat(delta_t * phi);
   arma::mat phi_hashtag = arma::kron(phi, I) + arma::kron(I, phi);
   arma::vec sigma_vec = arma::vectorise(sigma);
-  arma::vec psi_vec = arma::inv(phi_hashtag) *
-                      (arma::expmat(phi_hashtag * delta_t) - J) * sigma_vec;
+  // arma::vec psi_vec = arma::inv(phi_hashtag) * (arma::expmat(phi_hashtag *
+  // delta_t) - J) * sigma_vec;
+  arma::vec psi_vec = arma::solve(
+      phi_hashtag, (arma::expmat(phi_hashtag * delta_t) - J) * sigma_vec);
+  // arma::mat total_cov = arma::reshape(arma::inv(J - arma::kron(total, total))
+  // * psi_vec, phi.n_rows, phi.n_cols);
   arma::mat total_cov =
-      arma::reshape(arma::inv(J - arma::kron(total, total)) * psi_vec, p, p);
+      arma::reshape(arma::solve(J - arma::kron(total, total), psi_vec),
+                    phi.n_rows, phi.n_cols);
   arma::mat sd_row = arma::diagmat(arma::sqrt(total_cov.diag()));
   arma::mat sd_col_inv = arma::diagmat(1.0 / arma::sqrt(total_cov.diag()));
   arma::mat total_std = sd_row * total * sd_col_inv;

.setup/cpp/cTMed-med-std.cpp

@@ -18,10 +18,14 @@ Rcpp::NumericVector MedStd(const arma::mat& phi, const arma::mat& sigma,
   arma::mat total = arma::expmat(delta_t * phi);
   arma::mat phi_hashtag = arma::kron(phi, I) + arma::kron(I, phi);
   arma::vec sigma_vec = arma::vectorise(sigma);
-  arma::vec psi_vec = arma::inv(phi_hashtag) *
-                      (arma::expmat(phi_hashtag * delta_t) - J) * sigma_vec;
+  // arma::vec psi_vec = arma::inv(phi_hashtag) * (arma::expmat(phi_hashtag *
+  // delta_t) - J) * sigma_vec;
+  arma::vec psi_vec = arma::solve(
+      phi_hashtag, (arma::expmat(phi_hashtag * delta_t) - J) * sigma_vec);
+  // arma::mat total_cov = arma::reshape(arma::inv(J - arma::kron(total, total))
+  // * psi_vec, phi.n_rows, phi.n_cols);
   arma::mat total_cov =
-      arma::reshape(arma::inv(J - arma::kron(total, total)) * psi_vec,
+      arma::reshape(arma::solve(J - arma::kron(total, total), psi_vec),
                     phi.n_rows, phi.n_cols);
   arma::mat sd_row = arma::diagmat(arma::sqrt(total_cov.diag()));
   arma::mat sd_col_inv = arma::diagmat(1.0 / arma::sqrt(total_cov.diag()));

.setup/cpp/cTMed-total-std-delta-t.cpp

@@ -13,10 +13,14 @@ Rcpp::NumericVector TotalStdDeltaT(const arma::mat& phi, const arma::mat& sigma,
   arma::mat total = arma::expmat(delta_t * phi);
   arma::mat phi_hashtag = arma::kron(phi, I) + arma::kron(I, phi);
   arma::vec sigma_vec = arma::vectorise(sigma);
-  arma::vec psi_vec = arma::inv(phi_hashtag) *
-                      (arma::expmat(phi_hashtag * delta_t) - J) * sigma_vec;
+  // arma::vec psi_vec = arma::inv(phi_hashtag) * (arma::expmat(phi_hashtag *
+  // delta_t) - J) * sigma_vec;
+  arma::vec psi_vec = arma::solve(
+      phi_hashtag, (arma::expmat(phi_hashtag * delta_t) - J) * sigma_vec);
+  // arma::mat total_cov = arma::reshape(arma::inv(J - arma::kron(total, total))
+  // * psi_vec, phi.n_rows, phi.n_cols);
   arma::mat total_cov =
-      arma::reshape(arma::inv(J - arma::kron(total, total)) * psi_vec,
+      arma::reshape(arma::solve(J - arma::kron(total, total), psi_vec),
                     phi.n_rows, phi.n_cols);
   arma::mat sd_row = arma::diagmat(arma::sqrt(total_cov.diag()));
   arma::mat sd_col_inv = arma::diagmat(1.0 / arma::sqrt(total_cov.diag()));

.setup/cpp/cTMed-total-std-vec.cpp

@@ -7,7 +7,7 @@
 // [[Rcpp::depends(RcppArmadillo)]]
 // [[Rcpp::export(.TotalStdVec)]]
 arma::vec TotalStdVec(const arma::vec& v, const double& delta_t) {
-  arma::uword q = v.size();
+  arma::uword q = v.n_elem;
   arma::uword p = (-1 + std::sqrt(1 + 24 * q)) / 6;
   arma::mat phi = arma::mat(v.subvec(0, p * p - 1)).reshape(p, p);
   arma::vec sigma_vech = v.subvec(p * p, q - 1);
@@ -27,10 +27,15 @@ arma::vec TotalStdVec(const arma::vec& v, const double& delta_t) {
   arma::mat total = arma::expmat(delta_t * phi);
   arma::mat phi_hashtag = arma::kron(phi, I) + arma::kron(I, phi);
   arma::vec sigma_vec = arma::vectorise(sigma);
-  arma::vec psi_vec = arma::inv(phi_hashtag) *
-                      (arma::expmat(phi_hashtag * delta_t) - J) * sigma_vec;
+  // arma::vec psi_vec = arma::inv(phi_hashtag) * (arma::expmat(phi_hashtag *
+  // delta_t) - J) * sigma_vec;
+  arma::vec psi_vec = arma::solve(
+      phi_hashtag, (arma::expmat(phi_hashtag * delta_t) - J) * sigma_vec);
+  // arma::mat total_cov = arma::reshape(arma::inv(J - arma::kron(total, total))
+  // * psi_vec, phi.n_rows, phi.n_cols);
   arma::mat total_cov =
-      arma::reshape(arma::inv(J - arma::kron(total, total)) * psi_vec, p, p);
+      arma::reshape(arma::solve(J - arma::kron(total, total), psi_vec),
+                    phi.n_rows, phi.n_cols);
   arma::mat sd_row = arma::diagmat(arma::sqrt(total_cov.diag()));
   arma::mat sd_col_inv = arma::diagmat(1.0 / arma::sqrt(total_cov.diag()));
   arma::mat total_std = sd_row * total * sd_col_inv;

.setup/cpp/cTMed-total-std.cpp

@@ -13,10 +13,14 @@ arma::mat TotalStd(const arma::mat& phi, const arma::mat& sigma,
   arma::mat total = arma::expmat(delta_t * phi);
   arma::mat phi_hashtag = arma::kron(phi, I) + arma::kron(I, phi);
   arma::vec sigma_vec = arma::vectorise(sigma);
-  arma::vec psi_vec = arma::inv(phi_hashtag) *
-                      (arma::expmat(phi_hashtag * delta_t) - J) * sigma_vec;
+  // arma::vec psi_vec = arma::inv(phi_hashtag) * (arma::expmat(phi_hashtag *
+  // delta_t) - J) * sigma_vec;
+  arma::vec psi_vec = arma::solve(
+      phi_hashtag, (arma::expmat(phi_hashtag * delta_t) - J) * sigma_vec);
+  // arma::mat total_cov = arma::reshape(arma::inv(J - arma::kron(total, total))
+  // * psi_vec, phi.n_rows, phi.n_cols);
   arma::mat total_cov =
-      arma::reshape(arma::inv(J - arma::kron(total, total)) * psi_vec,
+      arma::reshape(arma::solve(J - arma::kron(total, total), psi_vec),
                     phi.n_rows, phi.n_cols);
   arma::mat sd_row = arma::diagmat(arma::sqrt(total_cov.diag()));
   arma::mat sd_col_inv = arma::diagmat(1.0 / arma::sqrt(total_cov.diag()));