Title: | Animal Movement Modelling using Hidden Markov Models |
---|---|
Description: | Provides tools for animal movement modelling using hidden Markov models. These include processing of tracking data, fitting hidden Markov models to movement data, visualization of data and fitted model, decoding of the state process, etc. <doi:10.1111/2041-210X.12578>. |
Authors: | Theo Michelot (aut, cre), Roland Langrock (aut, cre), Toby Patterson (aut, cre), Brett McClintock (ctb), Eric Rexstad (ctb) |
Maintainer: | Theo Michelot <[email protected]> |
License: | GPL-3 |
Version: | 1.9 |
Built: | 2024-11-05 04:43:58 UTC |
Source: | https://github.com/TheoMichelot/moveHMM |
Akaike information criterion of a moveHMM model.
## S3 method for class 'moveHMM' AIC(object, ..., k = 2)
## S3 method for class 'moveHMM' AIC(object, ..., k = 2)
object |
A |
... |
Optional additional |
k |
Penalty per parameter. Default: 2; for classical AIC. |
The AIC of the model(s) provided. If several models are provided, the AICs are output in ascending order.
# m is a moveHMM object (as returned by fitHMM), automatically loaded with the package m <- example$m AIC(m)
# m is a moveHMM object (as returned by fitHMM), automatically loaded with the package m <- example$m AIC(m)
Simulation-based computation of confidence intervals for the parameters of the angle distribution.
Used in CI
.
angleCI(m, alpha, nbSims = 10^6)
angleCI(m, alpha, nbSims = 10^6)
m |
A |
alpha |
Range of the confidence intervals. Default: 0.95 (i.e. 95% CIs). |
nbSims |
Number of simulations. Default: 10^6. |
A list of the following objects:
lower |
Lower bound of the confidence interval for the parameters of the angle distribution |
upper |
Upper bound of the confidence interval for the parameters of the angle distribution |
Computes the confidence intervals of the step length and turning angle parameters, as well as for the transition probabilities regression parameters.
CI(m, alpha = 0.95, nbSims = 10^6)
CI(m, alpha = 0.95, nbSims = 10^6)
m |
A |
alpha |
Range of the confidence intervals. Default: 0.95 (i.e. 95% CIs). |
nbSims |
Number of simulations in the computation of the CIs for the angle parameters. Default: 10^6. |
A list of the following objects:
stepPar |
Confidence intervals for the parameters of the step lengths distribution |
anglePar |
Confidence intervals for the parameters of the turning angles distribution |
beta |
Confidence intervals for the regression coefficients of the transition probabilities. |
# m is a moveHMM object (as returned by fitHMM), automatically loaded with the package m <- example$m CI(m)
# m is a moveHMM object (as returned by fitHMM), automatically loaded with the package m <- example$m CI(m)
Probability density function of the exponential distribution (written in C++)
dexp_rcpp(x, rate, foo = 0)
dexp_rcpp(x, rate, foo = 0)
x |
Vector of quantiles |
rate |
Rate |
foo |
Unused (for compatibility with template) |
Vector of densities
Probability density function of the gamma distribution (written in C++)
dgamma_rcpp(x, mu, sigma)
dgamma_rcpp(x, mu, sigma)
x |
Vector of quantiles |
mu |
Mean |
sigma |
Standard deviation |
Vector of densities
Probability density function of the log-normal distribution (written in C++)
dlnorm_rcpp(x, meanlog, sdlog)
dlnorm_rcpp(x, meanlog, sdlog)
x |
Vector of quantiles |
meanlog |
Mean of the distribution on the log-scale |
sdlog |
Standard deviation of the distribution on the log-scale |
Vector of densities
Probability density function of the Von Mises distribution, defined as a function of the modified Bessel function of order 0 (written in C++)
dvm_rcpp(x, mu, kappa)
dvm_rcpp(x, mu, kappa)
x |
Vector of quantiles |
mu |
Mean |
kappa |
Concentration |
Vector of densities
Probability density function of the Weibull distribution (written in C++)
dweibull_rcpp(x, shape, scale)
dweibull_rcpp(x, shape, scale)
x |
Vector of quantiles |
shape |
Shape |
scale |
Scale |
Vector of densities
Probability density function of the wrapped Cauchy distribution (written in C++)
dwrpcauchy_rcpp(x, mu, rho)
dwrpcauchy_rcpp(x, mu, rho)
x |
Vector of quantiles |
mu |
Mean |
rho |
Concentration |
Vector of densities
It is a data frame with the following columns:
ID
Track identifier
Easting
Easting coordinate of locations
Northing
Northing coordinate of locations
dist_water
Distance of elk to water (in metres)
elk_data
elk_data
This data is generated by the function exGen
, and used in the examples and tests of
other functions to keep them as short as possible.
example
example
It is a list of the following objects:
data
A moveData
object
m
A moveHMM
object
simPar
The parameters used to simulate data
par0
The initial parameters in the optimization to fit m
Generate the file data/example.RData
, used in other functions' examples and unit tests.
exGen()
exGen()
Fit an hidden Markov model to the data provided, using numerical optimization of the log-likelihood function.
fitHMM( data, nbStates, stepPar0, anglePar0 = NULL, beta0 = NULL, delta0 = NULL, formula = ~1, stepDist = c("gamma", "weibull", "lnorm", "exp"), angleDist = c("vm", "wrpcauchy", "none"), angleMean = NULL, stationary = FALSE, knownStates = NULL, verbose = 0, nlmPar = NULL, fit = TRUE )
fitHMM( data, nbStates, stepPar0, anglePar0 = NULL, beta0 = NULL, delta0 = NULL, formula = ~1, stepDist = c("gamma", "weibull", "lnorm", "exp"), angleDist = c("vm", "wrpcauchy", "none"), angleMean = NULL, stationary = FALSE, knownStates = NULL, verbose = 0, nlmPar = NULL, fit = TRUE )
data |
An object |
nbStates |
Number of states of the HMM. |
stepPar0 |
Vector of initial state-dependent step length distribution parameters.
The parameters should be in the order expected by the pdf of |
anglePar0 |
Vector of initial state-dependent turning angle distribution parameters.
The parameters should be in the order expected by the pdf of |
beta0 |
Initial matrix of regression coefficients for the transition probabilities (more
information in "Details").
Default: |
delta0 |
Initial value for the initial distribution of the HMM. Default: |
formula |
Regression formula for the covariates. Default: |
stepDist |
Name of the distribution of the step lengths (as a character string). Supported distributions are: gamma, weibull, lnorm, exp. Default: gamma. |
angleDist |
Name of the distribution of the turning angles (as a character string).
Supported distributions are: vm, wrpcauchy. Set to |
angleMean |
Vector of means of turning angles if not estimated (one for each state).
Default: |
stationary |
|
knownStates |
Vector of values of the state process which are known prior to fitting the model (if any). Default: NULL (states are not known). This should be a vector with length the number of rows of 'data'; each element should either be an integer (the value of the known states) or NA if the state is not known. |
verbose |
Determines the print level of the optimizer. The default value of 0 means that no printing occurs, a value of 1 means that the first and last iterations of the optimization are detailed, and a value of 2 means that each iteration of the optimization is detailed. |
nlmPar |
List of parameters to pass to the optimization function |
fit |
|
The matrix beta
of regression coefficients for the transition probabilities has
one row for the intercept, plus one row for each covariate, and one column for
each non-diagonal element of the transition probability matrix. For example, in a 3-state
HMM with 2 covariates, the matrix beta
has three rows (intercept + two covariates)
and six columns (six non-diagonal elements in the 3x3 transition probability matrix -
filled in row-wise).
In a covariate-free model (default), beta
has one row, for the intercept.
The choice of initial parameters is crucial to fit a model. The algorithm might not find the global optimum of the likelihood function if the initial parameters are poorly chosen.
A moveHMM
object, i.e. a list of:
mle |
The maximum likelihood estimates of the parameters of the model
(if the numerical algorithm has indeed identified the global maximum of the
likelihood function), which is a list of: |
data |
The movement data |
mod |
The object returned by the numerical optimizer |
conditions |
A few conditions used to fit the model ( |
rawCovs |
Raw covariate values, as found in the data (if any). Used in
|
knownStates |
Vector of states known a priori, as provided in input (if
any, |
nlmTime |
Computing time for optimisation, obtained with system.time |
Patterson T.A., Basson M., Bravington M.V., Gunn J.S. 2009. Classifying movement behaviour in relation to environmental conditions using hidden Markov models. Journal of Animal Ecology, 78 (6), 1113-1123.
Langrock R., King R., Matthiopoulos J., Thomas L., Fortin D., Morales J.M. 2012. Flexible and practical modeling of animal telemetry data: hidden Markov models and extensions. Ecology, 93 (11), 2336-2342.
### 1. simulate data # define all the arguments of simData nbAnimals <- 2 nbStates <- 2 nbCovs <- 2 mu<-c(15,50) sigma<-c(10,20) angleMean <- c(pi,0) kappa <- c(0.7,1.5) stepPar <- c(mu,sigma) anglePar <- c(angleMean,kappa) stepDist <- "gamma" angleDist <- "vm" zeroInflation <- FALSE obsPerAnimal <- c(50,100) data <- simData(nbAnimals=nbAnimals,nbStates=nbStates,stepDist=stepDist,angleDist=angleDist, stepPar=stepPar,anglePar=anglePar,nbCovs=nbCovs,zeroInflation=zeroInflation, obsPerAnimal=obsPerAnimal) ### 2. fit the model to the simulated data # define initial values for the parameters mu0 <- c(20,70) sigma0 <- c(10,30) kappa0 <- c(1,1) stepPar0 <- c(mu0,sigma0) # no zero-inflation, so no zero-mass included anglePar0 <- kappa0 # the angle mean is not estimated, so only the concentration parameter is needed formula <- ~cov1+cos(cov2) m <- fitHMM(data=data,nbStates=nbStates,stepPar0=stepPar0,anglePar0=anglePar0,formula=formula, stepDist=stepDist,angleDist=angleDist,angleMean=angleMean) print(m)
### 1. simulate data # define all the arguments of simData nbAnimals <- 2 nbStates <- 2 nbCovs <- 2 mu<-c(15,50) sigma<-c(10,20) angleMean <- c(pi,0) kappa <- c(0.7,1.5) stepPar <- c(mu,sigma) anglePar <- c(angleMean,kappa) stepDist <- "gamma" angleDist <- "vm" zeroInflation <- FALSE obsPerAnimal <- c(50,100) data <- simData(nbAnimals=nbAnimals,nbStates=nbStates,stepDist=stepDist,angleDist=angleDist, stepPar=stepPar,anglePar=anglePar,nbCovs=nbCovs,zeroInflation=zeroInflation, obsPerAnimal=obsPerAnimal) ### 2. fit the model to the simulated data # define initial values for the parameters mu0 <- c(20,70) sigma0 <- c(10,30) kappa0 <- c(1,1) stepPar0 <- c(mu0,sigma0) # no zero-inflation, so no zero-mass included anglePar0 <- kappa0 # the angle mean is not estimated, so only the concentration parameter is needed formula <- ~cov1+cos(cov2) m <- fitHMM(data=data,nbStates=nbStates,stepPar0=stepPar0,anglePar0=anglePar0,formula=formula, stepDist=stepDist,angleDist=angleDist,angleMean=angleMean) print(m)
Discrete colour palette for states
getPalette(nbStates)
getPalette(nbStates)
nbStates |
Number of states |
Vector of colours, of length nbStates.
Data to produce plots of fitted model
getPlotData(m, type, format = "wide", alpha = 0.95)
getPlotData(m, type, format = "wide", alpha = 0.95)
m |
Fitted HMM object, as output by fitHMM. |
type |
Type of plot, one of: "dist", "tpm", "stat" |
format |
Format of data, either "wide" (for base graphics) or "long" (for ggplot) |
alpha |
Level of confidence intervals. Default: 0.95, i.e., 95% confidence intervals |
If type = "dist", the function evaluates each state-dependent distribution over the range of observed variable (step length or turning angle), and weighs them by the proportion of time spent in each state (obtained from Viterbi state sequence).
If type = "tpm", the function returns transition probabilities estimated over a range of covariate values. Other covariates are fixed to their mean values.
Data frame (or list of data frames) containing data in a format that can easily be plotted. If type = "dist", the output is a list with two elements, "step" and "angle". If type = "tpm" or "stat", the output is a list with one element for each covariate. See details for more extensive description of output.
Data frame of the first three tracks from Michelot et al. (2016), with columns:
ID
Track identifier
x
Easting coordinate of locations
y
Northing coordinate of locations
slope
Terrain slope (in degrees)
temp
Air temperature (in degrees Celsius)
haggis_data
haggis_data
Check that an object is of class moveData
. Used in fitHMM
.
is.moveData(x)
is.moveData(x)
x |
An R object |
TRUE
if x
is of class moveData
, FALSE
otherwise.
Check that an object is of class moveHMM
. Used in CI
,
plotPR
, plotStates
, pseudoRes
, stateProbs
,
and viterbi
.
is.moveHMM(x)
is.moveHMM(x)
x |
An R object |
TRUE
if x
is of class moveHMM
, FALSE
otherwise.
Used in stateProbs
and pseudoRes
.
logAlpha(m)
logAlpha(m)
m |
A |
The matrix of forward log-probabilities.
## Not run: # m is a moveHMM object (as returned by fitHMM), automatically loaded with the package m <- example$m la <- logAlpha(m) ## End(Not run)
## Not run: # m is a moveHMM object (as returned by fitHMM), automatically loaded with the package m <- example$m la <- logAlpha(m) ## End(Not run)
Used in stateProbs
.
logBeta(m)
logBeta(m)
m |
A |
The matrix of backward log-probabilities.
## Not run: # m is a moveHMM object (as returned by fitHMM), automatically loaded with the package m <- example$m lb <- logBeta(m) ## End(Not run)
## Not run: # m is a moveHMM object (as returned by fitHMM), automatically loaded with the package m <- example$m lb <- logBeta(m) ## End(Not run)
moveData
objectsThis constructor is not intended to be used, except inside the function
prepData
. Refer to the documentation for that function.
moveData(data)
moveData(data)
data |
A data frame with columns: |
An object of class moveData
.
moveHMM
objectsThis constructor is not intended to be used, except inside the function
fitHMM
. Refer to the documentation for that function.
moveHMM(m)
moveHMM(m)
m |
A list of attributes of the fitted model: |
An object moveHMM
.
Scales each parameter from its natural interval to the set of real numbers, to allow for unconstrained optimization. Used during the optimization of the log-likelihood.
n2w(par, bounds, beta, delta = NULL, nbStates, estAngleMean)
n2w(par, bounds, beta, delta = NULL, nbStates, estAngleMean)
par |
Vector of state-dependent distributions parameters. |
bounds |
Matrix with 2 columns and as many rows as there are elements in |
beta |
Matrix of regression coefficients for the transition probabilities. |
delta |
Initial distribution. Default: |
nbStates |
The number of states of the HMM. |
estAngleMean |
|
A vector of unconstrained parameters.
## Not run: nbStates <- 3 par <- c(0.001,0.999,0.5,0.001,1500.3,7.1) bounds <- matrix(c(0,1, # bounds for first parameter 0,1, # bounds for second parameter 0,1, # ... 0,Inf, 0,Inf, 0,Inf), byrow=TRUE,ncol=2) beta <- matrix(rnorm(18),ncol=6,nrow=3) delta <- c(0.6,0.3,0.1) # vector of working parameters wpar <- n2w(par=par,bounds=bounds,beta=beta,delta=delta,nbStates=nbStates, estAngleMean=FALSE) ## End(Not run)
## Not run: nbStates <- 3 par <- c(0.001,0.999,0.5,0.001,1500.3,7.1) bounds <- matrix(c(0,1, # bounds for first parameter 0,1, # bounds for second parameter 0,1, # ... 0,Inf, 0,Inf, 0,Inf), byrow=TRUE,ncol=2) beta <- matrix(rnorm(18),ncol=6,nrow=3) delta <- c(0.6,0.3,0.1) # vector of working parameters wpar <- n2w(par=par,bounds=bounds,beta=beta,delta=delta,nbStates=nbStates, estAngleMean=FALSE) ## End(Not run)
Negative log-likelihood function
nLogLike( wpar, nbStates, bounds, parSize, data, stepDist = c("gamma", "weibull", "lnorm", "exp"), angleDist = c("vm", "wrpcauchy", "none"), angleMean = NULL, zeroInflation = FALSE, stationary = FALSE, knownStates = NULL )
nLogLike( wpar, nbStates, bounds, parSize, data, stepDist = c("gamma", "weibull", "lnorm", "exp"), angleDist = c("vm", "wrpcauchy", "none"), angleMean = NULL, zeroInflation = FALSE, stationary = FALSE, knownStates = NULL )
wpar |
Vector of working parameters. |
nbStates |
Number of states of the HMM. |
bounds |
Matrix with 2 columns and as many rows as there are elements in |
parSize |
Vector of two values: number of parameters of the step length distribution, number of parameters of the turning angle distribution. |
data |
An object |
stepDist |
Name of the distribution of the step lengths (as a character string). Supported distributions are: gamma, weibull, lnorm, exp. Default: gamma. |
angleDist |
Name of the distribution of the turning angles (as a character string).
Supported distributions are: vm, wrpcauchy. Set to |
angleMean |
Vector of means of turning angles if not estimated (one for each state).
Default: |
zeroInflation |
|
stationary |
|
knownStates |
Vector of values of the state process which are known prior to fitting the model (if any). Default: NULL (states are not known). This should be a vector with length the number of rows of 'data'; each element should either be an integer (the value of the known states) or NA if the state is not known. |
The negative log-likelihood of the parameters given the data.
## Not run: # data is a moveData object (as returned by prepData), automatically loaded with the package data <- example$data simPar <- example$simPar par0 <- example$par0 estAngleMean <- is.null(simPar$angleMean) bounds <- parDef(simPar$stepDist,simPar$angleDist,simPar$nbStates, estAngleMean,simPar$zeroInflation)$bounds parSize <- parDef(simPar$stepDist,simPar$angleDist,simPar$nbStates, estAngleMean,simPar$zeroInflation)$parSize par <- c(par0$stepPar0,par0$anglePar0) wpar <- n2w(par,bounds,par0$beta0,par0$delta0,simPar$nbStates,FALSE) l <- nLogLike(wpar=wpar,nbStates=simPar$nbStates,bounds=bounds,parSize=parSize,data=data, stepDist=simPar$stepDist,angleDist=simPar$angleDist,angleMean=simPar$angleMean, zeroInflation=simPar$zeroInflation) ## End(Not run)
## Not run: # data is a moveData object (as returned by prepData), automatically loaded with the package data <- example$data simPar <- example$simPar par0 <- example$par0 estAngleMean <- is.null(simPar$angleMean) bounds <- parDef(simPar$stepDist,simPar$angleDist,simPar$nbStates, estAngleMean,simPar$zeroInflation)$bounds parSize <- parDef(simPar$stepDist,simPar$angleDist,simPar$nbStates, estAngleMean,simPar$zeroInflation)$parSize par <- c(par0$stepPar0,par0$anglePar0) wpar <- n2w(par,bounds,par0$beta0,par0$delta0,simPar$nbStates,FALSE) l <- nLogLike(wpar=wpar,nbStates=simPar$nbStates,bounds=bounds,parSize=parSize,data=data, stepDist=simPar$stepDist,angleDist=simPar$angleDist,angleMean=simPar$angleMean, zeroInflation=simPar$zeroInflation) ## End(Not run)
Computation of the negative log-likelihood (forward algorithm - written in C++)
nLogLike_rcpp( nbStates, beta, covs, data, stepDist, angleDist, stepPar, anglePar, delta, aInd, zeroInflation, stationary, knownStates )
nLogLike_rcpp( nbStates, beta, covs, data, stepDist, angleDist, stepPar, anglePar, delta, aInd, zeroInflation, stationary, knownStates )
nbStates |
Number of states |
beta |
Matrix of regression coefficients for the transition probabilities |
covs |
Covariates |
data |
A |
stepDist |
The name of the step length distribution |
angleDist |
The name of the turning angle distribution |
stepPar |
State-dependent parameters of the step length distribution |
anglePar |
State-dependent parameters of the turning angle distribution |
delta |
Stationary distribution |
aInd |
Vector of indices of the rows at which the data switches to another animal |
zeroInflation |
|
stationary |
|
knownStates |
Vector of values of the state process which are known prior to fitting the model (if any). Default: NULL (states are not known). This should be a vector with length the number of rows of 'data'; each element should either be an integer (the value of the known states) or NA if the state is not known. |
Negative log-likelihood
Parameters definition
parDef(stepDist, angleDist, nbStates, estAngleMean, zeroInflation)
parDef(stepDist, angleDist, nbStates, estAngleMean, zeroInflation)
stepDist |
Name of the distribution of the step lengths. |
angleDist |
Name of the distribution of the turning angles.
Set to |
nbStates |
Number of states of the HMM. |
estAngleMean |
|
zeroInflation |
|
A list of:
parSize |
Vector of two values: number of parameters of the step length distribution, number of parameters of the turning angle distribution |
bounds |
Matrix with 2 columns and |
parNames |
Names of parameters of step distribution (the names of the parameters of the angle distribution are always the same). |
moveData
Plot moveData
## S3 method for class 'moveData' plot(x, animals = NULL, compact = FALSE, ask = TRUE, breaks = "Sturges", ...)
## S3 method for class 'moveData' plot(x, animals = NULL, compact = FALSE, ask = TRUE, breaks = "Sturges", ...)
x |
An object |
animals |
Vector of indices or IDs of animals for which information will be plotted.
Default: |
compact |
|
ask |
If |
breaks |
Histogram parameter. See |
... |
Currently unused. For compatibility with generic method. |
# data is a moveData object (as returned by prepData), automatically loaded with the package data <- example$data plot(data,compact=TRUE,breaks=20,ask=FALSE)
# data is a moveData object (as returned by prepData), automatically loaded with the package data <- example$data plot(data,compact=TRUE,breaks=20,ask=FALSE)
moveHMM
Plot the fitted step and angle densities over histograms of the data, transition probabilities as functions of the covariates, and maps of the animals' tracks colored by the decoded states.
## S3 method for class 'moveHMM' plot( x, animals = NULL, ask = TRUE, breaks = "Sturges", col = NULL, plotTracks = TRUE, plotCI = FALSE, alpha = 0.95, ... )
## S3 method for class 'moveHMM' plot( x, animals = NULL, ask = TRUE, breaks = "Sturges", col = NULL, plotTracks = TRUE, plotCI = FALSE, alpha = 0.95, ... )
x |
Object |
animals |
Vector of indices or IDs of animals for which information will be plotted.
Default: |
ask |
If |
breaks |
Histogram parameter. See |
col |
Vector or colors for the states (one color per state). |
plotTracks |
If |
plotCI |
If |
alpha |
Significance level of the confidence intervals if plotCI=TRUE. Default: 0.95 (i.e. 95% CIs). |
... |
Currently unused. For compatibility with generic method. |
The state-dependent densities are weighted by the frequency of each state in the most
probable state sequence (decoded with the function viterbi
). For example, if the
most probable state sequence indicates that one third of observations correspond to the first
state, and two thirds to the second state, the plots of the densities in the first state are
weighted by a factor 1/3, and in the second state by a factor 2/3.
# m is a moveHMM object (as returned by fitHMM), automatically loaded with the package m <- example$m plot(m,ask=TRUE,animals=1,breaks=20)
# m is a moveHMM object (as returned by fitHMM), automatically loaded with the package m <- example$m plot(m,ask=TRUE,animals=1,breaks=20)
Plots time series, qq-plots (against the standard normal distribution), and sample ACF functions of the pseudo-residuals
plotPR(m)
plotPR(m)
m |
A |
If some turning angles in the data are equal to pi, the corresponding pseudo-residuals will not be included. Indeed, given that the turning angles are defined on (-pi,pi], an angle of pi results in a pseudo-residual of +Inf (check Section 6.2 of reference for more information on the computation of pseudo-residuals).
If some steps are of length zero (i.e. if there is zero-inflation), the corresponding pseudo- residuals are shown as segments, because pseudo-residuals for discrete data are defined as segments (see Zucchini and MacDonald, 2009, Section 6.2).
Zucchini, W. and MacDonald, I.L. 2009. Hidden Markov Models for Time Series: An Introduction Using R. Chapman & Hall (London).
# m is a moveHMM object (as returned by fitHMM), automatically loaded with the package m <- example$m plotPR(m)
# m is a moveHMM object (as returned by fitHMM), automatically loaded with the package m <- example$m plotPR(m)
Plot tracking data on a satellite map. This function only works with longitude
and latitude values (not with UTM coordinates), and uses the package ggmap
to fetch a satellite image from Google. An Internet connection is required to use
this function.
plotSat( data, zoom = NULL, location = NULL, segments = TRUE, compact = TRUE, col = NULL, alpha = 1, size = 1, states = NULL, animals = NULL, ask = TRUE, return = FALSE )
plotSat( data, zoom = NULL, location = NULL, segments = TRUE, compact = TRUE, col = NULL, alpha = 1, size = 1, states = NULL, animals = NULL, ask = TRUE, return = FALSE )
data |
Data frame of the data, with necessary fields 'x' (longitude values) and 'y' (latitude values). |
zoom |
The zoom level, as defined for |
location |
Location of the center of the map to be plotted. |
segments |
|
compact |
|
col |
Palette of colours to use for the dots and segments. If not specified, uses default palette. |
alpha |
Transparency argument for |
size |
Size argument for |
states |
A sequence of integers, corresponding to the decoded states for these data (such that the observations are colored by states). |
animals |
Vector of indices or IDs of animals/tracks to be plotted.
Default: |
ask |
If |
return |
If |
If the plot displays the message "Sorry, we have no imagery here", try a lower level of zoom.
D. Kahle and H. Wickham. ggmap: Spatial Visualization with ggplot2. The R Journal, 5(1), 144-161. URL: http://journal.r-project.org/archive/2013-1/kahle-wickham.pdf
Plot the states and states probabilities.
plotStates(m, animals = NULL, ask = TRUE)
plotStates(m, animals = NULL, ask = TRUE)
m |
A |
animals |
Vector of indices or IDs of animals for which states will be plotted. |
ask |
If |
# m is a moveHMM object (as returned by fitHMM), automatically loaded with the package m <- example$m # plot states for first and second animals plotStates(m,animals=c(1,2))
# m is a moveHMM object (as returned by fitHMM), automatically loaded with the package m <- example$m # plot states for first and second animals plotStates(m,animals=c(1,2))
Plot stationary state probabilities
plotStationary(m, col = NULL, plotCI = FALSE, alpha = 0.95)
plotStationary(m, col = NULL, plotCI = FALSE, alpha = 0.95)
m |
An object |
col |
Vector or colors for the states (one color per state). |
plotCI |
Logical. Should 95% confidence intervals be plotted? (Default: FALSE) |
alpha |
Significance level of the confidence intervals if plotCI=TRUE. Default: 0.95 (i.e. 95% CIs). |
# m is a moveHMM object (as returned by fitHMM), automatically loaded with the package m <- example$m plotStationary(m)
# m is a moveHMM object (as returned by fitHMM), automatically loaded with the package m <- example$m plotStationary(m)
Predict stationary state probabilities
predictStationary( m, newData, beta = m$mle$beta, returnCI = FALSE, alpha = 0.95 )
predictStationary( m, newData, beta = m$mle$beta, returnCI = FALSE, alpha = 0.95 )
m |
Fitted moveHMM object, as returned by |
newData |
Data frame with columns for the covariates |
beta |
Optional matrix of regression coefficients for the transition
probability model. By default, uses estimates in |
returnCI |
Logical indicating whether confidence intervals should be returned. Default: FALSE. |
alpha |
Confidence level if returnCI = TRUE. Default: 0.95, i.e., 95% confidence intervals. |
List with elements 'mle', 'lci', and 'uci' (the last two only if returnCI = TRUE). Each element is a matrix of stationary state probabilities with one row for each row of newData and one column for each state.
Predict transition probabilities for new covariate values
predictTPM(m, newData, beta = m$mle$beta, returnCI = FALSE, alpha = 0.95)
predictTPM(m, newData, beta = m$mle$beta, returnCI = FALSE, alpha = 0.95)
m |
Fitted moveHMM object, as returned by |
newData |
Data frame with columns for the covariates |
beta |
Optional matrix of regression coefficients for the transition
probability model. By default, uses estimates in |
returnCI |
Logical indicating whether confidence intervals should be returned. Default: FALSE. |
alpha |
Confidence level if returnCI = TRUE. Default: 0.95, i.e., 95% confidence intervals. |
List with elements 'mle', 'lci', and 'uci' (the last two only if returnCI = TRUE). Each element is an array, where each layer is a transition probability matrix corresponding to a row of newData.
Preprocessing of the tracking data
prepData( trackData, type = c("LL", "UTM"), coordNames = c("x", "y"), LLangle = NULL )
prepData( trackData, type = c("LL", "UTM"), coordNames = c("x", "y"), LLangle = NULL )
trackData |
A dataframe of the tracking data, including at least coordinates
(either longitude/latitude values or cartesian coordinates), and optionnaly a field |
type |
|
coordNames |
Names of the columns of coordinates in the data frame. Default: |
LLangle |
Logical. If TRUE, the turning angle is calculated with |
An object moveData
, i.e. a dataframe of:
ID |
The ID(s) of the observed animal(s) |
step |
The step lengths - in kilometers if longitude/latitude provided, and in the metrics of the data otherwise |
angle |
The turning angles (if any) - in radians |
x |
Either Easting or longitude (or e.g. depth for 1D data) |
y |
Either Northing or latitude (all zero if 1D data) |
... |
Covariates (if any) |
coord1 <- c(1,2,3,4,5,6,7,8,9,10) coord2 <- c(1,1,1,2,2,2,1,1,1,2) trackData <- data.frame(coord1=coord1,coord2=coord2) d <- prepData(trackData,type='UTM',coordNames=c("coord1","coord2"))
coord1 <- c(1,2,3,4,5,6,7,8,9,10) coord2 <- c(1,1,1,2,2,2,1,1,1,2) trackData <- data.frame(coord1=coord1,coord2=coord2) d <- prepData(trackData,type='UTM',coordNames=c("coord1","coord2"))
moveHMM
Print moveHMM
## S3 method for class 'moveHMM' print(x, ...)
## S3 method for class 'moveHMM' print(x, ...)
x |
A |
... |
Currently unused. For compatibility with generic method. |
# m is a moveHMM object (as returned by fitHMM), automatically loaded with the package m <- example$m print(m)
# m is a moveHMM object (as returned by fitHMM), automatically loaded with the package m <- example$m print(m)
The pseudo-residuals of a moveHMM model, as described in Zucchini and McDonad (2009).
pseudoRes(m)
pseudoRes(m)
m |
A |
If some turning angles in the data are equal to pi, the corresponding pseudo-residuals will not be included. Indeed, given that the turning angles are defined on (-pi,pi], an angle of pi results in a pseudo-residual of +Inf (check Section 6.2 of reference for more information on the computation of pseudo-residuals).
A list of:
stepRes |
The pseudo-residuals for the step lengths |
angleRes |
The pseudo-residuals for the turning angles |
Zucchini, W. and MacDonald, I.L. 2009. Hidden Markov Models for Time Series: An Introduction Using R. Chapman & Hall (London).
# m is a moveHMM object (as returned by fitHMM), automatically loaded with the package m <- example$m res <- pseudoRes(m) qqnorm(res$stepRes) qqnorm(res$angleRes)
# m is a moveHMM object (as returned by fitHMM), automatically loaded with the package m <- example$m res <- pseudoRes(m) qqnorm(res$stepRes) qqnorm(res$angleRes)
Simulates movement data from an HMM.
simData( nbAnimals = 1, nbStates = 2, stepDist = c("gamma", "weibull", "lnorm", "exp"), angleDist = c("vm", "wrpcauchy", "none"), stepPar = NULL, anglePar = NULL, beta = NULL, covs = NULL, nbCovs = 0, zeroInflation = FALSE, obsPerAnimal = c(500, 1500), model = NULL, states = FALSE )
simData( nbAnimals = 1, nbStates = 2, stepDist = c("gamma", "weibull", "lnorm", "exp"), angleDist = c("vm", "wrpcauchy", "none"), stepPar = NULL, anglePar = NULL, beta = NULL, covs = NULL, nbCovs = 0, zeroInflation = FALSE, obsPerAnimal = c(500, 1500), model = NULL, states = FALSE )
nbAnimals |
Number of observed individuals to simulate. |
nbStates |
Number of behavioural states to simulate. |
stepDist |
Name of the distribution of the step lengths (as a character string). Supported distributions are: gamma, weibull, lnorm, exp. Default: gamma. |
angleDist |
Name of the distribution of the turning angles (as a character string).
Supported distributions are: vm, wrpcauchy. Set to |
stepPar |
Parameters of the step length distribution. |
anglePar |
Parameters of the turning angle distribution. |
beta |
Matrix of regression parameters for the transition probabilities (more information in "Details"). |
covs |
Covariate values to include in the model, as a dataframe. Default: |
nbCovs |
Number of covariates to simulate (0 by default). Does not need to be specified of
|
zeroInflation |
|
obsPerAnimal |
Either the number of the number of observations per animal (if single value),
or the bounds of the number of observations per animal (if vector of two values). In the latter case,
the numbers of obervations generated for each animal are uniformously picked from this interval.
Default: |
model |
A moveHMM object. This option can be used to simulate from a fitted model. Default: NULL. Note that, if this argument is specified, most other arguments will be ignored – except for nbAnimals, obsPerAnimal, covs (if covariate values different from those in the data should be specified), and states. |
states |
|
The matrix beta
of regression coefficients for the transition probabilities has
one row for the intercept, plus one row for each covariate, and one column for
each non-diagonal element of the transition probability matrix. For example, in a 3-state
HMM with 2 covariates, the matrix beta
has three rows (intercept + two covariates)
and six columns (six non-diagonal elements in the 3x3 transition probability matrix - filled in
row-wise).
In a covariate-free model (default), beta
has one row, for the intercept.
If the length of covariate values passed (either through 'covs', or 'model') is not the same as the number of observations suggested by 'nbAnimals' and 'obsPerAnimal', then the series of covariates is either shortened (removing last values - if too long) or extended (starting over from the first values - if too short).
An object moveData, i.e. a dataframe of:
ID |
The ID(s) of the observed animal(s) |
step |
The step lengths |
angle |
The turning angles (if any) |
x |
Either easting or longitude |
y |
Either northing or latitude |
... |
Covariates (if any) |
# 1. Pass a fitted model to simulate from # (m is a moveHMM object - as returned by fitHMM - automatically loaded with the package) # We keep the default nbAnimals=1. m <- example$m obsPerAnimal=c(50,100) data <- simData(model=m,obsPerAnimal=obsPerAnimal) # 2. Pass the parameters of the model to simulate from stepPar <- c(1,10,1,5,0.2,0.3) # mean1, mean2, sd1, sd2, z1, z2 anglePar <- c(pi,0,0.5,2) # mean1, mean2, k1, k2 stepDist <- "gamma" angleDist <- "vm" data <- simData(nbAnimals=5,nbStates=2,stepDist=stepDist,angleDist=angleDist,stepPar=stepPar, anglePar=anglePar,nbCovs=2,zeroInflation=TRUE,obsPerAnimal=obsPerAnimal) stepPar <- c(1,10,1,5) # mean1, mean2, sd1, sd2 anglePar <- c(pi,0,0.5,0.7) # mean1, mean2, k1, k2 stepDist <- "weibull" angleDist <- "wrpcauchy" data <- simData(nbAnimals=5,nbStates=2,stepDist=stepDist,angleDist=angleDist,stepPar=stepPar, anglePar=anglePar,obsPerAnimal=obsPerAnimal) # step length only and zero-inflation stepPar <- c(1,10,1,5,0.2,0.3) # mean1, mean2, sd1, sd2, z1, z2 stepDist <- "gamma" data <- simData(nbAnimals=5,nbStates=2,stepDist=stepDist,angleDist="none",stepPar=stepPar, nbCovs=2,zeroInflation=TRUE,obsPerAnimal=obsPerAnimal) # include covariates # (note that it is useless to specify "nbCovs", which is determined # by the number of columns of "cov") cov <- data.frame(temp=rnorm(500,20,5)) stepPar <- c(1,10,1,5) # mean1, mean2, sd1, sd2 anglePar <- c(pi,0,0.5,2) # mean1, mean2, k1, k2 stepDist <- "gamma" angleDist <- "vm" data <- simData(nbAnimals=5,nbStates=2,stepDist=stepDist,angleDist=angleDist,stepPar=stepPar, anglePar=anglePar,covs=cov)
# 1. Pass a fitted model to simulate from # (m is a moveHMM object - as returned by fitHMM - automatically loaded with the package) # We keep the default nbAnimals=1. m <- example$m obsPerAnimal=c(50,100) data <- simData(model=m,obsPerAnimal=obsPerAnimal) # 2. Pass the parameters of the model to simulate from stepPar <- c(1,10,1,5,0.2,0.3) # mean1, mean2, sd1, sd2, z1, z2 anglePar <- c(pi,0,0.5,2) # mean1, mean2, k1, k2 stepDist <- "gamma" angleDist <- "vm" data <- simData(nbAnimals=5,nbStates=2,stepDist=stepDist,angleDist=angleDist,stepPar=stepPar, anglePar=anglePar,nbCovs=2,zeroInflation=TRUE,obsPerAnimal=obsPerAnimal) stepPar <- c(1,10,1,5) # mean1, mean2, sd1, sd2 anglePar <- c(pi,0,0.5,0.7) # mean1, mean2, k1, k2 stepDist <- "weibull" angleDist <- "wrpcauchy" data <- simData(nbAnimals=5,nbStates=2,stepDist=stepDist,angleDist=angleDist,stepPar=stepPar, anglePar=anglePar,obsPerAnimal=obsPerAnimal) # step length only and zero-inflation stepPar <- c(1,10,1,5,0.2,0.3) # mean1, mean2, sd1, sd2, z1, z2 stepDist <- "gamma" data <- simData(nbAnimals=5,nbStates=2,stepDist=stepDist,angleDist="none",stepPar=stepPar, nbCovs=2,zeroInflation=TRUE,obsPerAnimal=obsPerAnimal) # include covariates # (note that it is useless to specify "nbCovs", which is determined # by the number of columns of "cov") cov <- data.frame(temp=rnorm(500,20,5)) stepPar <- c(1,10,1,5) # mean1, mean2, sd1, sd2 anglePar <- c(pi,0,0.5,2) # mean1, mean2, k1, k2 stepDist <- "gamma" angleDist <- "vm" data <- simData(nbAnimals=5,nbStates=2,stepDist=stepDist,angleDist=angleDist,stepPar=stepPar, anglePar=anglePar,covs=cov)
Defines a new ID variable, which changes where there is a long gap in
the data. This is sometimes useful for preprocessing data prior to using
prepData
.
splitAtGaps(data, maxGap = 60, shortestTrack = 0, units = "mins")
splitAtGaps(data, maxGap = 60, shortestTrack = 0, units = "mins")
data |
Data frame with (at least) columns for "ID" and "time" |
maxGap |
Longest allowed gap, in minutes by default (but see argument "units"). Track will be split at longer gaps. |
shortestTrack |
Shortest track to keep after splitting, in minutes. Shorter tracks will be removed from the output data set. |
units |
Character string, e.g., "mins" (default), "secs", "hours", "days". Time units used for the maxGap and shortestTrack arguments. |
Data frame with identical structure as input, where ID column has been replaced by new ID for split tracks. Old ID still accessible as ID_old column
For a given model, computes the probability of the process being in the different states at each time point.
stateProbs(m)
stateProbs(m)
m |
A |
The matrix of state probabilities, with element [i,j] the probability of being in state j in observation i.
Zucchini, W. and MacDonald, I.L. 2009. Hidden Markov Models for Time Series: An Introduction Using R. Chapman & Hall (London).
# m is a moveHMM object (as returned by fitHMM), automatically loaded with the package m <- example$m sp <- stateProbs(m)
# m is a moveHMM object (as returned by fitHMM), automatically loaded with the package m <- example$m sp <- stateProbs(m)
Calculates the stationary probabilities of each state, for given covariate values.
stationary(m, covs, beta = m$mle$beta)
stationary(m, covs, beta = m$mle$beta)
m |
Fitted model (as output by |
covs |
Either a data frame or a design matrix of covariates. |
beta |
Optional matrix of regression coefficients for the transition
probability model. By default, uses estimates in |
Matrix of stationary state probabilities. Each row corresponds to a row of covs, and each column corresponds to a state.
# m is a moveHMM object (as returned by fitHMM), automatically loaded with the package m <- example$m # data frame of covariates stationary(m, covs = data.frame(cov1 = 0, cov2 = 0)) # design matrix (each column corresponds to row of m$mle$beta) stationary(m, covs = matrix(c(1,0,cos(0)),1,3))
# m is a moveHMM object (as returned by fitHMM), automatically loaded with the package m <- example$m # data frame of covariates stationary(m, covs = data.frame(cov1 = 0, cov2 = 0)) # design matrix (each column corresponds to row of m$mle$beta) stationary(m, covs = matrix(c(1,0,cos(0)),1,3))
moveData
Summary moveData
## S3 method for class 'moveData' summary(object, details = TRUE, ...)
## S3 method for class 'moveData' summary(object, details = TRUE, ...)
object |
A |
details |
|
... |
Currently unused. For compatibility with generic method. |
# m is a moveData object (as returned by prepData), automatically loaded with the package data <- example$data summary(data)
# m is a moveData object (as returned by prepData), automatically loaded with the package data <- example$data summary(data)
Computation of the transition probability matrix, as a function of the covariates and the regression
parameters. Written in C++. Used in fitHMM
, logAlpha
, logBeta
,
plot.moveHMM
, pseudoRes
, and viterbi
.
trMatrix_rcpp(nbStates, beta, covs)
trMatrix_rcpp(nbStates, beta, covs)
nbStates |
Number of states |
beta |
Matrix of regression parameters |
covs |
Matrix of covariate values |
Three dimensional array trMat
, such that trMat[,,t]
is the transition matrix at
time t.
Used in prepData
.
turnAngle(x, y, z, LLangle)
turnAngle(x, y, z, LLangle)
x |
First point |
y |
Second point |
z |
Third point |
LLangle |
Logical. If TRUE, the turning angle is calculated with
|
The angle between vectors (x,y) and (y,z)
## Not run: x <- c(0,0) y <- c(4,6) z <- c(10,7) turnAngle(x,y,z,LLangle=FALSE) ## End(Not run)
## Not run: x <- c(0,0) y <- c(4,6) z <- c(10,7) turnAngle(x,y,z,LLangle=FALSE) ## End(Not run)
For a given model, reconstructs the most probable states sequence, using the Viterbi algorithm.
viterbi(m, newdata = NULL)
viterbi(m, newdata = NULL)
m |
An object |
newdata |
An object |
The sequence of most probable states.
Zucchini, W. and MacDonald, I.L. 2009. Hidden Markov Models for Time Series: An Introduction Using R. Chapman & Hall (London).
# m is a moveHMM object (as returned by fitHMM), automatically loaded with the package m <- example$m # reconstruction of states sequence states <- viterbi(m)
# m is a moveHMM object (as returned by fitHMM), automatically loaded with the package m <- example$m # reconstruction of states sequence states <- viterbi(m)
Scales each parameter from the set of real numbers, back to its natural interval. Used during the optimization of the log-likelihood.
w2n(wpar, bounds, parSize, nbStates, nbCovs, estAngleMean, stationary)
w2n(wpar, bounds, parSize, nbStates, nbCovs, estAngleMean, stationary)
wpar |
Vector of state-dependent distributions unconstrained parameters. |
bounds |
Matrix with 2 columns and as many rows as there are elements in |
parSize |
Vector of two values: number of parameters of the step length distribution, number of parameters of the turning angle distribution. |
nbStates |
The number of states of the HMM. |
nbCovs |
The number of covariates. |
estAngleMean |
|
stationary |
|
A list of:
stepPar |
Matrix of natural parameters of the step length distribution |
anglePar |
Matrix of natural parameters of the turning angle distribution |
beta |
Matrix of regression coefficients of the transition probabilities |
delta |
Initial distribution |
## Not run: nbStates <- 3 nbCovs <- 2 par <- c(0.001,0.999,0.5,0.001,1500.3,7.1) parSize <- c(1,1) bounds <- matrix(c(0,1,0,1,0,1, 0,Inf,0,Inf,0,Inf), byrow=TRUE,ncol=2) beta <- matrix(rnorm(18),ncol=6,nrow=3) delta <- c(0.6,0.3,0.1) wpar <- n2w(par,bounds,beta,delta,nbStates,FALSE) print(w2n(wpar,bounds,parSize,nbStates,nbCovs,estAngleMean=FALSE,stationary=FALSE)) ## End(Not run)
## Not run: nbStates <- 3 nbCovs <- 2 par <- c(0.001,0.999,0.5,0.001,1500.3,7.1) parSize <- c(1,1) bounds <- matrix(c(0,1,0,1,0,1, 0,Inf,0,Inf,0,Inf), byrow=TRUE,ncol=2) beta <- matrix(rnorm(18),ncol=6,nrow=3) delta <- c(0.6,0.3,0.1) wpar <- n2w(par,bounds,beta,delta,nbStates,FALSE) print(w2n(wpar,bounds,parSize,nbStates,nbCovs,estAngleMean=FALSE,stationary=FALSE)) ## End(Not run)