Become an expert in R — Interactive courses, Cheat Sheets, certificates and more!
Get Started for Free

log_rss

Calculate log-RSS for a fitted model


Description

Calculate log-RSS(x1, x2) for a fitted RSF or (i)SSF

Usage

log_rss(object, ...)

## S3 method for class 'glm'
log_rss(object, x1, x2, ci = NA, ci_level = 0.95, n_boot = 1000, ...)

## S3 method for class 'fit_clogit'
log_rss(object, x1, x2, ci = NA, ci_level = 0.95, n_boot = 1000, ...)

Arguments

object

[fit_logit, fit_clogit]
A fitted RSF or (i)SSF model.

...

Further arguments, none implemented.

x1

[data.frame]
A data.frame representing the habitat values at location x_1. Must contain all fitted covariates as expected by predict().

x2

[data.frame]
A 1-row data.frame representing the single hypothetical location of x_2. Must contain all fitted covariates as expected by predict().

ci

[character]
Method for estimating confidence intervals around log-RSS. NA skips calculating CIs. Character string "se" uses standard error method and "boot" uses empirical bootstrap method.

ci_level

[numeric]
Level for confidence interval. Defaults to 0.95 for a 95% confidence interval.

n_boot

[integer]
Number of bootstrap samples to estimate confidence intervals. Ignored if ci != "boot".

Details

This function assumes that the user would like to compare relative selection strengths from at least one proposed location (x1) to exactly one reference location (x2).

The objects object$model, x1, and x2 will be passed to predict(). Therefore, the columns of x1 and x2 must match the terms in the model formula exactly.

Value

Returns a list of class log_rss.

Author(s)

Brian J. Smith

References

Avgar, T., Lele, S.R., Keim, J.L., and Boyce, M.S.. (2017). Relative Selection Strength: Quantifying effect size in habitat- and step-selection inference. Ecology and Evolution, 7, 5322–5330.

See Also

See Avgar et al. 2017 for details about relative selection strength.

Default plotting method available: plot.log_rss()

Examples

## Not run: 
# RSF -------------------------------------------------------
# Fit an RSF, then calculate log-RSS to visualize results.

# Load packages
library(ggplot2)

#Load data
data("amt_fisher")

# Prepare data for RSF
rsf_data <- amt_fisher %>%
  filter(name == "Lupe") %>%
  make_track(x_, y_, t_) %>%
  random_points() %>%
  extract_covariates(amt_fisher_covar$elevation) %>%
  extract_covariates(amt_fisher_covar$popden) %>%
  extract_covariates(amt_fisher_covar$landuse) %>%
  mutate(lu = factor(landuse))

# Fit RSF
m1 <- rsf_data %>%
  fit_rsf(case_ ~ lu + elevation + popden)

# Calculate log-RSS
# data.frame of x1s
x1 <- data.frame(lu = factor(50, levels = levels(rsf_data$lu)),
                 elevation = seq(90, 120, length.out = 100),
                 popden = mean(rsf_data$popden))
# data.frame of x2 (note factor levels should be same as model data)
x2 <- data.frame(lu = factor(50, levels = levels(rsf_data$lu)),
                 elevation = mean(rsf_data$elevation),
                 popden = mean(rsf_data$popden))
# Calculate (use se method for confidence interval)
logRSS <- log_rss(object = m1, x1 = x1, x2 = x2, ci = "se")

# Plot
ggplot(logRSS$df, aes(x = elevation_x1, y = log_rss)) +
  geom_hline(yintercept = 0, linetype = "dashed", color = "gray") +
  geom_ribbon(aes(ymin = lwr, ymax = upr), fill = "gray80") +
  geom_line() +
  xlab(expression("Elevation " * (x[1]))) +
  ylab("log-RSS") +
  ggtitle(expression("log-RSS" * (x[1] * ", " * x[2]))) +
  theme_bw()

# SSF -------------------------------------------------------
# Fit an SSF, then calculate log-RSS to visualize results.

 #Prepare data for SSF
ssf_data <- deer %>%
  steps_by_burst() %>%
  random_steps(n = 15) %>%
  extract_covariates(sh_forest) %>%
  mutate(forest = factor(sh.forest, levels = 1:2,
                    labels = c("forest", "non-forest")),
  cos_ta = cos(ta_),
  log_sl = log(sl_))

# Fit an SSF (note model = TRUE necessary for predict() to work)
m2 <- ssf_data %>%
  fit_clogit(case_ ~ forest + strata(step_id_), model = TRUE)

# Calculate log-RSS
# data.frame of x1s
x1 <- data.frame(forest = factor(c("forest", "non-forest")))
# data.frame of x2
x2 <- data.frame(forest = factor("forest", levels = levels(ssf_data$forest)))
# Calculate
logRSS <- log_rss(object = m2, x1 = x1, x2 = x2, ci = "se")

# Plot
ggplot(logRSS$df, aes(x = forest_x1, y = log_rss)) +
  geom_hline(yintercept = 0, linetype = "dashed", color = "gray") +
  geom_errorbar(aes(ymin = lwr, ymax = upr), width = 0.25) +
  geom_point(size = 3) +
  xlab(expression("Forest Cover " * (x[1]))) +
  ylab("log-RSS") +
  ggtitle(expression("log-RSS" * (x[1] * ", " * x[2]))) +
  theme_bw()

## End(Not run)

amt

Animal Movement Tools

v0.1.4
GPL-3
Authors
Johannes Signer [aut, cre], Brian Smith [ctb], Bjoern Reineking [ctb], Ulrike Schlaegel [ctb], John Fieberg [ctb], Scott LaPoint [dtc]
Initial release

We don't support your browser anymore

Please choose more modern alternatives, such as Google Chrome or Mozilla Firefox.