Revision 86019c746fd23b26ad2f3d811953cf3623e18981 authored by Wenchao Ma on 28 March 2018, 07:18:50 UTC, committed by cran-robot on 28 March 2018, 07:18:50 UTC
1 parent a55a4dc
NgRg.cpp
#include <RcppArmadillo.h>
// [[Rcpp::depends(RcppArmadillo)]]
using namespace Rcpp;
// [[Rcpp::export]]
List NgRg(arma::mat & mlogPost,
arma::mat & mX,
arma::mat & mloc,
arma::vec weights){
int J = mX.n_cols;
//int L = mlogPost.n_cols;
int Ljmax = mloc.max();
arma::mat mIndmiss = arma::zeros<arma::mat>(arma::size(mX));
mIndmiss.elem( arma::find_finite(mX) ).ones(); //missing - 0;nonmissing - 1
arma::mat Ng = arma::zeros<arma::mat>(J,Ljmax);
arma::mat Rg = arma::zeros<arma::mat>(J,Ljmax);
arma::mat mPost = exp ( mlogPost ); //N x L
mPost.each_col() %= weights;
--mloc;
for (int j=0;j<J;++j){ //for each item
int Kjmax = mloc.row(j).max()+1;
//missing values in item j are removed from posterior and X
arma::mat mPostmiss=mPost.rows(arma::find(mIndmiss.col(j)==1));
arma::mat mXmiss=mX.rows(arma::find(mIndmiss.col(j)==1));
for (int k=0;k<Kjmax;++k){
arma::mat tmp = mPostmiss.cols(arma::find(mloc.row(j)==k));
Ng(j,k) = accu(tmp);
Rg(j,k) = accu((mXmiss.col(j)*arma::ones<arma::rowvec>(tmp.n_cols)) % tmp);
}
}
return Rcpp::List::create(Rcpp::Named("Ng") = Ng,
Rcpp::Named("Rg") = Rg);
}
Computing file changes ...