In this vignette, we illustrate how UHC plots can be used to validate step-selection functions fit to animal telemetry data. We will use the moose data set from Fieberg et al. (in review), which is also available as part of the uhcplots package. We can access the data and explore the first 6 records as follows:
library(uhcplots) data("moose12687") head(moose12687)
The data set consists of 689 GPS locations from 2013 and from 2014 from a single collared moose. The data were subsampled to achieve a consistent 4.25-hour fix rate. A total of 10 available locations were generated for each used location by randomly selecting 10 step lengths and 10 turn angles to project the animal forward in time from the previous location (see Street, G., J. Fieberg, A.R. Rodgers, M. Carstensen, R. Moen, S.A. Moore, S.K. Windels, and J.D. Forester. In press. Habitat functional response mitigates reduced foraging opportunity: implications for animal fitness and space use. Landscape Ecology for a more detailed description). The group of (1 used, 10 available) locaitons at each time step forms a stratum.
We defined resource availability at used and available locations as the proportional cover of 4 land cover types within a 50 m radius buffer: deciduous forest (decid50), mixedwood forest (mixed50), coniferous forest (conif50) and treed wetlands (treedwet50). Lastly, the variable step records the distance between each location from the previous location (after diving through by 1000).
The steps for generating a UHC plot are identical to those from our first example, except that we:
We will use the uhcsimstrat function rather than the uhcsim function to generate our predicted distributions of habitat covariates at the used locations in the test data set.
# load the survival library for the clogit function library(survival) library(dplyr)
The data set includes observations from 2 years (2013, 2014). We will treat the 2013 data as training data & 2014 as test data.
table(moose12687$year) mdat.train<-subset(moose12687, year==2013) mdat.test<-subset(moose12687, year==2014)
Here we explore UHC plots associated with the model containing deciduous forest (decid50), mixedwood forest (mixed50), coniferous forest (conif50) and treed wetlands (treedwet50). We also include step (for a justification,k see Forester et al. (2009) and Avgar et al. (2016)).
# set up textplot1 <- (expression(y %~% decid50 + mixed50 + conif50 + treedwet50 + step + strata(stratum))) form1a <- (presence ~ decid50 + mixed50 + conif50 + treedwet50 + step + strata(stratum)) form2a <- ~decid50 + mixed50 + conif50 + treedwet50 + step -1
We begin by fitting a step-selection function to the data for this animal for the year 2013 (training data).
ssf.train.full <- clogit(form1a, data=mdat.train) summary(ssf.train.full)
We then:
Create simulation envelopes for the environmental characteristics at the observed locations in the test data using the uhcsimstrat function
Get density estimates for the environmental characteristics associated with the observed locations in the test data and also associated with the randomly chosen locations generated by the uhcsimstrat function. To do this, we will again use the uhcdenscalc function
design.mat.test.full <- model.matrix(form2a, data=mdat.test) z <- model.matrix(~decid50 + mixed50 + conif50 + treedwet50 + step -1, data = mdat.test)[,-5] xchoice.full <- uhcsimstrat(nsims = 1000, xmat = design.mat.test.full, stratum = mdat.test$stratum, fit_ssf = ssf.train.full, z = z) denshats.decid.full <- uhcdenscalc(rand_sims = xchoice.full[,,1], dat = z[mdat.test$presence==1,1], avail = z[mdat.test$presence==0,1]) denshats.mixed.full <- uhcdenscalc(rand_sims=xchoice.full[,,2], dat=z[mdat.test$presence==1,2], avail=z[mdat.test$presence==0,2]) denshats.conif.full <- uhcdenscalc(rand_sims=xchoice.full[,,3], dat=z[mdat.test$presence==1,3], avail=z[mdat.test$presence==0,3]) denshats.treedwt.full <- uhcdenscalc(rand_sims=xchoice.full[,,4], dat=z[mdat.test$presence==1,4], avail=z[mdat.test$presence==0,4])
Then, we create UHC plots again using the uhcdensplot function illustrating:
the distribution of covariates at the available points in the test data (red dashed lines)
par(mfrow=c(1,4), mar=c(4,2,2,2), oma=c(3, 4, 7, 0), bty="L") # First, the density plots uhcdensplot(densdat = denshats.decid.full$densdat, densrand = denshats.decid.full$densrand, includeAvail = TRUE, densavail = denshats.decid.full$densavail) mtext(outer=F, side=2, line=3, "Density") mtext(outer=F, side=3, line=1, "Deciduous", cex=1) uhcdensplot(densdat = denshats.mixed.full$densdat, densrand = denshats.mixed.full$densrand, includeAvail = TRUE, densavail = denshats.mixed.full$densavail) mtext(outer=F, side=3, line=1, "Mixedwood", cex=1) uhcdensplot(densdat = denshats.conif.full$densdat, densrand = denshats.conif.full$densrand, includeAvail = TRUE, densavail = denshats.conif.full$densavail) mtext(outer=F, side=3, line=1, "Conifer", cex=1) uhcdensplot(densdat = denshats.treedwt.full$densdat, densrand = denshats.treedwt.full$densrand, includeAvail = TRUE, densavail = denshats.treedwt.full$densavail) mtext(outer=F, side=3, line=1, "Treed Wetlands", cex=1) mtext(outer=T, side=3, line=3, textplot1, cex=1)
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.