The hardware and bandwidth for this mirror is donated by METANET, the Webhosting and Full Service-Cloud Provider.
If you wish to report a bug, or if you are interested in having us mirror your free-software or open-source project, please feel free to contact us at mirror[@]metanet.ch.
Move*
objects
Move*
objectsThe move package provides a series of functions to import, visualize and analyse animal movement data. With the move package movement data associated with individual animals or a group of animals can be imported as a single object. The minimum requirements are timestamps, coordinates, and a unique ID per animal. Further attributes associated with the individuals (e.g. sex, age, species, etc), the locations (e.g. temperature, behavior, battery voltage, etc.), coordinates (projection) and the study (e.g. citation, license terms) can be also included.
Movebank (www.movebank.org) is a free online database of animal tracking data, where data owners can manage their data and have the option to share it with colleagues or the public. Movebank users retain ownership of their data and choose whether and with whom to share it. If the public or a registered user has permission to see a study, the study can be downloaded as a .csv file and imported directly into R. Those with Movebank accounts can also log in, browse and download data they have access to from Movebank directly within the Move package (see the ”browsing Movebank” vignette for more information). Besides importing data from Movebank, it is also possible to use personal data sets, as long as they contain the minimum required information: timestamps, coordinates, and animal IDs.
A series of functions allow you to plot, summarize, and analyse the
movement data. Individual Move
objects can be stacked to a
MoveStack
object, which includes a series of animals and the
additional data that are associated with them. This allows you to work
with multiple animals at the same time. Tracks can be also bursted by a
specific characteristic, e.g.specific behaviors (migratory,
non-migratory, resting behavior), which facilitates analyzing segments
of the trajectory belonging to a specific behavior or variable. This
package also provides functions to calculate the dynamic Brownian Bridge
Movement Model (dBBMM) and the utilization distribution(UD) of a
trajectory.
This vignette gives an overview of the main functions of the
move
library. For more functions and more details about
functions we recommend viewing the corresponding help files.
The import of tracking data into R can happen in three ways:
getMovebankData()
,
getMovebankLocationData()
,
getMovebankNonLocationData()
).move()
function.move()
function.In all cases multiple individuals can be imported simultaneously.
See the “browsing Movebank” vignette for detailed explanation.
Files from Movebank can be
imported with move(x="path")
, where “path” is a character
string that locates the file on the hard drive. The files can also be
compressed csv or .zip files downloaded with the EnvDATA tool. It may look
like:
Any data set containing movement data (position and time information) can be imported as a move object. It is mandatory to indicate the columns that contain the coordinates and the timestamp column in the correct format. The projection of the coordinates should be indicated, the data frame that includes all of the imported data can be added (if it contains additional information associated to the locations) and if the data correspond to a single individual, its name/ID can be indicated as a character, but if the data correspond to several individuals the column with the animal names/IDs has to be indicated.
Any data set containing movement data (position and time information) can be imported as a move object. You must define the columns that contain the coordinates, the projection of the coordinates, and the format of the timestamp column. If the data correspond to a single individual, its name/ID can be indicated as a character, but if the data correspond to several individuals, the column with the animal names/IDs must be indicated. If the data frame contains additional columns, these additional attributes can be included.
In the following example a move object is created from a custom file read in as a data frame:
data <- read.csv(system.file("extdata","leroy.csv.gz",package="move"))
leroy <- move(x=data$location.long, y=data$location.lat,
time=as.POSIXct(data$timestamp, format="%Y-%m-%d %H:%M:%OS", tz="UTC"),
proj=CRS("+proj=longlat +ellps=WGS84 +datum=WGS84"),
data=data, animal=data$individual.local.identifier, sensor=data$sensor)
By entering the object an overview is provided. This overview
indicates the object class (Move
), number of coordinates,
the extent of the coordinates, the coordinates reference system
(projection), the number of columns of the imported data.frame, and the
names of the columns of the data.frame, as well as the first and last
timestamp and the duration of the observation.
## class : Move
## features : 919
## extent : -73.93067, -73.84366, 42.70898, 42.7687 (xmin, xmax, ymin, ymax)
## crs : +proj=longlat +datum=WGS84 +no_defs
## variables : 17
## names : timestamp, eobs.battery.voltage, eobs.horizontal.accuracy.estimate, eobs.key.bin.checksum, eobs.speed.accuracy.estimate, eobs.start.timestamp, eobs.status, eobs.temperature, eobs.type.of.fix, eobs.used.time.to.get.fix, ground.speed, heading, height.above.ellipsoid, utm.easting, utm.northing, ...
## min values : 1234354605, 3596, 3.07, 3258904, 0.27, 2009-02-11 12:14:59.000, A, 13, 3, 4, 0.01, 0, -169.6, 587507.837877134, 4729143.16566605, ...
## max values : 1236158219.998, 3666, 97.02, 4291715164, 33.04, 2009-03-04 09:15:01.000, A, 35, 3, 119, 31.71, 359.79, 349, 594679.382402579, 4735720.47868847, ...
## timestamps : 2009-02-11 12:16:45 ... 2009-03-04 09:16:59.998 Time difference of 21 days (start ... end, duration)
## sensors : gps
## indiv. data : eobs.fix.battery.voltage, manually.marked.outlier, visible, sensor.type, individual.taxon.canonical.name, tag.local.identifier, individual.local.identifier, study.name, study.timezone
## indiv. value: 0 NA true gps Martes pennanti 74 Leroy Urban fisher GPS tracking Eastern Standard Time
## unused rec. : 1071
## study name : Urban fisher GPS tracking
## date created: 2023-01-20 16:56:43
The following objects containing movement data can also be read in by
the move()
function:
ltraj
from the library
adehabitatLTtelemetry
or a list of telemetry
objects from the library ctmmtrack
from the library
bcpatrack_xyt
from the library
amtbinClstPath
or binClstStck
from the library EMbCdata.frame
containing data downloaded
from Movebank webpage or with getMovebankLocationDataIf the data are from Movebank, the function will automatically
recognize that there are multiple individuals contained in the file. If
the custom dataset includes more than one animal, you must indicate the
column with the animal IDs in the argument animal
. When
multiple individuals are imported, the move()
function
automatically creates a stack of Move
objects called
MoveStack
. A MoveStack
can also be created
from a list of multiple Move
or MoveStack
objects with the function moveStack()
. In this case all
objects have to be in the same projection, and their timestamps have to
be in the same time zone. Most functions of the Move package are capable
of working with MoveStacks
.
In the following example a MoveStack
is created from two
Move
objects:
ricky<-move(system.file("extdata","ricky.csv.gz", package="move"))
data(leroy)
## if argument forceTz is not stated, the timestamp is converted to the computer timezone
leroyP<-spTransform(leroy, proj4string(ricky))
myStack <- moveStack(list(leroyP, ricky),forceTz="UTC")
In the overview the information of all individuals is combined.
## class : MoveStack
## features : 9877
## extent : -73.94032, -73.84366, 42.70898, 42.851 (xmin, xmax, ymin, ymax)
## crs : +proj=longlat +datum=WGS84 +no_defs
## variables : 20
## names : timestamp, eobs.battery.voltage, eobs.horizontal.accuracy.estimate, eobs.key.bin.checksum, eobs.speed.accuracy.estimate, eobs.start.timestamp, eobs.status, eobs.temperature, eobs.type.of.fix, eobs.used.time.to.get.fix, ground.speed, heading, height.above.ellipsoid, utm.easting, utm.northing, ...
## min values : 1234354605, 3449, 2.05, 524362, 0.22, 2009-02-11 12:14:59.000, A, -12, 3, 3, 0, 0, -315.6, 586587.171471437, 4729143.16566605, ...
## max values : 1270056686, 3740, 97.02, 4294919150, 48.67, 2010-03-31 17:30:02.000, A, 35, 3, 120, 43.12, 359.79, 436.6, 594679.382402579, 4744838.73355732, ...
## timestamps : 2009-02-11 12:16:45 ... 2010-03-31 17:31:26 Time difference of 413 days (start ... end, duration)
## sensors : gps
## indiv. data : eobs.fix.battery.voltage, manually.marked.outlier, visible, sensor.type, individual.taxon.canonical.name, tag.local.identifier, individual.local.identifier, study.name, study.timezone, behavioural.classification
## min ID Data : NA, NA, true, gps, Martes pennanti, 74, Leroy, Urban fisher GPS tracking, Eastern Standard Time, NA
## max ID Data : NA, NA, true, gps, Martes pennanti, 1016, Ricky.T, Urban fisher GPS tracking, Eastern Standard Time, NA
## individuals : Leroy, Ricky.T
## unused rec. : 2959
## date created: 2024-11-20 16:55:07
It is also possible to break down a MoveStack
into a
list of Move
objects using the split()
function. This is often useful to do calculations with functions that
require lists as an input as e.g. lapply()
. This list of
Move
objects can be transformed easily back into a
MoveStack
with the function moveStack()
.
Given that an animal cannot be at two places at the same time,
Move*
objects account for this and cannot be created if
there are duplicated timestamps present in the data. There are several
options to remove these duplicated timestamps:
If the study is from Movebank and it contains duplicated timestamps
you can set the argument removeDuplicatedTimestamps=TRUE
when reading in the .csv file with the move()
function.
This will retain the first of multiple records with the same animal ID
and timestamp, and remove any subsequent duplicates. In some cases, one
duplicate record might contain more complete information or a better
location estimate than the other. In case you want to control which of
the duplicate timestamps are kept and which are deleted, we recommend to
download the data as a .csv file from Movebank or to use the function
getMovebankLocationData()
, find the duplicates using
e.g. getDuplicatedTimestamps()
, decide which of the
duplicated timestamp to retain, and than create a move/moveStack object
with the function {move()
. Another option is to edit the
records in movebank and mark the appropriate records as outliers.
One option is to remove the duplicated timestamps with the function
duplicated()
, nevertheless in this case, depending on the
settings, the first or the last location will be retained. In some
cases, one duplicate record might contain more complete information or a
better location estimate than the other. In case you want to control
which of the duplicate timestamps are kept and which are deleted, you
can find the duplicates using
e.g. getDuplicatedTimestamps()
, decide which of the
duplicated timestamp to retain, and than create a move/moveStack object
with the function move()
.
Move*
objectsThe Move*
objects contains a series of slots containing
the timestamps, coordinates and additional information associated with
the individuals and the locations. To get an overview of all slots use
str()
:
There are a set of functions to extract the information contained in
Move*
objects:
## [1] "2009-02-11 12:16:45 UTC" "2009-02-11 12:31:38 UTC"
## [3] "2009-02-11 12:45:48 UTC"
## For a MoveStack the output is a continuous vector of timestamps of all individuals:
timestamps(myStack)[1:3]
## [1] "2009-02-11 12:16:45 UTC" "2009-02-11 12:31:38 UTC"
## [3] "2009-02-11 12:45:48 UTC"
## location.long location.lat
## 45 -73.89880 42.74370
## 46 -73.89872 42.74369
## 47 -73.89869 42.74364
## For a MoveStack the output is one matrix containing the coordinates of all individuals:
coordinates(myStack)[1:3,]
## location.long location.lat
## 45 -73.89880 42.74370
## 46 -73.89872 42.74369
## 47 -73.89869 42.74364
Move*
object with non-Movebank data:## [1] "+proj=longlat +datum=WGS84 +no_defs"
## class : Extent
## xmin : -73.93067
## xmax : -73.84366
## ymin : 42.70898
## ymax : 42.7687
## min max
## location.long -73.93067 -73.84366
## location.lat 42.70898 42.76870
## eobs.fix.battery.voltage manually.marked.outlier visible sensor.type
## Leroy 0 NA true gps
## individual.taxon.canonical.name tag.local.identifier
## Leroy Martes pennanti 74
## individual.local.identifier study.name
## Leroy Leroy Urban fisher GPS tracking
## study.timezone
## Leroy Eastern Standard Time
## eobs.fix.battery.voltage manually.marked.outlier visible sensor.type
## Leroy 0 NA true gps
## Ricky.T NA NA true gps
## individual.taxon.canonical.name tag.local.identifier
## Leroy Martes pennanti 74
## Ricky.T Martes pennanti 1016
## individual.local.identifier study.name
## Leroy Leroy Urban fisher GPS tracking
## Ricky.T Ricky.T Urban fisher GPS tracking
## study.timezone behavioural.classification
## Leroy Eastern Standard Time NA
## Ricky.T Eastern Standard Time NA
## timestamp eobs.battery.voltage eobs.horizontal.accuracy.estimate
## 45 2009-02-11 12:16:45 3615 14.85
## 46 2009-02-11 12:31:38 3623 5.38
## 47 2009-02-11 12:45:48 3627 5.38
## eobs.key.bin.checksum eobs.speed.accuracy.estimate eobs.start.timestamp
## 45 2992317972 5.65 2009-02-11 12:14:59.000
## 46 1723246055 4.69 2009-02-11 12:30:01.000
## 47 1910450098 4.19 2009-02-11 12:45:01.000
## eobs.status eobs.temperature eobs.type.of.fix eobs.used.time.to.get.fix
## 45 A 20 3 106
## 46 A 26 3 97
## 47 A 24 3 47
## ground.speed heading height.above.ellipsoid utm.easting utm.northing
## 45 2.10 125.17 79.3 590130.0 4732942
## 46 0.51 3.28 94.2 590136.3 4732940
## 47 0.16 91.10 82.5 590138.5 4732935
## utm.zone study.local.timestamp
## 45 18N 2009-02-11 07:16:45
## 46 18N 2009-02-11 07:31:38
## 47 18N 2009-02-11 07:45:48
## [1] "Leroy"
## [1] "Leroy" "Ricky.T"
## [1] 919
## Leroy Ricky.T
## 919 8958
## [1] 1
## [1] 2
Move*
object with
non-Movebank data:## [1] "gps"
## timestamp location.long location.lat eobs.battery.voltage
## 1 2009-02-11 00:00:01 NA NA 3642
## 2 2009-02-11 00:15:01 NA NA 3635
## 3 2009-02-11 00:30:01 NA NA 3632
## eobs.horizontal.accuracy.estimate eobs.key.bin.checksum
## 1 NA 3717396718
## 2 NA 3757665008
## 3 NA 30534636
## eobs.speed.accuracy.estimate eobs.start.timestamp eobs.status
## 1 NA 2009-02-11 00:00:01.000 D
## 2 NA 2009-02-11 00:15:01.000 D
## 3 NA 2009-02-11 00:30:01.000 D
## eobs.temperature eobs.type.of.fix eobs.used.time.to.get.fix ground.speed
## 1 27 0 120 NA
## 2 23 0 120 NA
## 3 21 0 120 NA
## heading height.above.ellipsoid utm.easting utm.northing utm.zone
## 1 NA NA NA NA
## 2 NA NA NA NA
## 3 NA NA NA NA
## study.local.timestamp sensor timestamps
## 1 2009-02-10 19:00:01 gps 2009-02-11 00:00:01
## 2 2009-02-10 19:15:01 gps 2009-02-11 00:15:01
## 3 2009-02-10 19:30:01 gps 2009-02-11 00:30:01
## [1] "2009-02-11 00:00:01 UTC" "2009-02-11 00:15:01 UTC"
## [3] "2009-02-11 00:30:01 UTC"
## timestamp location.long location.lat eobs.battery.voltage
## 1 2009-02-11 00:00:01 NA NA 3642
## 2 2009-02-11 00:15:01 NA NA 3635
## 3 2009-02-11 00:30:01 NA NA 3632
## eobs.horizontal.accuracy.estimate eobs.key.bin.checksum
## 1 NA 3717396718
## 2 NA 3757665008
## 3 NA 30534636
## eobs.speed.accuracy.estimate eobs.start.timestamp eobs.status
## 1 NA 2009-02-11 00:00:01.000 D
## 2 NA 2009-02-11 00:15:01.000 D
## 3 NA 2009-02-11 00:30:01.000 D
## eobs.temperature eobs.type.of.fix eobs.used.time.to.get.fix ground.speed
## 1 27 0 120 NA
## 2 23 0 120 NA
## 3 21 0 120 NA
## heading height.above.ellipsoid utm.easting utm.northing utm.zone
## 1 NA NA NA NA
## 2 NA NA NA NA
## 3 NA NA NA NA
## study.local.timestamp
## 1 2009-02-10 19:00:01
## 2 2009-02-10 19:15:01
## 3 2009-02-10 19:30:01
## [1] "gps"
Move*
was created by directly downloading the
data from Movebank via the getMovebankData()
function,
these slots will contain information if provided by the data owner:coatis_bci <- getMovebankData(study="Coatis on BCI Panama (data from Powell et al. 2017)")
## the study name:
coatis_bci@study
# [1] "Coatis on BCI Panama (data from Powell et al. 2017)"
## how to cite the study:
citations(coatis_bci)
# [1] "Powell RA, Ellwood S, Kays R, Maran T (2017) Stink or swim: techniques to meet the challenges for the study and conservation of small critters that hide, swim, or climb, and may otherwise make themselves unpleasant. In Macdonald DW, Newman C, Harrington LA, eds, Biology and Conservation of Musteloids. Oxford University Press, Oxford. p 216–230. doi:10.1093/oso/9780198759805.003.0008"
## license terms of the study
licenseTerms(coatis_bci)
# [1] "These data have been published by the Movebank Data Repository with DOI "http://dx.doi.org/10.5441/001/1.41076dq1"
MoveStack
objectThe MoveStack
object is very similar to the
Move
object but contains an additional slot
(@trackId
) which ensures the correct association of the
information belonging to each individual.
Obtain a vector with the individual names/IDs which each location is associated:
## [1] Leroy Leroy Leroy
## Levels: Leroy Ricky.T
MoveStack
objectExtract one specific individual:
Extract the first individual:
Extract several specific individuals (Note: individuals have to be
listed in the same order as they are in the moveStack
):
Extract several individuals:
MoveStack
objectDelete one specific individual:
Delete the first individual:
Delete several specific individuals (Note: individuals have to be
listed in the same order as they are in the moveStack
):
Delete several individuals:
Move*
objectSubset a Move
object for specific locations:
Subset a MoveStack
object for specific locations:
## select the locations 1-50 from a movestack. WARNING: this will just select the 50 first locations in order of occurrence, in this case they correspond to the first individual of the movestack
myStackSub <- myStack[1:50]
## to select locations 1-50 from each individual
myStackSubs <- moveStack(lapply(split(myStack), function(x){x[1:50]}))
Subset a Move*
object to a specific time range, for
example:
## selecting a specific day
leroyOneDay <- leroy[as.Date(timestamps(leroy))==as.Date("2009-02-25")]
## selecting a range of days
leroy3Days <- leroy[as.Date(timestamps(leroy))%in%c(as.Date("2009-02-25"):as.Date("2009-02-28"))]
## selecting a specific month
myStackMarch <- myStack[month(timestamps(myStack))==3]
GPS data are often collected in UTC, but for some analysis, it is
useful to transform the timestamps into the local time zone. This can be
done with the function with_tz()
of the library
lubridate
:
## [1] "2009-02-11 12:16:45 UTC"
leroyLocal <- leroy
timestamps(leroyLocal) <- with_tz(timestamps(leroyLocal), tz="America/New_York")
leroyLocal@timestamps[1]
## [1] "2009-02-11 07:16:45 EST"
When the Move of MoveStack object was created a projection was
declared, therefore these objects can be reprojected into any projection
with the function spTransform()
:
## [1] "+proj=longlat +datum=WGS84 +no_defs"
leroyProj <- spTransform(leroy, CRSobj="+proj=utm +zone=18 +ellps=GRS80 +datum=NAD83 +units=m +no_defs")
projection(leroyProj)
## [1] "+proj=utm +zone=18 +datum=NAD83 +units=m +no_defs"
The functions plot()
, points()
,
lines()
can be used with any Move*
object and
all graphical parameters can be passed along.
plot(leroy, xlab="Longitude", ylab="Latitude", type="l", pch=16, lwd=0.5)
points(leroy, pch=20, cex=0.5)
plot(myStack, xlab="Longitude", ylab="Latitude",type="b", pch=16, cex=0.5, col=c("blue","magenta")[myStack@trackId])
Also ggplot()
can be used, but in this case the
Move*
object has to be transformed into a
data.frame
.
library(ggplot2)
myStackDF <- as.data.frame(myStack)
ggplot(data = myStackDF, aes(x = location.long, y = location.lat, color = trackId)) +
geom_path() + geom_point(size = 0.5) + theme_bw() + coord_cartesian()
The tracks can also be plotted on a background map like open street map:
ggmap
library(ggmap)
require(mapproj)
leroyDF <- as.data.frame(leroy)
register_stadiamaps(your_stadia_key)
m <- get_map(bbox(extent(leroy)*1.1),maptype="stamen_terrain", source="stadia", zoom=12)
ggmap(m)+geom_path(data=leroyDF, aes(x=location.long, y=location.lat))
leaflet
: for an example see the vignette
“leafletPlot”.There are a series of functions that extract temporal (e.g. time lag
between locations) and spatial (e.g. distance between locations)
information from Move*
objects.
## [1] 14.88333 14.18330 14.45007 15.04998 14.90000
## List of 2
## $ Leroy : num [1:918] 14.9 14.2 14.5 15 14.9 ...
## $ Ricky.T: num [1:8957] 27.55 2.55 1.883 2.217 0.583 ...
## [1] 6.382461 5.606722 12.671263 9.651801 11.733921
## List of 2
## $ Leroy : num [1:918] 6.38 5.61 12.67 9.65 11.73 ...
## $ Ricky.T: num [1:8957] 153.35 12.62 7.81 10.23 13.18 ...
## [1] 0.007147213 0.006588408 0.014615000 0.010688607 0.013125191
## List of 2
## $ Leroy : num [1:918] 0.00715 0.00659 0.01461 0.01069 0.01313 ...
## $ Ricky.T: num [1:8957] 0.0928 0.0825 0.0691 0.0769 0.3766 ...
## [1] 101.44448 157.41345 26.47859 -133.22121 -108.37639
## List of 2
## $ Leroy : num [1:918] 101.4 157.4 26.5 -133.2 -108.4 ...
## $ Ricky.T: num [1:8957] 115.7 70.2 -18.2 -78.9 76.5 ...
## [1] 55.96892 -130.93488 -159.69985 24.84488 -179.16947
## List of 2
## $ Leroy : num [1:917] 56 -130.9 -159.7 24.8 -179.2 ...
## $ Ricky.T: num [1:8956] -45.5 -88.4 -60.7 155.4 -108.9 ...
It can be interesting to compare different parts of the track with
each other. For example, how do data points differ between winter and
summer, or between behaviors like migrating and non-migrating, or
between day and night. To indicate which point of the data set belongs
to which category, a track is ’bursted’ with the function
burst()
. A track is bursted by supplying a vector with the
specific behavior/variable associated with each location. The length of
this vector has to be one shorter than the number of coordinates. The
returned object belongs to the class MoveBurst
.
In the following example, the trajectory of ‘leroy’ is bursted according to day and night.
library(solartime)
elev<-computeSunPosition(timestamps(leroy), coordinates(leroy)[,2],coordinates(leroy)[,1])[,'elevation']
DayNight <- c("Night",'Day')[1+(elev>0)]
## assigning to each segment if it is during daytime or night
leroy.burst <- burst(x=leroy, f=DayNight[-n.locs(leroy)])
The MoveBurst
object has a very similar structure to a
Move
object. The MoveBurst
object has an
additional slot @burstId
which contains the information of
the assignment of each location to each behavior/category (see
str(leroy.burst)
). In the overview of the
MoveBurst
there is an additional information feature under
“bursts” with the information how many locations fall within each
category. MoveBurst
objects cannot be stacked into a
MoveStack
.
The resulting MoveBurst
can be plotted highlighting the
different categories in different colors. Note that R by default only
uses 8 different colors, if there are more than 8 categories colors will
be recycled, in this case we recommended specifying the colors.
plot(leroy.burst, type="l", col=c("red", "black"), asp=1)
legend("bottomleft",legend=c("day","night"), col=c("red", "black"), lty=1)
The plotBursts()
function plots a circle at the midpoint
of each burst segment (consecutive coordinates that belong to a single
burst). The size of the cycles can have different meanings, depending on
what is calculated in the sizeFUN
argument (for details see
?plotBursts
). The default size refers to the relative time
of the burst segment compared to the whole track.
### in the default the size of the cicles is relative to the total time spent within each burst segment
plotBursts(leroy.burst,breaks=5, col=c("red", "black"), pch=19, add=F,main="Size of points: total time spent in burst segment", asp=1)
legend("bottomleft",legend=c("day","night"), col=c("red", "black"), pch=19)
## here, e.g., the size of the circles is relative to the total distance moved within each burst segment
plotBursts(object=leroy.burst, breaks=5, sizeFUN=function(x){sum(distance(x))},col=c("red", "black"), pch=19, add=F,main="Size of points: total distance moved in burst segment", asp=1)
legend("bottomleft",legend=c("day","night"), col=c("red", "black"), pch=19)
A MoveBurst
object can be split into a list of
Move
objects, one per each burst segment:
leroyB.split <- split(leroy.burst)
par(mfrow=c(1,4))
plot(leroyB.split[[1]], type="b", pch=20, main=paste0(names(leroyB.split[1])," 1"))
plot(leroyB.split[[2]], type="b", pch=20, main=paste0(names(leroyB.split[2])," 1"))
plot(leroyB.split[[3]], type="b", pch=20, main=paste0(names(leroyB.split[3])," 2"))
plot(leroyB.split[[4]], type="b", pch=20, main=paste0(names(leroyB.split[4])," 2"))
The corridor()
function identifies sections of the track
that suggest corridor use behavior (i.e. parallel and fast movements).
For details see LaPoint et al. (2013) Animal Behavior, Cost-based
Corridor Models, and Real Corridors. Landscape Ecology. https://doi.org/10.1007/s10980-013-9910-0.
The speed threshold and how parallel the segments should be can be adjusted in the function. The proportion of speeds that are high enough to be a valid corridor point is be defined in the argument “speedProp”. The default value selects speeds that are greater than 75% of all speeds. The proportion of the circular variances low enough to be a valid corridor point is defined in the argument “circProp”. Low values of the circular variance indicate that the segments are (near) parallel. The default value selects variances that are lower than 25 % of all variances.
LeroyCorr <- corridor(leroy)
plot(LeroyCorr, type="l", xlab="Longitude", ylab="Latitude", col=c("black", "grey"), lwd=c(1,2))
legend("bottomleft", c("Corridor", "Non-corridor"), col=c("black", "grey"), lty=c(1,1), bty="n")
The resulting object is of class MoveBurst
which in the
@data
slot contains new attributes associated to the
locations, where “segMid_x” and “segMid_y” are the coordinates of the
midpoint of each segment, and “speed”, “azimuth”, “pseudoAzimuth” and
“circVar” are the variables used to identify the segments that suggest
corridor behavior.
The dBBMM is a method to calculate the occurrence distribution of a given track, i.e. it estimates where the animal was located during the observed (tracking) period. It calculates the probability landscape for the transition between any two known consecutive locations given the amount of time it had available (assuming a conditional random (Brownian) motion between locations). It is “dynamic” because it allows the variance to change along the trajectory, using the behavioral change point analysis in a sliding window. For details see Kranstauber et. all (2012) A dynamic Brownian bridge movement model to estimate utilization distributions for heterogeneous animal movement. Journal of Animal Ecology. https://doi.org/10.1111/j.1365-2656.2012.01955.x
The Move*
object used to calculate the dBBMM needs to be
projected and, the location error, margin, window size and raster
options need to be specified. Note that in most cases different values
of window size and the margin do not have a great effect on the results.
For details on the different ways to specify the raster options see the
‘Details’ section of ?brownian.bridge.dyn
. In the following
examples, the calculation is done on a raster with a pixel size of
100m:
leroyRed <- leroy[1:200] # reducing dataset for example
leroy.prj <- spTransform(leroyRed, center=TRUE) # center=T: the center of the coordinate system is the center of the track. Units are in meters
dBB.leroy <- brownian.bridge.dyn(leroy.prj, ext=.85, raster=100, location.error=20)
The resulting object will be of class DBBMM
,
DBBMMStack
, or DBBMMBurstStack
depending on
the input class. All of them contain a raster or a rasterStack were all
pixels of each raster sum up to 1 (besides the
DBBMMBurstStack
, for details see ‘Calculating a dBBMM from
a MoveBurst’ section)
Plotting the result provides limited visual information, as the values correspond to the probability the animal was present in a given pixel during the observed period. These results can be used to estimate probability of the animal being present in a certain area during the tracking period (see ‘Utilization distribution’ section), or to calculate the probability of the animal being present in a certain environment during the tracking period by overlaying the resulting raster over a categorical raster (e.g. of land cover) and summing up the dBBMM pixel values that fall within each category.
From the Dynamic Brownian Bridge object, the UD can be calculated, which represents the minimum area in which an animal has some specified probability of being located.
UDleroy <- getVolumeUD(dBB.leroy)
par(mfrow=c(1,2))
plot(UDleroy, main="UD")
## also a contour can be added
plot(UDleroy, main="UD and contour lines")
contour(UDleroy, levels=c(0.5, 0.95), add=TRUE, lwd=c(0.5, 0.5), lty=c(2,1))
The obtained UD can be subset to a specific probability:
par(mfrow=c(1,3))
## mantaining the lower probabilities
ud95 <- UDleroy
ud95[ud95>.95] <- NA
plot(ud95, main="UD95")
## or extracting the area with a given probability, where cells that belong to the given probability will get the value 1 while the others get 0
ud95 <- UDleroy<=.95
plot(ud95, main="UD95")
ud50 <- UDleroy<=.5
plot(ud50, main="UD50")
MoveBurst
The dBBMM can also be calculated for a MoveBurst
, for
all burst types, or only for a subset of them. The result will provide
one raster per burst segment in an object of class
DBBMMBurstStack
. By default it is calculated for all burst
types:
leroyBRed <- leroy.burst[1:155] # reducing dataset for example
leroyB.prj <- spTransform(leroyBRed, center=TRUE) # center=T: the center of the coordinate system is the center of the track. Units are in meters
dBB.leroyB <- brownian.bridge.dyn(leroyB.prj, ext=1.25, raster=100, location.error=20)
plot(dBB.leroyB)
The UD can also be calculated for each raster layer of the resulting
DBBMMBurstStack
. In this case to calculate the UD, each
raster layer needs to be standardized so the sum of all pixels is 1.
This can be done with the function UDStack()
To calculate the UD for each burst type, the layers of each burst
type have to be summed up and converted to the .UD
class.
par(mfrow=c(1,2))
## select all layers corresponding to "Day", sum them up, and convert them to class .UD
daylayers <- grep("Day", names(dBB.leroyB))
rasterDay <- sum(dBB.leroyB[[daylayers]])
rasterDayUD <- as(rasterDay,".UD")
UDleroy.day <- getVolumeUD(rasterDayUD)
plot(UDleroy.day, main="UD of 'Day'")
## same for layer corresponding to "Night"
nightlayers <- grep("Night", names(dBB.leroyB))
rasterNight <- sum(dBB.leroyB[[nightlayers]])
rasterNightUD <- as(rasterNight,".UD")
UDleroy.night <- getVolumeUD(rasterNightUD)
plot(UDleroy.night, main="UD of 'Night'")
To calculate the dBBMM for a specific or a subset of burst types, these have to be stated in the argument “burstType”. In this case the dBBMM is calculated only taking into account the burst segments “Night”:
dBB.leroyB.night <- brownian.bridge.dyn(leroyB.prj, ext=.75, raster=100, location.error=20, burstType="Night")
plot(dBB.leroyB.night,nr=1)
From the resulting object, the UD per layer, or the total UD can be calculated as above.
Below some common issues that may arise while calculating the dBBMM and possible solutions to them.
There can be several reasons why the calculation is taking very long:
By default the dBBMM will take the shortest time lag (in minutes) of the trajectory divided by 15 as the time interval taken for every integration step. If the shortest time lag is e.g. 15min, than the trajectory will be divided into steps of 1min. But if for some reason there is one time lag of 1sec (0.0167mins), the trajectory will be divided into steps of 0.001 mins, increasing the calculation time immensely.
Using the argument time.step
the integration step time
can be set manually. As a general rule of thumb one can take the
intended GPS schedule in minutes divided by 15.
## check the timeLag of data. If there are timelags shorter than the intended scedule, use the argument "timestep"
rickyRed <- ricky[1:100]
summary(timeLag(rickyRed,"mins"))
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 0.200 1.875 2.117 9.422 6.883 60.700
ricky.prj <- spTransform(rickyRed, center=TRUE)
ts <- median(timeLag(rickyRed,"mins"))
BB.ricky <- brownian.bridge.dyn(ricky.prj, ext=.45, dimSize=100, location.error=20,time.step=ts/15)
For trajectories collected at a very high resolution, e.g. 1Hz, or
trajectories containing high resolution bursts, the
time.step
argument can be set to 0.017 min (1sec), as there
is no need to estimate where the animal was at a shorter time interval,
because it is known where it was every 1sec.
In any case giving higher values to the time.step
argument will speed up the calculation time.
This problem is mostly due to data that contains large chunks of missing data. During these large time gaps, the uncertainty of where the animal may have been is very large, or in other words the animal could have been anywhere.
Here is an example:
## creating a gappy data set
leroyWithGap <- leroy[-c(50:500,550:850)]
leroyWithGap_p <- spTransform(leroyWithGap, center=TRUE)
## calculate the dBBMM with the default extent gives an error that it is too small
dbb <- brownian.bridge.dyn(leroyWithGap_p, raster=100, location.error=20)
## Error in .local(object, raster, location.error = location.error, ext = ext, : Lower x grid not large enough, consider extending the raster in that direction or enlarging the ext argument
## making the extent bigger seems to solve the problem
dbb <- brownian.bridge.dyn(leroyWithGap_p, raster=100, location.error=20, ext=4)
## but than the UD is not very informative
ud <- getVolumeUD(dbb)
par(mfrow=c(1,2))
plot(ud, main="UD")
contour(ud, levels=c(0.5, 0.95), add=TRUE, lwd=c(0.5, 0.5), lty=c(2,1))
plot(ud, main="UD with locations")
points(leroyWithGap_p, col="red", cex=.5, pch=20)
contour(ud, levels=c(0.5, 0.95), add=TRUE, lwd=c(0.5, 0.5), lty=c(2,1))
And the solution of such a case:
The solution is to remove the variance of the segments corresponding
to the large time gaps. For this first the variance is calculated with
the brownian.motion.variance.dyn()
function, then the
segments corresponding to the large time gaps are set to FALSE, and
finally the dBBMM is calculated.
## calculate the dynamic brownian motion variance of the gappy track
dbbv <- brownian.motion.variance.dyn(leroyWithGap_p, location.error=20, window.size=31, margin=11)
## the intended GPS fix rate of leroy was 15min, so we will ignore for example all segments that have a larger time lag than 5hours. The 'dBMvariance' object resulting from the function above, contains the slot '@interest' in which those segments marked as FALSE won't be included in the calculation of the dBBMM. Therefore we set all segments with time lag larger than 300mins to false
dbbv@interest[timeLag(leroyWithGap_p,"mins")>300] <- FALSE
## then we use the 'dBMvariance' object to calculate the dBBMM
dbb.corrected <- brownian.bridge.dyn(dbbv, raster=100, ext=.45,location.error=20)
## now the UD makes more sense
ud.corrected <- getVolumeUD(dbb.corrected)
par(mfrow=c(1,2))
plot(ud.corrected, main="UD")
contour(ud.corrected, levels=c(0.5, 0.95), add=TRUE, lwd=c(0.5, 0.5), lty=c(2,1))
plot(ud.corrected, main="UD with locations")
points(leroyWithGap_p, col="red", cex=.5, pch=20)
contour(ud.corrected, levels=c(0.5, 0.95), add=TRUE, lwd=c(0.5, 0.5), lty=c(2,1))
The earth mover’s distance function emd()
quantifies
similarity between utilization distributions by calculating the effort
it takes to shape one utilization distribution landscape into
another.
## e.g. compare the utilization distributions of the day and night UDs of Leroy
dBB.leroyB <- brownian.bridge.dyn(leroyB.prj, ext=1.5, raster=500, location.error=20)
leroyBud <- UDStack(dBB.leroyB)
## to optimize the calculation, the cells outside of the 99.99% UD contour are removed by setting them to zero.
values(leroyBud)[values(getVolumeUD(leroyBud))>.999999]<-0
## the rasters have to be standardized so the pixels sum up to 1
leroyBud_2<-(leroyBud/cellStats(leroyBud,sum))
emd(leroyBud_2)
## Day.1 Night.1 Day.2
## Night.1 3316.270
## Day.2 4531.383 1804.934
## Night.2 2499.574 1095.517 2247.205
The function interpolateTime()
allows to interpolate
trajectories. The interpolation can be done to obtain a regular track
with a specific number of locations, to obtain a track with location at
a specific regular time interval, or to obtain positions at given
timestamps. The new locations can be interpolated using an euclidean,
great circle or rhumb line function.
## providing the number of locations. In this case the interpolated trajectory will have 500 locations with a regular timelag
interp500p <- interpolateTime(leroyRed, time=500, spaceMethod='greatcircle')
plot(leroyRed, col="red",pch=20, main="By number of locations")
points(interp500p)
lines(leroyRed, col="red")
legend("bottomleft", c("True locations", "Interpolated locations"), col=c("red", "black"), pch=c(20,1))
summary(timeLag(interp500p, "mins"))
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 11.54 11.54 11.54 11.54 11.54 11.54
## providing a time interval. In this case the interpolated trajectory will have a location every 5mins
interp5min <- interpolateTime(leroyRed, time=as.difftime(5, units="mins"), spaceMethod='greatcircle')
plot(leroyRed, col="red",pch=20, main="By time interval")
points(interp5min)
lines(leroyRed, col="red")
legend("bottomleft", c("True locations", "Interpolated locations"), col=c("red", "black"), pch=c(20,1))
## providing a vector of timestamps
ts <- as.POSIXct(c("2009-02-12 06:00:00", "2009-02-12 12:00:00", "2009-02-12 23:00:00",
"2009-02-14 06:00:00", "2009-02-14 12:00:00", "2009-02-14 23:00:00"),
format="%Y-%m-%d %H:%M:%S", tz="UTC")
interpCusom <- interpolateTime(leroyRed, time=ts, spaceMethod='greatcircle')
plot(leroyRed, col="red",pch=20, main="By custom timestamps")
points(interpCusom)
lines(leroyRed, col="red")
legend("bottomleft", c("True locations", "Interpolated locations"), col=c("red", "black"), pch=c(20,1))
In the previous case (with interpolateTime()
), tracks
are regularized by “making up” new locations. With the function
thinTrackTime()
trajectories can be thinned, and somewhat
regularized by selecting locations that are separated at a specific time
interval. Trajectories can be also thinned by selecting locations that
are separated by a specific distance with the function
thinDistanceAlongTrack()
.
The function searches for consecutive segments with a cumulative sum
of the time corresponding to interval and tolerance values. From each
selected chunk of the track, only the first and last location are kept
in the new object and this new segment is labelled with “selected”. The
segments labelled as “notSelected” are those parts of the track that did
not fulfill the indicated interval, because e.g. there is a large time
gap without locations in the data. Consecutive segments that are larger
than the indicated interval get condensed into one “notSelected”
segment. The returned object is a MoveBurst
.
## selecting those segments that have a time interval of 45 mins plus/minus 5 mins. The original data have a fix every ~15mins.
thintime <- thinTrackTime(leroyRed, interval = as.difftime(45, units='mins'),
tolerance = as.difftime(5, units='mins'))
Looking at the correspondence between time lag and burstID, we can see that “notSelected” bursts correspond to segments that have a larger timelag than ~45mins. These “notSelected” bursts can correspond to multiple consecutive segments that have a larger timelag than ~45min, or single large time gaps that are present in the original data.
## TL burst
## 15 45.05000 selected
## 16 45.21663 selected
## 17 44.55000 selected
## 18 45.00003 selected
## 19 45.03332 selected
## 20 420.68335 notSelected
## 21 45.01665 selected
## 22 44.88337 selected
## 23 45.36665 selected
## 24 779.53335 notSelected
## 25 45.11663 selected
The function searches for consecutive segments with a cumulative sum
of distances traveled corresponding to interval and tolerance values.
From each selected chunk of the track, only the first and last location
are kept in the new object, this new segment is labelled with
“selected”. The segments labelled as “notSelected” are those segments
that are larger than the distance indicated. These long segments are
present in the original data. The returned object is a
MoveBurst
.
Note that in the case of thinDistanceAlongTrack()
, the
distances between the locations in the new object do not represent the
distance that the animal actually traveled, as the intermediate
locations have been removed.
Move*
objectsA Move*
object can be transformed into:
Data.frame
SpatialPointsDataFrame
ltraj
objectStoring and loading as .RData
## save the Move* object as a RData file
save(leroy, file="C:/User/Documents/leroy.RData")
## load an Move* object saved as a RData file
load("C:/User/Documents/leroy.RData")
Storing as text file
## save as a text file
leroyDF <- as.data.frame(leroy)
write.table(leroyDF, file="C:/User/Documents/leroyDF.csv", sep=",", row.names = FALSE)
Storing as ESRI shapefile file
## save as a shape file
library(rgdal)
writeOGR(leroy, "C:/User/Documents/leroySHP/", layer="leroy", driver="ESRI Shapefile")
Storing as kml or kmz
## kml or kmz of movestack ##
library(plotKML)
## open a file to write the content
kml_open('myStack.kml')
## write the movement data individual-wise
for(i in levels(trackId(myStack))){
kml_layer(as(myStack[[i]],'SpatialLines'))
}
## close the file
kml_close('myStack.kml')
## export KML using writeOGR ##
for(i in 1:nrow(myStack@idData)){
writeOGR(as(myStack[[i]], "SpatialPointsDataFrame"),
paste(row.names(myStack@idData)[i], ".kml", sep=""),
row.names(myStack@idData)[i], driver="KML")
writeOGR(as(myStack[[i]], "SpatialLinesDataFrame"),
paste(row.names(myStack@idData)[i], "-track.kml", sep=""),
row.names(myStack@idData)[i], driver="KML")
}
These binaries (installable software) and packages are in development.
They may not be fully stable and should be used with caution. We make no claims about them.