Commit 01d0925d authored by hilaire drouineau's avatar hilaire drouineau
Browse files

try a method to decrease autocorrelation

parent c8e491dc
......@@ -5,7 +5,7 @@
#' @useDynLib RCaN
NULL
fitCaN <- function(N, A, b, C, v, L, x0, thin) {
.Call('_RCaN_fitCaN', PACKAGE = 'RCaN', N, A, b, C, v, L, x0, thin)
fitCaN <- function(N, A, b, C, v, L, x0, thin, test) {
.Call('_RCaN_fitCaN', PACKAGE = 'RCaN', N, A, b, C, v, L, x0, thin, test)
}
......@@ -7,6 +7,7 @@
#' @param nchain the number of mcmc chains
#' @param ncore number of cores to use
#' @param thin thinning interval
#' @param test if TRUE try a new method aiming at reducing autocorrelation
#' @return a \code{\link[coda]{mcmc.list}}
#' @export
#'
......@@ -41,7 +42,8 @@ fitmyCaNmod <- function(myCaNmod,
N,
nchain = 1,
ncore = 1,
thin = 1) {
thin = 1,
test = FALSE) {
ncore <- min(min(detectCores() - 1, ncore), nchain)
`%myinfix%` <- `%do%`
......@@ -91,7 +93,8 @@ fitmyCaNmod <- function(myCaNmod,
as.matrix(myCaNmod$C),
myCaNmod$v,
as.matrix(myCaNmod$L),
x0
x0,
test
)
names(res) <- c("F", "B")
res$F <- res$F[, -seq_len(length(myCaNmod$species))]
......
......@@ -4,7 +4,8 @@
\alias{fitmyCaNmod}
\title{fitmyCaNmod}
\usage{
fitmyCaNmod(myCaNmod, N, nchain = 1, ncore = 1, thin = 1)
fitmyCaNmod(myCaNmod, N, nchain = 1, ncore = 1, thin = 1,
test = FALSE)
}
\arguments{
\item{myCaNmod}{a CaNmod object with following elements}
......@@ -16,6 +17,8 @@ fitmyCaNmod(myCaNmod, N, nchain = 1, ncore = 1, thin = 1)
\item{ncore}{number of cores to use}
\item{thin}{thinning interval}
\item{test}{if TRUE try a new method aiming at reducing autocorrelation}
}
\value{
a \code{\link[coda]{mcmc.list}}
......
......@@ -7,8 +7,8 @@
using namespace Rcpp;
// fitCaN
List fitCaN(const int N, const Eigen::MatrixXd& A, const Eigen::VectorXd& b, const Eigen::MatrixXd& C, const Eigen::VectorXd& v, const Eigen::MatrixXd& L, const Eigen::VectorXd& x0, const int thin);
RcppExport SEXP _RCaN_fitCaN(SEXP NSEXP, SEXP ASEXP, SEXP bSEXP, SEXP CSEXP, SEXP vSEXP, SEXP LSEXP, SEXP x0SEXP, SEXP thinSEXP) {
List fitCaN(const int N, const Eigen::MatrixXd& A, const Eigen::VectorXd& b, const Eigen::MatrixXd& C, const Eigen::VectorXd& v, const Eigen::MatrixXd& L, const Eigen::VectorXd& x0, const int thin, const bool test);
RcppExport SEXP _RCaN_fitCaN(SEXP NSEXP, SEXP ASEXP, SEXP bSEXP, SEXP CSEXP, SEXP vSEXP, SEXP LSEXP, SEXP x0SEXP, SEXP thinSEXP, SEXP testSEXP) {
BEGIN_RCPP
Rcpp::RObject rcpp_result_gen;
Rcpp::RNGScope rcpp_rngScope_gen;
......@@ -20,13 +20,14 @@ BEGIN_RCPP
Rcpp::traits::input_parameter< const Eigen::MatrixXd& >::type L(LSEXP);
Rcpp::traits::input_parameter< const Eigen::VectorXd& >::type x0(x0SEXP);
Rcpp::traits::input_parameter< const int >::type thin(thinSEXP);
rcpp_result_gen = Rcpp::wrap(fitCaN(N, A, b, C, v, L, x0, thin));
Rcpp::traits::input_parameter< const bool >::type test(testSEXP);
rcpp_result_gen = Rcpp::wrap(fitCaN(N, A, b, C, v, L, x0, thin, test));
return rcpp_result_gen;
END_RCPP
}
static const R_CallMethodDef CallEntries[] = {
{"_RCaN_fitCaN", (DL_FUNC) &_RCaN_fitCaN, 8},
{"_RCaN_fitCaN", (DL_FUNC) &_RCaN_fitCaN, 9},
{NULL, NULL, 0}
};
......
......@@ -18,15 +18,16 @@ using Eigen::VectorXd;
List fitCaN(const int N, const Eigen::MatrixXd &A ,const Eigen::VectorXd &b,
const Eigen::MatrixXd &C ,const Eigen::VectorXd &v,
const Eigen::MatrixXd &L,
const Eigen::VectorXd &x0, const int thin) {
const Eigen::VectorXd &x0, const int thin,
const bool test) {
int p=A.cols();
int m2=C.rows();
MatrixXd F(N, p);
MatrixXd B(N, L.rows());
if(m2>0){ //there are equality constraints
F=cpgsR::cpgsEquality(N, A, b, C, v, x0, thin);
F=cpgsR::cpgsEquality(N, A, b, C, v, x0, thin, test);
} else{
F=cpgsR::cpgs(N, A, b, x0, thin);
F=cpgsR::cpgs(N, A, b, x0, thin, test);
}
for (int i=0;i<N;++i){
B.row(i)=L*F.row(i).transpose();
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment