Kernel Smoothing in Space and Time of the Animals' Use of Space

Share:

Description

These functions estimate the utilization distribution (UD) in space and time of animals monitored using radio-telemetry, using the product kernel estimator advocated by Keating and Cherry (2009).

Note that this approach has also been useful for the analysis of recoveries in programs involving ringed birds (Calenge et al. in prep., see section examples below).

kernelkc estimate the UD of several animals from an object of class ltraj.
kernelkcbase estimate one UD from a data frame with three columns indicating the spatial coordinates and associated timing.
exwc allows to search for the best value of the time smoothing parameter in the case where the time is considered as a circular variable (see details).
getvolumeUDk and getvolumeUDs provide utilities for home-range estimation.
getverticeshrk and getverticeshrs store the home range contours.

Usage

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
kernelkc(tr, h, tcalc, t0, grid = 40, circular = FALSE,
         cycle = 24 * 3600, same4all = FALSE,
         byburst = FALSE, extent = 0.5)

kernelkcbase(xyt, h, tcalc, t0, grid=40, circular=FALSE,
             cycle=24*3600, extent=0.5)

getvolumeUDk(kc)
getvolumeUDs(asc)

getverticeshrk(kc, lev = 95)
getverticeshrs(ascv, lev = 95)

exwc(hv)

Arguments

tr

an object of class ltraj

xyt

a data frame with three columns indicating the x and y coordinates, as well as the timing of the relocations.

kc

an object of class kernelkco returned by the function kernelkc

asc

an object of class asc returned by the function kernelkcbase

ascv

an object of class asc returned by the function getvolumeUDs

lev

the percentage level for home range contour estimation.

h

a numeric vector with three elements indicating the value of the smoothing parameters: the first and second elements are the smoothing parameters of the X and Y coordinates respectively, the third element is the smoothing parameter for the time dimension. If circular=TRUE it should be a smoothing parameter in the interval 0-1 (see details). If circular=FALSE this smoothing parameter should be given in seconds.

tcalc

the time at which the UD is to be estimated

t0

if circular=TRUE, this parameter indicates the time at which the time cycle begins (see examples).

grid

a number giving the size of the grid on which the UD should be estimated. Alternatively, this parameter may be an object of class asc, or a list of objects of class asc, with named elements corresponding to each level of the burst/id

circular

logical. Indicates whether the time should be considered as a circular variable (e.g., the 31th december 2007 is considered to be one day before the 1st january 2007) or not (e.g., the 31th december 2007 is considered to be one year after the 1st january 2007).

cycle

if circular=TRUE, the duration of the time cycle. for kernelkc, it should be given in seconds, and for kernelkcbase, in the units of the data (the units of the third column of xyt).

same4all

logical. If TRUE, the same grid is used for all levels of id/burst. If FALSE, one grid per id/burst is used.

byburst

logical. Indicates whether one UD should be estimated by burst of tr, or whether the data should be pooled across all bursts of each value of id in tr

extent

a value indicating the extent of the grid used for the estimation (the extent of the grid on the abscissa is equal to (min(xy[,1]) + extent * diff(range(xy[,1])))).

hv

a value of smoothing parameter for the time dimension.

...

additional arguments to be passed to the function contour.

Details

Keating and Cherry (in press) advocated the estimation of the UD in time and space using the product kernel estimator. These functions implement exactly this methodology.

For the spatial coordinates, the implemented kernel function is the biweight kernel.

Two possible approaches are possible to manage the time in the estimation process: (i) the time may be considered as a linear variable (e.g., the 31th december 2007 is considered to be one day before the 1st january 2007), or (ii) the time may be considered as a circular variable (e.g., the 31th december 2007 is considered to be one year after the 1st january 2007).

If the time is considered as a linear variable, the kernel function used in the estimation process is the biweight kernel. If the time is considered as a circular variable, the implemented kernel is the wrapped Cauchy distribution (as in the article of Keating and Cherry). In this latter case, the smoothing parameter should be chosen in the interval 0-1, with a value of 1 corresponding to a stronger smoothing.

These functions can only be used on objects of class "ltraj", but the estimation of the UD in time and space is also possible with other types of data (see the help page of kernelkcbase). Note that both kernelkc and kernelkcbase return conditional probability density function (pdf), i.e. the pdf to relocate an animal at a place, given that it has been relocated at time tcalc (i.e. the volume under the UD estimated at time tcalc is equal to 1 whatever tcalc).

The functions getvolumeUDk and getvolumeUDs modifies the UDs, so that the contour of the UD displayed by the function contour corresponds to the different percentage levels of home-range estimation (see examples).

The function exwc draws a graph of the wrapped Cauchy distribution for the chosen h parameter (for circular time), so that it is possible to make one's mind concerning the weight that can be given to the neighbouring points of a given time point. Note that although Keating and Cherry (in press) advocated the use of an automatic algorithm to select "optimal" values for the smoothing parameter, it is not implemented in adehabitat. Indeed, different smoothing parameters may allow to identify patterns at different scales, and we encourage the user to try several values before subjectively choosing the value which allows to more clearly identify the patterns of the UD.

Value

kernelkc returns a list of class "kernelkco" containing objects of class asc, mapping one estimate of the UD per burst or id (depending on the value of the parameter byburst).

getvolumekc returns a list of class "kernelkcov" containing objects of class asc, mapping the volume under the UD for each burst or id.

kernelkcbase returns an object of class "asc" mapping the estimated UD. The volume under the UD can be computed using the function getvolumeUDs.

getverticeshrk stores the home range contour as objects of class area in a list of class kver, with one component per burst/id.

getverticeshrs stores the home range contour as an object of class area.

Author(s)

Clement Calenge clement.calenge@oncfs.gouv.fr

References

Keating, K. and Cherry, S. (in press) Modeling utilization distributions in space and time. Ecology, in press.

Calenge, C., Guillemain, M., Gauthier-Clerc, M. and Simon, G. in prep. A new exploratory approach to the study of the spatio-temporal distribution of ring recoveries - the example of Teal (Anas crecca) ringed in Camargue, Southern France.

See Also

as.ltraj for additional information on objects of class ltraj, kernelUD for the "classical" kernel home range estimates.

Examples

  1
  2
  3
  4
  5
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
## Not run: 


################################################
##
## Illustrates the analysis of recoveries of
## ringed data

data(teal)
head(teal)

## compute the sequence of dates at which the
## probability density function (pdf) of recoveries is to be estimated

vv <- seq(min(teal$date), max(teal$date), length=50)
head(vv)

## The package "maps" should be installed for the example below
library(maps)


re <- lapply(1:length(vv), function(i) {

             ## Estimate the pdf. We choose a smoothing parameter of
             ## 2 degrees of lat-long for X and Y coordinates,
             ## and of 2 months for the time
             uu <- kernelkcbase(teal, c(2.5,2.5,2*30*24*3600), tcalc =
                                vv[i], grid=100, extent=0.1)

             ## now, we show the result
             ## potentially, we could type
             ##
             ## jpeg(paste("prdefu", i, ".jpg", sep=""))
             ##
             ## to store the figures in a file, and then to build a
             ## movie with the resulting files:
             ##

             image(uu, main=vv[i])

             ## highlight the area on which there is a probability
             ## equal to 0.95 to recover a bird
             hh <- getvolumeUDs(uu)

             contour(hh, levels=95, col="red",
                     drawlabels=FALSE, add=TRUE, lwd=2)

             ## The map:
             map(xlim=c(-20,70), ylim=c(30,80), add=TRUE)

             ## and if we had typed jpeg(...) before, we have to type
             ## dev.off()
             ## to close the device. When we have finished this loop
             ## We could combine the resulting files with imagemagick
             ## (windows) or mencoder (linux)
             })




################################################
##
## Illustrates how to explore the UD in time and
## space with the bear data

data(bear)

## compute the sequence of dates at which the UD is to be
## estimated
vv <- seq(min(bear[[1]]$date), max(bear[[1]]$date), length=50)
head(vv)

## estimates the UD at each time point
re <- lapply(1:length(vv), function(i) {

             ## estimate the UD. We choose a smoothing parameter of
             ## 1000 meters for X and Y coordinates, and of 72 hours
             ## for the time (after a visual exploration)
             uu <- kernelkc(bear, h = c(1000,1000,72*3600),
                            tcalc= vv[i])

             ## now, we show the result
             ## potentially, we could type
             ##
             ## jpeg(paste("UD", i, ".jpg", sep=""))
             ##
             ## to store the figures in a file, and then to build a
             ## movie with the resulting files:
             ##
             image(uu[[1]], main=vv[i])

             ## highlight the 95 percent home range
             hh <- getvolumeUDk(uu)

             contour(hh[[1]], levels=95, col="red",
                     drawlabels=FALSE, add=TRUE)

             ## and if we had typed jpeg(...) before, we have to type
             ## dev.off()
             ## to close the device. When we have finished this loop
             ## We could combine the resulting files with imagemagick
             ## (windows) or mencoder (linux)
             })


## Or, just show the home range:
re <- lapply(1:length(vv), function(i) {

             uu <- kernelkc(bear, h = c(1000,1000,72*3600),
                            tcalc= vv[i])

             pc <- getverticeshrk(uu, lev=95)
             plot(pc, xlim=c(510000, 530000),
                  ylim=c(6810000, 6825000), main=vv[i])
             })








##################################################
##
## Example with several wild boars (linear time)

## load wild boar data
data(puechcirc)


## keep only the first two circuits:
puech <- puechcirc[1:2]


## Now load te map of the elevation
data(puechabon)
elevation <- getkasc(puechabon$kasc, 1)


## compute the time point at which the UD is to be estimated
vv <- seq(min(puechcirc[[2]]$date), max(puechcirc[[2]]$date),
          length=50)


## The estimate the UD
re <- lapply(1:length(vv),
             function(i) {

                 ## We choose a smoothing parameter of 300 meters for
                 ## the x and y coordinates and of one hour for the time
                 ## (but try to play with these smoothing parameters)

                 uu <- kernelkc(puechcirc, h=c(300,300,3600),
                                tcalc = vv[i], same4all=TRUE,
                                extent=0.1)

                 ## show the elevation
                 image(elevation, main=vv[i],
                       xlim=c(698000,704000),
                       ylim=c(3156000,3160000))

                 ## and the UD, with contour lines
                 colo <- c("red","blue")
                 lapply(1:length(uu), function(i) {
                    contour(uu[[i]], col=colo[i], add=TRUE)
                 })

                 ## the blue contour lines show the UD of the mother and
                 ## the red ones correspond to her son. Adult wild boars
                 ## are known to be more "shy" that the youger ones.
                 ## Here, the low elevation corresponds to crop area
                 ## (vineyards). The young boar is the first and the
                 ## last in the crops
             })








##################################################
##
## Example with the bear, to illustrate (circular time)

data(bear)

## We consider a time cycle of 24 hours.
## the following vector contains the time points on the
## time circle at which the UD is to be estimated (note that
## the time is given in seconds)
vv <- seq(0, 24*3600-1, length=40)


## for each time point:
re <- lapply(1:length(vv),
             function(i) {

                 ## Estimation of the UD for the bear. We choose
                 ## a smoothing parameter of 1000 meters for the spatial
                 ## coordinates and a smoothing parameter equal to 0.2
                 ## for the time. We set the beginning of the time
                 ## cycle at midnight (no particular reason, just to
                 ## illustrate the function). So we pass, as t0, any
                 ## object of class POSIXct corresponding t a date at
                 ## this hour, for example the 12/25/2012 at 00H00
                 t0 <- as.POSIXct("2012-12-25 00:00")
                 uu <- kernelkc(bear, h=c(1000,1000,0.2), cycle=24*3600,
                                 tcalc=vv[i], t0=t0, circular=TRUE)

                 ## shows the results
                 ## first compute the hour for the title
                 hour <- paste(floor(vv[i]/3600), "hours",
                               floor((vv[i]%%3600)/60), "min")

                 ## compute the 95% home range
                 pc <- getverticeshrk(uu, lev=95)
                 plot(pc, xlim=c(510000, 530000),
                      ylim=c(6810000, 6825000), main=hour)

                 ## compute the 50% home range
                 pc <- getverticeshrk(uu, lev=50)
                 plot(pc, add=TRUE, colpol="blue")
                 
             })
 ## Now, each home range computed at a given time point corresponds to
 ## the area used by the animal at this time period.  We may for example
 ## try to identify the main difference in habitat composition of the
 ## home-range between different time, to identify differences in
 ## habitat use between different time of the day. We do not do it here
 ## (lack of example data)










##################################################
##
## Example of the use of the function kernelkcbase and
## related functions


## load the data
data(puechabon)
locs <- puechabon$locs

## keeps only the wild boar Chou
locs <- locs[locs$Name=="Jean",]

## compute the number of days since the beginning
## of the monitoring
dd <- cumsum(c(0, diff(strptime(locs$Date, "%y%m%d"))))
dd

## compute xyt. Note that t is here the number of
## days since the beginning of the monitoring (it
## is not an object of class POSIXt, but it may be)
xyt <- data.frame(locs[,c("X", "Y")], dd)

## Now compute the time points at which the UD is to be estimated:
vv <- 1:61

## and finally, show the UD changed with time:
re <- lapply(1:length(vv),
             function(i) {
                 ud <- kernelkcbase(xyt, h=c(300,300,20), 
                                    tcalc=vv[i], grid=100)
                 image(ud, main=vv[i])
                 contour(getvolumeUDs(ud), level=95,
                         col="red", lwd=2, add=TRUE)

                 ## Just to slow down the process
                 Sys.sleep(0.2)
                 })




## End(Not run)

Want to suggest features or report bugs for rdrr.io? Use the GitHub issue tracker.