-
Notifications
You must be signed in to change notification settings - Fork 4
Commit
This commit does not belong to any branch on this repository, and may belong to a fork outside of the repository.
* small updates * Update vignette.Rmd * attempt at residuals for family="fixed" * fix docs * update to allow DHARMa using predictive samples * update citation * fix bug in loo_residuals * fix CI * add track_progress argument to loo_residuals * Update vignette.Rmd * add upper and lower bounds for optimizer * fix residuals demo * Update vignette.Rmd * update plots * fix last commit * Update dsem.R * fit error in DHARMa vignette example * add `stepwise_selection(.)` * fix notes and warnings * add `equations_to_text` * update vignette to show `convert_equations` * update DFA analysis using lag-and-equations syntax * fix bug in docs * update make_dsem_ram to use additive-kronecker * try to fix CI error * fix attempt-2 * attempt-3 * change DFA vignette back * next attempt * another fix attempt * update citation and pkgdown yaml * adding cAIC ... ... which only works when using models with measurement errors * merging add_priors to dev branch * fix merge * add test-priors CI check * test fix for R-CMD-check error * add dsemRTMB * add reduced-rank GMRF to dsemRTMB * adding log_prior argument to dsemRTMB * Update test_dsemRTMB.R * finish initial dsemRTMB * get dsemRTMB working as package * allow vector-valued priors * more simulator scratch-scripts * adding obj$simulator(.) slot * fix dsemRTMB to work without loading locally * add run_time slot (to compare TMB and RTMB) * fix error in sparseMatrix stuff * updates for merge * try making CI happy * another attempt at R_CMD_check * another try * another try to R_CMD_check * cautious update * make check_win_devel(.) happy * add CI for fixed and mapped parameters and fix issue in dsemRTMB with these * make check_win_devel(.) happy * Update make_matrices.R --------- Co-authored-by: Jim Thorson <James.T.Thorson@gmail.com>
- Loading branch information
1 parent
298cbb7
commit 0876fa8
Showing
43 changed files
with
2,765 additions
and
116 deletions.
There are no files selected for viewing
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| Original file line number | Diff line number | Diff line change |
|---|---|---|
| @@ -0,0 +1,122 @@ | ||
|
|
||
| #' @title Calculate conditional AIC | ||
| #' | ||
| #' @description | ||
| #' Calculates the conditional Akaike Information criterion (cAIC). | ||
| #' | ||
| #' @param object Output from \code{\link{dsem}} | ||
| #' @param what Whether to return the cAIC or the effective degrees of freedom | ||
| #' (EDF) for each group of random effects. | ||
| #' | ||
| #' @details | ||
| #' cAIC is designed to optimize the expected out-of-sample predictive | ||
| #' performance for new data that share the same random effects as the | ||
| #' in-sample (fitted) data, e.g., spatial interpolation. In this sense, | ||
| #' it should be a fast approximation to optimizing the model structure | ||
| #' based on k-fold crossvalidation. | ||
| #' By contrast, \code{AIC} calculates the | ||
| #' marginal Akaike Information Criterion, which is designed to optimize | ||
| #' expected predictive performance for new data that have new random effects, | ||
| #' e.g., extrapolation, or inference about generative parameters. | ||
| #' | ||
| #' cAIC also calculates as a byproduct the effective degrees of freedom, | ||
| #' i.e., the number of fixed effects that would have an equivalent impact on | ||
| #' model flexibility as a given random effect. | ||
| #' | ||
| #' Both cAIC and EDF are calculated using Eq. 6 of Zheng Cadigan Thorson 2024. | ||
| #' | ||
| #' Note that, for models that include profiled fixed effects, these profiles | ||
| #' are turned off. | ||
| #' | ||
| #' @return | ||
| #' Either the cAIC, or the effective degrees of freedom (EDF) by group | ||
| #' of random effects | ||
| #' | ||
| #' @references | ||
| #' | ||
| #' **Deriving the general approximation to cAIC used here** | ||
| #' | ||
| #' Zheng, N., Cadigan, N., & Thorson, J. T. (2024). | ||
| #' A note on numerical evaluation of conditional Akaike information for | ||
| #' nonlinear mixed-effects models (arXiv:2411.14185). arXiv. | ||
| #' \doi{10.48550/arXiv.2411.14185} | ||
| #' | ||
| #' **The utility of EDF to diagnose hierarchical model behavior** | ||
| #' | ||
| #' Thorson, J. T. (2024). Measuring complexity for hierarchical | ||
| #' models using effective degrees of freedom. Ecology, | ||
| #' 105(7), e4327 \doi{10.1002/ecy.4327} | ||
| #' | ||
| #' @export | ||
| cAIC <- | ||
| function( object, | ||
| what = c("cAIC","EDF") ){ | ||
|
|
||
| what = match.arg(what) | ||
| data = object$tmb_inputs$data | ||
|
|
||
| # Error checks | ||
| if(any(is.na(object$tmb_inputs$map$x_tj))){ | ||
| stop("cAIC is not implemented when fixing states at data using family=`fixed`") | ||
| } | ||
|
|
||
| # Turn on all GMRF parameters | ||
| map = object$tmb_inputs$map | ||
| map$x_tj = factor(seq_len(prod(dim(data$y_tj)))) | ||
|
|
||
| # Make sure profile = NULL | ||
| #if( is.null(object$internal$control$profile) ){ | ||
| obj = object$obj | ||
| #}else{ | ||
| obj = TMB::MakeADFun( data = data, | ||
| parameters = object$internal$parhat, | ||
| random = object$tmb_inputs$random, | ||
| map = map, | ||
| profile = NULL, | ||
| DLL="dsem", | ||
| silent = TRUE ) | ||
| #} | ||
|
|
||
| # Weights = 0 is equivalent to data = NA | ||
| data$y_tj[] = NA | ||
| # Make obj_new | ||
| obj_new = TMB::MakeADFun( data = data, | ||
| parameters = object$internal$parhat, | ||
| map = map, | ||
| random = object$tmb_inputs$random, | ||
| DLL = "dsem", | ||
| profile = NULL ) | ||
|
|
||
| # | ||
| par = obj$env$parList() | ||
| parDataMode <- obj$env$last.par | ||
| indx = obj$env$lrandom() | ||
| q = sum(indx) | ||
| p = length(object$opt$par) | ||
|
|
||
| ## use - for Hess because model returns negative loglikelihood; | ||
| Hess_new = -Matrix::Matrix(obj_new$env$f(parDataMode,order=1,type="ADGrad"),sparse = TRUE) | ||
| #cov_Psi_inv = -Hess_new[indx,indx]; ## this is the marginal prec mat of REs; | ||
| Hess_new = Hess_new[indx,indx] | ||
|
|
||
| ## Joint hessian etc | ||
| Hess = -Matrix::Matrix(obj$env$f(parDataMode,order=1,type="ADGrad"),sparse = TRUE) | ||
| Hess = Hess[indx,indx] | ||
| #negEDF = diag(as.matrix(solve(ddlj.r)) %*% ddlr.r) | ||
| negEDF = Matrix::diag(Matrix::solve(Hess, Hess_new)) | ||
| # | ||
|
|
||
| if(what=="cAIC"){ | ||
| jnll = obj$env$f(parDataMode) | ||
| cnll = jnll - obj_new$env$f(parDataMode) | ||
| cAIC = 2*cnll + 2*(p+q) - 2*sum(negEDF) | ||
| return(cAIC) | ||
| } | ||
| if(what=="EDF"){ | ||
| #Sdims = object$tmb_inputs$tmb_data$Sdims | ||
| #group = rep.int( seq_along(Sdims), times=Sdims ) | ||
| #names(negEDF) = names(obj$env$last.par)[indx] | ||
| EDF = length(negEDF) - sum(negEDF) | ||
| return(EDF) | ||
| } | ||
| } |
Oops, something went wrong.