Nothing
"kernelbb" <- function(ltr, sig1, sig2, grid = 40,
same4all=FALSE, byburst=FALSE, extent = 0.5,
nalpha=25)
{
## verifications
if (!inherits(ltr, "ltraj"))
stop("tr should be of class \"ltraj\"")
x <- do.call("rbind",ltr)[,c("x","y")]
bu <- burst(ltr)
x$burst <- unlist(lapply(1:length(bu), function(i) {
rep(bu[i], nrow(ltr[[i]]))
}))
idd <- id(ltr)
x$id <- unlist(lapply(1:length(idd), function(i) {
rep(idd[i], nrow(ltr[[i]]))
}))
x$date <- unlist(lapply(ltr, function(y) y$date))
## Bases
sorties <- list()
gr <- grid
x <- x[!is.na(x$x),]
x <- x[!is.na(x$y),]
xy<-x[,c("x","y")]
sig12<-sig1^2
sig22<-sig2^2
h<-c(sig1, sig2)
names(h)<-c("sig1","sig2")
fac<-x$burst
if (!byburst)
fac<-x$id
fac<-factor(fac)
lixy<-split(x,fac)
if (is.list(grid)) {
if (is.null(names(grid)))
stop("when grid is a list, it should have named elements")
nn <- names(grid)
lev <- levels(fac)
if (length(lev) != length(nn))
stop("the length of the grid list should be equal to the number of levels of id/burst")
if (!all(lev%in%nn))
stop("some levels of id/bursts do not have corresponding grids")
}
if (same4all) {
if (inherits(grid, "SpatialPoints"))
stop("when same4all is TRUE, grid should be a number")
grid <- .makegridUD(xy, grid, extent)
}
for (i in 1:length(lixy)) {
dft<-lixy[[i]]
df<-dft[,c("x","y")]
if (!is.list(gr)) {
if (!inherits(gr, "SpatialPoints")) {
if (length(as.vector(gr)) == 1) {
if (!is.numeric(gr))
stop("non convenient grid")
if (!same4all) {
grid <- .makegridUD(xy, gr, extent)
}
}
}
} else {
grid <- gr[[names(lixy)[i]]]
}
gridded(grid) <- TRUE
fullgrid(grid) <- TRUE
pfs <- proj4string(grid)
grrw <- gridparameters(grid)
if (nrow(grrw) > 2)
stop("grid should be defined in two dimensions")
if (is.loaded("adehabitatMA")) {
opteps <- adehabitatMA::adeoptions()$epsilon
} else {
opteps <- 1e-08
}
if ((grrw[1, 2] - grrw[2, 2])> opteps)
stop("the cellsize should be the same in x and y directions")
xyg <- coordinates(grid)
xg <- unique(xyg[,1])
yg <- unique(xyg[,2])
date<-as.double(dft$date)-min(as.double(dft$date))
toto<-.C("kernelbbc", double(length(xg)*length(yg)), as.double(xg),
as.double(yg), as.integer(length(yg)),as.integer(length(xg)),
as.integer(nrow(df)), as.double(sig12), as.double (sig22),
as.double(df$x), as.double(df$y), as.double(date),
as.integer(nalpha),
PACKAGE="adehabitatHR")
too <- c(t(matrix(toto[[1]], nrow=length(xg), byrow=TRUE)))
UD <- data.frame(ud=too)
coordinates(UD) <- expand.grid(yg,xg)[,2:1]
gridded(UD) <- TRUE
UD <- new("estUD", UD)
slot(UD, "h") <- list(values=c(sig12=sig12, sig22=sig22),
meth="BB-specified")
slot(UD, "vol") <- FALSE
if (!is.na(pfs))
proj4string(UD) <- CRS(pfs)
sorties[[names(lixy)[i]]] <- UD
}
class(sorties) <- "estUDm"
if (length(sorties)==1)
sorties <- sorties[[1]]
return(sorties)
}
liker <- function (tr, rangesig1, sig2, le=1000,
byburst = FALSE, plotit=TRUE)
{
locchecktz <- function(x)
{
if (!inherits(x, "ltraj"))
stop("x should be of class \"ltraj\"")
tz <- lapply(x, function(y) attr(y$date, "tzone"))
atz <- all(sapply(1:length(tz), function(i) identical(tz[[i]],tz[[1]])))
if (!atz)
stop("multiple time zones not allowed in objects of class ltraj")
if (is.null(tz[[1]]))
return("")
if (tz[[1]]=="")
return("")
return(tz[[1]])
}
locltraj2traj <- function(x)
{
if (!inherits(x, "ltraj"))
stop("x should be of class \"ltraj\"")
if (!attr(x,"typeII"))
stop("x should be of type II (time recorded")
id <- factor(unlist(lapply(x, function(y)
id <- rep(attr(y,"id"), nrow(y)))))
burst <- factor(unlist(lapply(x, function(y)
id <- rep(attr(y,"burst"), nrow(y)))))
if (attr(x,"typeII"))
tz <- locchecktz(x)
if (!is.null(infolocs(x)))
infol <- do.call("rbind", infolocs(x))
res <- do.call("rbind", x)
res <- cbind(id,burst,res)
if (!is.null(infolocs(x)))
res <- cbind(res, infol)
class(res) <- c("traj","data.frame")
return(res)
}
x <- locltraj2traj(tr)
if (!inherits(x, "traj"))
stop("tr should be of class \"ltraj\"")
sorties <- list()
gr <- grid
x <- x[!is.na(x$x), ]
x <- x[!is.na(x$y), ]
xy <- x[, c("x", "y")]
sig22 <- sig2^2
fac <- x$burst
if (!byburst)
fac <- x$id
fac <- factor(fac)
lixy <- split(x, fac)
if (plotit)
par(mfrow=n2mfrow(length(lixy)))
so <- list()
for (i in 1:length(lixy)) {
dft <- lixy[[i]]
df <- dft[, c("x", "y")]
vsig <- seq(rangesig1[1], rangesig1[2],
length=le)
date <- as.double(dft$date) - min(as.double(dft$date))
huhu <- .C("CVL", as.double(t(as.matrix(df))),
as.double(date), as.integer(nrow(df)),
double(length(vsig)), as.double(vsig^2),
as.integer(length(vsig)),
as.double(sig22),
PACKAGE="adehabitatHR")
sig1=vsig[which.max(huhu[[4]])]
if (plotit)
plot(vsig, huhu[[4]], ty="l", main=paste(names(lixy)[i],
"\nMax for:",round(sig1,4)),
xlab="sig1", ylab="log-likelihood")
so[[i]] <- list(sig1=sig1, sig2=sig2,
cv=data.frame(sigma1=vsig, LogLikelihood=huhu[[4]]))
}
names(so) <- names(lixy)
class(so) <- "liker"
return(so)
}
print.liker <- function(x, ...)
{
if (!inherits(x, "liker"))
stop("x should be of class liker")
cat("*****************************************\n\n")
cat("Maximization of the log-likelihood for parameter\n")
cat("sig1 of brownian bridge\n\n")
lapply(1:length(x), function(i) {
cat(paste(names(x)[i],": Sig1 =", round(x[[i]]$sig1,4),
"Sig2 =", x[[i]]$sig2,"\n"))
})
}
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.