Skip to content

Commit

Permalink
Merge branch 'kellijohnson-NOAA-check-plot-aniso'
Browse files Browse the repository at this point in the history
  • Loading branch information
seananderson committed Apr 24, 2023
2 parents b7527ae + 57ac4ca commit f5700b2
Show file tree
Hide file tree
Showing 3 changed files with 30 additions and 4 deletions.
22 changes: 20 additions & 2 deletions R/plot.R
Original file line number Diff line number Diff line change
Expand Up @@ -15,11 +15,14 @@
#' anisotropy. The ellipses are centered at coordinates of zero in the space of
#' the X-Y coordinates being modeled. The ellipses show the spatial and/or
#' spatiotemporal range (distance at which correlation is effectively
#' independent) in any direction from zero. Uses \pkg{ggplot2}.
#' independent) in any direction from zero. Uses \pkg{ggplot2}. If anisotropy
#' was turned off when fitting the model, `NULL` is returned instead of a
#' {ggplot2} object.
#'
#' `plot_anisotropy2()`: A plot of eigenvectors illustrating the estimated
#' anisotropy. A list of the plotted data is invisibly returned. Uses base
#' graphics.
#' graphics. If anisotropy was turned off when fitting the model, `NULL` is
#' returned instead of a plot object.
#' @references Code adapted from VAST R package
#' @importFrom rlang .data
#' @examplesIf inla_installed() && ggplot2_installed()
Expand All @@ -40,6 +43,7 @@
#' @rdname plot_anisotropy
plot_anisotropy <- function(object, return_data = FALSE) {
stopifnot(inherits(object, "sdmTMB"))
if (!check_for_H(object)) return(NULL)
report <- object$tmb_obj$report(object$tmb_obj$env$last.par.best)
delta <- isTRUE(object$family$delta)

Expand Down Expand Up @@ -156,6 +160,7 @@ plot_anisotropy <- function(object, return_data = FALSE) {
#' @rdname plot_anisotropy
plot_anisotropy2 <- function(object, model = 1) {
stopifnot(inherits(object, "sdmTMB"))
if (!check_for_H(object)) return(NULL)
report <- object$tmb_obj$report(object$tmb_obj$env$last.par.best)
if (model == 1) eig <- eigen(report$H)
if (model == 2) eig <- eigen(report$H2)
Expand Down Expand Up @@ -327,3 +332,16 @@ plot_smooth <- function(object, select = 1, n = 100, level = 0.95,
return(g)
}
}

check_for_H <- function(obj) {
H <- any(grepl(
pattern = "ln_H_input",
x = names(obj$sd_report$par.fixed),
ignore.case = TRUE
))
if (!H) {
cli::cli_inform("`anisotropy = FALSE` in `sdmTMB()`; no anisotropy figure is available.")
# FIXME in the future plot the isotropic covariance instead of NULL?
}
H
}
7 changes: 5 additions & 2 deletions man/plot_anisotropy.Rd

Some generated files are not rendered by default. Learn more about how customized files appear on GitHub.

5 changes: 5 additions & 0 deletions tests/testthat/test-2-print-anisotropy.R
Original file line number Diff line number Diff line change
Expand Up @@ -13,6 +13,8 @@ test_that("Print anisotropy prints correctly", {
)

expect_output(print(fit1), regexp = "range: 33.23")
expect_null(plot_anisotropy(fit1))
expect_null(plot_anisotropy2(fit1))

# -------------------
# Anisotropy with spatial only
Expand All @@ -26,6 +28,8 @@ test_that("Print anisotropy prints correctly", {
anisotropy = TRUE
)

plot_anisotropy(fit_sp_only)
plot_anisotropy2(fit_sp_only)
expect_output(print(fit_sp_only), regexp = "\\(spatial\\): 6.1 to 86.0 at 126")

# Anisotropy with only spatiotemporal random field
Expand Down Expand Up @@ -92,3 +96,4 @@ test_that("Print anisotropy prints correctly", {
expect_output(cat(print_anisotropy(fit_dg_not_shared, m = 2)), regexp = "\\(spatial\\): 0")
expect_output(cat(print_anisotropy(fit_dg_not_shared, m = 2)), regexp = "\\(spatiotemporal\\): 9")
})

0 comments on commit f5700b2

Please sign in to comment.