[R-sig-eco] simm.levy in adehabitat
Clément Calenge
clement.calenge at oncfs.gouv.fr
Thu Oct 2 10:57:04 CEST 2008
Hi Tim,
>> p.sett
>>
> *********** List of class ltraj ***********
> Type of the traject: Type II (time recorded)
> Regular traject. Time lag between two locs: 900 seconds
> Characteristics of the bursts:
> id burst nb.reloc NAs date.begin date.end
> 1 p1801 1988038937 6 1 2006-08-09 09:00:00 2006-08-09 10:15:00
> 2 p1801 1988038938 3 0 2006-08-10 02:00:00 2006-08-10 02:30:00
> 3 p1801 1988038939 9 1 2006-08-11 03:00:00 2006-08-11 05:00:00
> 4 p1801 1988038940 17 5 2006-08-12 02:00:00 2006-08-12 06:00:00
> 5 p1801 1988038941 23 10 2006-08-13 02:00:00 2006-08-13 07:30:00
> ...
> 237 p1801 1988039304 126 110 2007-08-10 01:15:00 2007-08-11 08:30:00
> 238 p1801 1988039305 1 0 2007-08-11 02:45:00 2007-08-11 02:45:00
>
> Ideally what I'd like to do is define the simulation using the
> characteristics of a given 'burst' from my ltraj object. I'd appreciate some
> pointers on the syntax of calling the simulation using the the
> characteristics of each row in p.sett. Sample code follows:
>
>
>> sim.levy<-simm.levy(date=1:126, id=p.sett$id=="1801",
>>
> burst=p.sett$burst=="1988039304", x0=c(174.48255,-36.819058), typeII=T)
> Error in as.ltraj(data.frame(co, si), date, id, burst, typeII = typeII) :
> id should be of the same length as xy, or of length 1
>
> How do I define the arguments of the function using the elements of a given
> burst?
>
First note that the "burst" and "id" are stored as attributes of each
element of the list of class ltraj, and not as elements of the list
themselves. It means that p.sett$burst and p.sett$id do not exist. To
get the burst or the id of the trajectories of the list, you can simply type
burst(p.sett)
id(p.sett)
See the help page of as.ltraj() for trajectory management with
adehabitat. Now, it seems that you want to simulate a levy walk at the
same times as your actual trajectory, beginning at the same point. In
this case, you should use the following code:
res <- lapply(1:length(p.sett), function(i) {
simm.levy(date = p.sett[[i]]$date, id = id(p.sett)[i],
x0 = c(p.sett[[i]]$x[1], p.sett[[i]]$y[1]),
burst = burst(p.sett)[i], typeII = T, mu=2)
})
result <- do.call("c", res)
plot(result)
But note that you may have to change the parameter mu, depending on what
you want to do. See the help page of simm.levy.
BTW, there is a bug in the function concerning the parameter x0. If you
want to use parameter x0 (required only if you want to match the
trajectory with maps), you should use the following fixed function,
which will be included in the next version of adehabitat:
simm.levy <- function (date = 1:500, mu = 2, l0 = 1, x0 = c(0, 0), id =
"A1",
burst = id, typeII = TRUE)
{
if (typeII)
class(date) <- c("POSIX", "POSIXct")
n <- length(date)
dt <- c(diff(unclass(date)))
if (all(dt - dt[1] > 1e-07))
stop("the time lag between relocations should be constant")
ang <- runif(n - 2, -pi, pi)
v = dt * (l0 * (runif(n - 1)^(1/(1 - mu))))
ang = cumsum(c(runif(1, 0, 2 * pi), ang))
si = c(x0[2], x0[2]+cumsum(v * sin(ang)))
co = c(x0[1], x0[1]+cumsum(v * cos(ang)))
res <- as.ltraj(data.frame(co, si), date, id, burst, typeII = typeII)
return(res)
}
Hope this helps,
Clément.
--
Clément CALENGE
Office national de la chasse et de la faune sauvage
Saint Benoist - 78610 Auffargis
tel. (33) 01.30.46.54.14
More information about the R-sig-ecology
mailing list