diff options
author | Ken Kellner <ken@kenkellner.com> | 2023-07-06 11:14:37 -0400 |
---|---|---|
committer | Ken Kellner <ken@kenkellner.com> | 2023-07-06 11:14:52 -0400 |
commit | cd4783a0ee6aa9e17c7fe6e9901e780b4bf9407c (patch) | |
tree | 10eac548b7027049b77d340a3966e30097761f29 | |
parent | 5c9f793a19bbd8d79a36c1ebd9a3f38df88dcd14 (diff) |
Finish cpp updates and bump version. Fixes #258.
-rw-r--r-- | DESCRIPTION | 2 | ||||
-rw-r--r-- | NEWS.md | 6 | ||||
-rw-r--r-- | R/RcppExports.R | 8 | ||||
-rw-r--r-- | R/occuMulti.R | 14 | ||||
-rw-r--r-- | src/RcppExports.cpp | 57 | ||||
-rw-r--r-- | src/nll_occuMulti.cpp | 73 | ||||
-rw-r--r-- | src/nll_occuMulti.h | 11 | ||||
-rw-r--r-- | tests/testthat/test_occuMulti.R | 5 |
8 files changed, 106 insertions, 70 deletions
diff --git a/DESCRIPTION b/DESCRIPTION index 9800624..ba0d28e 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -1,5 +1,5 @@ Package: unmarked -Version: 1.3.1.9001 +Version: 1.3.2 Date: 2023-07-06 Type: Package Title: Models for Data from Unmarked Animals @@ -1,4 +1,8 @@ -# unmarked 1.3.1 +# unmarked 1.3.2 + +* Modernize some Cpp code to pass new LTO checks + +#unmarked 1.3.1 * Remove log.grad function to pass CRAN checks diff --git a/R/RcppExports.R b/R/RcppExports.R index c6aa6d2..79ab967 100644 --- a/R/RcppExports.R +++ b/R/RcppExports.R @@ -53,6 +53,14 @@ nll_occuMS <- function(beta, y, dm_state, dm_phi, dm_det, sind, pind, dind, prm, .Call(`_unmarked_nll_occuMS`, beta, y, dm_state, dm_phi, dm_det, sind, pind, dind, prm, S, T, J, N, naflag, guide) } +nll_occuMulti_loglik <- function(fStart, fStop, dmF, dmOcc, beta, dmDet, dStart, dStop, y, yStart, yStop, Iy0, z, fixed0) { + .Call(`_unmarked_nll_occuMulti_loglik`, fStart, fStop, dmF, dmOcc, beta, dmDet, dStart, dStop, y, yStart, yStop, Iy0, z, fixed0) +} + +nll_occuMulti <- function(fStart, fStop, dmF, dmOcc, beta, dmDet, dStart, dStop, y, yStart, yStop, Iy0, z, fixed0, penalty) { + .Call(`_unmarked_nll_occuMulti`, fStart, fStop, dmF, dmOcc, beta, dmDet, dStart, dStop, y, yStart, yStop, Iy0, z, fixed0, penalty) +} + nll_occuPEN <- function(y, X, V, beta_psi, beta_p, nd, knownOcc, navec, X_offset, V_offset, penalty) { .Call(`_unmarked_nll_occuPEN`, y, X, V, beta_psi, beta_p, nd, knownOcc, navec, X_offset, V_offset, penalty) } diff --git a/R/occuMulti.R b/R/occuMulti.R index 8b20dd1..7d31d47 100644 --- a/R/occuMulti.R +++ b/R/occuMulti.R @@ -88,10 +88,10 @@ occuMulti <- function(detformulas, stateformulas, data, maxOrder, #Likelihood function in C---------------------------------------------------- nll_C <- function(params) { - .Call("nll_occuMulti", + nll_occuMulti( fStart-1, fStop-1, t_dmF, dmOcc, params, dmDet, dStart-1, dStop-1, - y, yStart-1, yStop-1, Iy0, as.matrix(z), fixed0, penalty, 0, - PACKAGE = "unmarked") + y, yStart-1, yStop-1, Iy0, as.matrix(z), fixed0, penalty + ) } #---------------------------------------------------------------------------- @@ -169,13 +169,11 @@ occuMultiLogLik <- function(fit, data){ dmF <- Matrix::Matrix(dm$dmF, sparse=TRUE) t_dmF <- Matrix::t(dmF) - out <- .Call("nll_occuMulti", + out <- nll_occuMulti_loglik( dm$fStart-1, dm$fStop-1, t_dmF, dm$dmOcc, coef(fit), dm$dmDet, dm$dStart-1, dm$dStop-1, dm$y, - dm$yStart-1, dm$yStop-1, dm$Iy0, as.matrix(dm$z), dm$fixed0, 0, - #return site likelihoods - 1, - PACKAGE = "unmarked") + dm$yStart-1, dm$yStop-1, dm$Iy0, as.matrix(dm$z), dm$fixed0 + ) as.vector(out) diff --git a/src/RcppExports.cpp b/src/RcppExports.cpp index c69911a..6d89809 100644 --- a/src/RcppExports.cpp +++ b/src/RcppExports.cpp @@ -363,6 +363,55 @@ BEGIN_RCPP return rcpp_result_gen; END_RCPP } +// nll_occuMulti_loglik +arma::vec nll_occuMulti_loglik(Rcpp::IntegerVector fStart, Rcpp::IntegerVector fStop, arma::sp_mat dmF, Rcpp::List dmOcc, arma::colvec beta, Rcpp::List dmDet, Rcpp::IntegerVector dStart, Rcpp::IntegerVector dStop, arma::mat y, Rcpp::IntegerVector yStart, Rcpp::IntegerVector yStop, arma::mat Iy0, arma::mat z, Rcpp::LogicalVector fixed0); +RcppExport SEXP _unmarked_nll_occuMulti_loglik(SEXP fStartSEXP, SEXP fStopSEXP, SEXP dmFSEXP, SEXP dmOccSEXP, SEXP betaSEXP, SEXP dmDetSEXP, SEXP dStartSEXP, SEXP dStopSEXP, SEXP ySEXP, SEXP yStartSEXP, SEXP yStopSEXP, SEXP Iy0SEXP, SEXP zSEXP, SEXP fixed0SEXP) { +BEGIN_RCPP + Rcpp::RObject rcpp_result_gen; + Rcpp::RNGScope rcpp_rngScope_gen; + Rcpp::traits::input_parameter< Rcpp::IntegerVector >::type fStart(fStartSEXP); + Rcpp::traits::input_parameter< Rcpp::IntegerVector >::type fStop(fStopSEXP); + Rcpp::traits::input_parameter< arma::sp_mat >::type dmF(dmFSEXP); + Rcpp::traits::input_parameter< Rcpp::List >::type dmOcc(dmOccSEXP); + Rcpp::traits::input_parameter< arma::colvec >::type beta(betaSEXP); + Rcpp::traits::input_parameter< Rcpp::List >::type dmDet(dmDetSEXP); + Rcpp::traits::input_parameter< Rcpp::IntegerVector >::type dStart(dStartSEXP); + Rcpp::traits::input_parameter< Rcpp::IntegerVector >::type dStop(dStopSEXP); + Rcpp::traits::input_parameter< arma::mat >::type y(ySEXP); + Rcpp::traits::input_parameter< Rcpp::IntegerVector >::type yStart(yStartSEXP); + Rcpp::traits::input_parameter< Rcpp::IntegerVector >::type yStop(yStopSEXP); + Rcpp::traits::input_parameter< arma::mat >::type Iy0(Iy0SEXP); + Rcpp::traits::input_parameter< arma::mat >::type z(zSEXP); + Rcpp::traits::input_parameter< Rcpp::LogicalVector >::type fixed0(fixed0SEXP); + rcpp_result_gen = Rcpp::wrap(nll_occuMulti_loglik(fStart, fStop, dmF, dmOcc, beta, dmDet, dStart, dStop, y, yStart, yStop, Iy0, z, fixed0)); + return rcpp_result_gen; +END_RCPP +} +// nll_occuMulti +double nll_occuMulti(Rcpp::IntegerVector fStart, Rcpp::IntegerVector fStop, arma::sp_mat dmF, Rcpp::List dmOcc, arma::colvec beta, Rcpp::List dmDet, Rcpp::IntegerVector dStart, Rcpp::IntegerVector dStop, arma::mat y, Rcpp::IntegerVector yStart, Rcpp::IntegerVector yStop, arma::mat Iy0, arma::mat z, Rcpp::LogicalVector fixed0, double penalty); +RcppExport SEXP _unmarked_nll_occuMulti(SEXP fStartSEXP, SEXP fStopSEXP, SEXP dmFSEXP, SEXP dmOccSEXP, SEXP betaSEXP, SEXP dmDetSEXP, SEXP dStartSEXP, SEXP dStopSEXP, SEXP ySEXP, SEXP yStartSEXP, SEXP yStopSEXP, SEXP Iy0SEXP, SEXP zSEXP, SEXP fixed0SEXP, SEXP penaltySEXP) { +BEGIN_RCPP + Rcpp::RObject rcpp_result_gen; + Rcpp::RNGScope rcpp_rngScope_gen; + Rcpp::traits::input_parameter< Rcpp::IntegerVector >::type fStart(fStartSEXP); + Rcpp::traits::input_parameter< Rcpp::IntegerVector >::type fStop(fStopSEXP); + Rcpp::traits::input_parameter< arma::sp_mat >::type dmF(dmFSEXP); + Rcpp::traits::input_parameter< Rcpp::List >::type dmOcc(dmOccSEXP); + Rcpp::traits::input_parameter< arma::colvec >::type beta(betaSEXP); + Rcpp::traits::input_parameter< Rcpp::List >::type dmDet(dmDetSEXP); + Rcpp::traits::input_parameter< Rcpp::IntegerVector >::type dStart(dStartSEXP); + Rcpp::traits::input_parameter< Rcpp::IntegerVector >::type dStop(dStopSEXP); + Rcpp::traits::input_parameter< arma::mat >::type y(ySEXP); + Rcpp::traits::input_parameter< Rcpp::IntegerVector >::type yStart(yStartSEXP); + Rcpp::traits::input_parameter< Rcpp::IntegerVector >::type yStop(yStopSEXP); + Rcpp::traits::input_parameter< arma::mat >::type Iy0(Iy0SEXP); + Rcpp::traits::input_parameter< arma::mat >::type z(zSEXP); + Rcpp::traits::input_parameter< Rcpp::LogicalVector >::type fixed0(fixed0SEXP); + Rcpp::traits::input_parameter< double >::type penalty(penaltySEXP); + rcpp_result_gen = Rcpp::wrap(nll_occuMulti(fStart, fStop, dmF, dmOcc, beta, dmDet, dStart, dStop, y, yStart, yStop, Iy0, z, fixed0, penalty)); + return rcpp_result_gen; +END_RCPP +} // nll_occuPEN double nll_occuPEN(arma::icolvec y, arma::mat X, arma::mat V, arma::colvec beta_psi, arma::colvec beta_p, Rcpp::IntegerVector nd, Rcpp::LogicalVector knownOcc, Rcpp::LogicalVector navec, arma::colvec X_offset, arma::colvec V_offset, double penalty); RcppExport SEXP _unmarked_nll_occuPEN(SEXP ySEXP, SEXP XSEXP, SEXP VSEXP, SEXP beta_psiSEXP, SEXP beta_pSEXP, SEXP ndSEXP, SEXP knownOccSEXP, SEXP navecSEXP, SEXP X_offsetSEXP, SEXP V_offsetSEXP, SEXP penaltySEXP) { @@ -500,7 +549,6 @@ END_RCPP RcppExport SEXP getDetVecs(void *, void *, void *, void *, void *); RcppExport SEXP getSingleDetVec(void *, void *, void *); -RcppExport SEXP nll_occuMulti(void *, void *, void *, void *, void *, void *, void *, void *, void *, void *, void *, void *, void *, void *, void *, void *); static const R_CallMethodDef CallEntries[] = { {"_unmarked_get_lik_trans", (DL_FUNC) &_unmarked_get_lik_trans, 2}, @@ -516,14 +564,15 @@ static const R_CallMethodDef CallEntries[] = { {"_unmarked_nll_nmixTTD", (DL_FUNC) &_unmarked_nll_nmixTTD, 13}, {"_unmarked_nll_occu", (DL_FUNC) &_unmarked_nll_occu, 11}, {"_unmarked_nll_occuMS", (DL_FUNC) &_unmarked_nll_occuMS, 15}, + {"_unmarked_nll_occuMulti_loglik", (DL_FUNC) &_unmarked_nll_occuMulti_loglik, 14}, + {"_unmarked_nll_occuMulti", (DL_FUNC) &_unmarked_nll_occuMulti, 15}, {"_unmarked_nll_occuPEN", (DL_FUNC) &_unmarked_nll_occuPEN, 11}, {"_unmarked_nll_occuRN", (DL_FUNC) &_unmarked_nll_occuRN, 10}, {"_unmarked_nll_occuTTD", (DL_FUNC) &_unmarked_nll_occuTTD, 17}, {"_unmarked_nll_pcount", (DL_FUNC) &_unmarked_nll_pcount, 11}, {"_unmarked_nll_pcountOpen", (DL_FUNC) &_unmarked_nll_pcountOpen, 35}, - {"getDetVecs", (DL_FUNC) &getDetVecs, 5}, - {"getSingleDetVec", (DL_FUNC) &getSingleDetVec, 3}, - {"nll_occuMulti", (DL_FUNC) &nll_occuMulti, 16}, + {"getDetVecs", (DL_FUNC) &getDetVecs, 5}, + {"getSingleDetVec", (DL_FUNC) &getSingleDetVec, 3}, {NULL, NULL, 0} }; diff --git a/src/nll_occuMulti.cpp b/src/nll_occuMulti.cpp index 3a44c7d..5310097 100644 --- a/src/nll_occuMulti.cpp +++ b/src/nll_occuMulti.cpp @@ -1,46 +1,21 @@ -#include "nll_occuMulti.h" +#include <RcppArmadillo.h> using namespace Rcpp; using namespace arma; -SEXP nll_occuMulti( SEXP fStartR, SEXP fStopR, SEXP dmFr, SEXP dmOccR, - SEXP betaR, SEXP dmDetR, SEXP dStartR, SEXP dStopR, SEXP yR, SEXP yStartR, - SEXP yStopR, SEXP Iy0r, SEXP zR, SEXP fixed0r, SEXP penaltyR, - SEXP returnLLr){ - - //Inputs - IntegerVector fStart(fStartR); - IntegerVector fStop(fStopR); - - //if Matrix is a dependency - sp_mat dmF = as<sp_mat>(dmFr); //already transposed - - //if Matrix not a dependency - //sp_mat dmF( as<mat>(dmFr) ); - - int nF = dmF.n_rows; - List dmOcc(dmOccR); - - colvec beta = as<colvec>(betaR); - LogicalVector fixed0(fixed0r); - - List dmDet(dmDetR); - IntegerVector dStart(dStartR); - IntegerVector dStop(dStopR); - - mat y = as<mat>(yR); +// [[Rcpp::export]] +arma::vec nll_occuMulti_loglik(Rcpp::IntegerVector fStart, Rcpp::IntegerVector fStop, + arma::sp_mat dmF, Rcpp::List dmOcc, + arma::colvec beta, Rcpp::List dmDet, + Rcpp::IntegerVector dStart, Rcpp::IntegerVector dStop, + arma::mat y, Rcpp::IntegerVector yStart, Rcpp::IntegerVector yStop, + arma::mat Iy0, arma::mat z, Rcpp::LogicalVector fixed0){ + + int nF = dmF.n_rows; //dmF is already transposed int S = y.n_cols; int J = y.n_rows; - IntegerVector yStart(yStartR); - IntegerVector yStop(yStopR); int N = yStart.size(); - mat Iy0 = as<mat>(Iy0r); - - mat z = as<mat>(zR); - double penalty = as<double>(penaltyR); - int returnLL = as<int>(returnLLr); - //psi calculation int index = 0; mat f(N, nF); @@ -54,7 +29,7 @@ SEXP nll_occuMulti( SEXP fStartR, SEXP fStopR, SEXP dmFr, SEXP dmOccR, } } - mat psi = exp( f * dmF ); + mat psi = exp( f * dmF ); for(unsigned int i = 0; i < psi.n_rows; i++){ psi.row(i) = psi.row(i) / sum( psi.row(i) ); } @@ -72,15 +47,15 @@ SEXP nll_occuMulti( SEXP fStartR, SEXP fStopR, SEXP dmFr, SEXP dmOccR, vec logLik(N); for(int i = 0; i < N; i++){ - + mat ysub = y.rows(yStart[i], yStop[i]); mat psub = p.rows(yStart[i], yStop[i]); rowvec cdp(S); for(int j = 0; j < S; j++){ - cdp(j) = exp( sum( ysub.col(j) % log(psub.col(j)) ) + + cdp(j) = exp( sum( ysub.col(j) % log(psub.col(j)) ) + sum( (1 - ysub.col(j)) % log( 1 - psub.col(j)) ) ); } - + rowvec prdProbY(M); for(int j = 0; j < M; j++){ prdProbY(j) = prod( z.row(j) % cdp + (1 - z.row(j)) % Iy0.row(i) ); @@ -89,11 +64,21 @@ SEXP nll_occuMulti( SEXP fStartR, SEXP fStopR, SEXP dmFr, SEXP dmOccR, logLik(i) = log( sum( psi.row(i) % prdProbY ) ); } - - if(returnLL){ - return(wrap(logLik)); - } + + return logLik; +} + +// [[Rcpp::export]] +double nll_occuMulti(Rcpp::IntegerVector fStart, Rcpp::IntegerVector fStop, + arma::sp_mat dmF, Rcpp::List dmOcc, + arma::colvec beta, Rcpp::List dmDet, + Rcpp::IntegerVector dStart, Rcpp::IntegerVector dStop, + arma::mat y, Rcpp::IntegerVector yStart, Rcpp::IntegerVector yStop, + arma::mat Iy0, arma::mat z, Rcpp::LogicalVector fixed0, double penalty){ + + vec logLik = nll_occuMulti_loglik(fStart, fStop, dmF, dmOcc, beta, dmDet, + dStart, dStop, y, yStart, yStop, Iy0, z, fixed0); double pen = penalty * 0.5 * accu(pow(beta, 2)); - return(wrap(-1.0 * (sum(logLik) - pen))); + return -1.0 * (sum(logLik) - pen); } diff --git a/src/nll_occuMulti.h b/src/nll_occuMulti.h deleted file mode 100644 index 47cb53e..0000000 --- a/src/nll_occuMulti.h +++ /dev/null @@ -1,11 +0,0 @@ -#ifndef _unmarked_NLL_OCCUMULTI_H -#define _unmarked_NLL_OCCUMULTI_H - -#include <RcppArmadillo.h> - -RcppExport SEXP nll_occuMulti( SEXP fStartR, SEXP fStopR, SEXP dmFr, SEXP dmOccR, - SEXP betaR, SEXP dmDetR, SEXP dStartR, SEXP dStopR, SEXP yR, SEXP yStartR, - SEXP yStopR, SEXP Iy0r, SEXP zR, SEXP fixed0r, SEXP penaltyR, - SEXP returnLLr) ; - -#endif diff --git a/tests/testthat/test_occuMulti.R b/tests/testthat/test_occuMulti.R index d37cbc5..6ce7921 100644 --- a/tests/testthat/test_occuMulti.R +++ b/tests/testthat/test_occuMulti.R @@ -536,7 +536,10 @@ test_that("occuMulti penalized likelihood works",{ set.seed(123) opt_fit <- optimizePenalty(fm, penalties=c(0,1), boot=2) expect_equal(opt_fit@call$penalty, 1) - + + ll <- unmarked:::occuMultiLogLik(fm, fm@data) + expect_equal(length(ll), numSites(fm@data)) + expect_equal(-sum(ll), fm@negLogLike) }) |