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
|
# Calculates the autocovariate to be used in autonormal, autopoisson
# or autologistic regression. Three distance-weighting schemes are available
#
# z is the response variable
# xy is the matrix of coordinates
# nbs is "neighbourhood size", selected by user; default is 1
# type defines the weighting scheme:
# "one" gives equal weight to all data points in the neighbourhood;
# "inverse" (the default) weights by inverse Euclidean distance;
# "inverse.squared" weights by the square of "inverse"
#
# by Carsten F. Dormann, 04.11.2005, carsten.dormann@ufz.de
# Re-implementation allowing list representation
# Roger Bivand 28.11.2005
# Upgrade to sp classes February 2007, longlat sanity June 2010
# style changed from "W" to "B" 2015-01-27; see Bardos et al. (2015) for details
autocov_dist <- function(z, xy, nbs=1, type="inverse", zero.policy=NULL,
style="B", longlat=NULL) {
if (is.null(zero.policy))
zero.policy <- get("zeroPolicy", envir = .spdepOptions)
stopifnot(is.logical(zero.policy))
stopifnot(is.vector(z))
if (type=="one") expo <- 0
if (type=="inverse") expo <- 1
if (type=="inverse.squared") expo <- 2
if (inherits(xy, "SpatialPoints")) {
if ((is.null(longlat) || !is.logical(longlat))
&& !is.na(is.projected(xy)) && !is.projected(xy)) {
longlat <- TRUE
} else longlat <- FALSE
xy <- coordinates(xy)
} else if (is.null(longlat) || !is.logical(longlat)) longlat <- FALSE
stopifnot(ncol(xy) == 2)
if (longlat) {
bb <- bbox(xy)
if (!.ll_sanity(bb))
warning("Coordinates are not geographical: longlat argument wrong")
}
nb <- dnearneigh(xy, 0, nbs, longlat=longlat)
if (any(card(nb) == 0)) warning(paste("With value", nbs,
"some points have no neighbours"))
nbd <- nbdists(nb, xy, longlat=longlat)
if (expo == 0) lw <- nb2listw(nb, style=style, zero.policy=zero.policy)
else {
gl <- lapply(nbd, function(x) 1/(x^expo))
lw <- nb2listw(nb, glist=gl, style=style, zero.policy=zero.policy)
}
lag(lw, z, zero.policy=zero.policy)
}
.ll_sanity <- function(bb) {
TOL <- get_ll_TOL()
tol <- .Machine$double.eps ^ TOL
W <- bb[1,1] < -180 &&
!isTRUE(all.equal((bb[1, 1] - -180), 0, tolerance = tol))
if (W) attr(W, "out") <- bb[1,1]
E <- bb[1,2] > 360 &&
!isTRUE(all.equal((bb[1, 2] - 360), 0, tolerance = tol))
if (E) attr(E, "out") <- bb[1,2]
S<- bb[2,1] < -90 &&
!isTRUE(all.equal((bb[2, 1] - -90), 0, tolerance = tol))
if (S) attr(S, "out") <- bb[2,1]
N <- bb[2,2] > 90 &&
!isTRUE(all.equal((bb[2, 2] - 90), 0, tolerance = tol))
if (N) attr(N, "out") <- bb[2,2]
res <- !(any(W || E || S || N))
attr(res, "details") <- list(W, E, S, N)
res
}
|
