updates

Author: drizopoulos

src/JMbayes2_Funs.h

@@ -158,11 +158,22 @@ field<mat> mat2field (const mat &b, const field<uvec> &ind_RE) {
   uword n = ind_RE.n_elem;
   field<mat> out(n);
   for (uword i = 0; i < n; i++) {
-    out.at(i) = b.cols(ind_RE.at(i));
+      uword first_col = ind_RE.at(i).front();
+      uword last_col  = ind_RE.at(i).back();
+      out.at(i) = b.cols(first_col, last_col);
   }
   return out;
 }
 
+inline void mat2field_inplace (field<mat> &out, const mat &b, const field<uvec> &ind_RE) {
+    uword n = ind_RE.n_elem;
+    for (uword i = 0; i < n; ++i) {
+        uword first_col = ind_RE.at(i).front();
+        uword last_col  = ind_RE.at(i).back();
+        out.at(i) = b.cols(first_col, last_col);
+    }
+}
+
 field<vec> vec2field (const vec &betas, const field<uvec> &ind_FE) {
   uword n = ind_FE.n_elem;
   field<vec> out(n);

src/JMbayes2_RE.h

@@ -50,10 +50,15 @@ void update_b (field<mat> &b, mat &b_mat, field<vec> &eta,
     double log_det_V_Sigma = -arma::sum(arma::log(L.diag())) - arma::sum(arma::log(sds));
     double other_terms = -(double)nRE / 2.0 * log2pi + log_det_V_Sigma;
     vec denominator_b = logLik_long + logLik_surv + logLik_re;
+    field<mat> proposed_b(ind_RE.n_elem);
+    for (uword i = 0; i < ind_RE.n_elem; ++i) {
+        proposed_b.at(i).set_size(b_mat.n_rows, ind_RE.at(i).n_elem);
+    }
     for (uword j = 0; j < nRE; ++j) {
         vec old_b_j = b_mat.col(j);
         b_mat.col(j) += scale_b.col(j) % arma::randn<arma::vec>(n);
-        field<mat> proposed_b = mat2field(b_mat, ind_RE);
+        //field<mat> proposed_b = mat2field(b_mat, ind_RE);
+        mat2field_inplace(proposed_b, b_mat, ind_RE);
         field<vec> eta_proposed = linpred_mixed(X, betas, Z, proposed_b, idL);
         vec logLik_long_proposed = log_long(y, eta_proposed, sigmas, extra_parms,
                                             families, links, ids, unq_ids, n);