1
- # Copyright 2001-9, 2021 by Roger Bivand
1
+ # Copyright 2001-24 by Roger Bivand
2
2
#
3
3
4
4
@@ -7,6 +7,7 @@ geary <- function(x, listw, n, n1, S0, zero.policy=attr(listw, "zero.policy")) {
7
7
zero.policy <- get(" zeroPolicy" , envir = .spdepOptions )
8
8
stopifnot(is.logical(zero.policy ))
9
9
stopifnot(is.vector(x ))
10
+ stopifnot(all(is.finite(x )))
10
11
# z <- scale(x, scale=FALSE)
11
12
z <- scale(x )
12
13
zz <- sum(z ^ 2 )
@@ -34,16 +35,24 @@ geary.intern <- function(x, listw, n, zero.policy=attr(listw, "zero.policy"), ty
34
35
}
35
36
36
37
geary.test <- function (x , listw , randomisation = TRUE , zero.policy = attr(listw , " zero.policy" ),
37
- alternative = " greater" , spChk = NULL , adjust.n = TRUE ) {
38
+ alternative = " greater" , spChk = NULL , adjust.n = TRUE , na.action = na.fail ) {
38
39
if (is.null(zero.policy ))
39
40
zero.policy <- get(" zeroPolicy" , envir = .spdepOptions )
40
41
stopifnot(is.logical(zero.policy ))
41
42
alternative <- match.arg(alternative , c(" less" , " greater" , " two.sided" ))
42
- if (! inherits(listw , " listw" )) stop(paste(deparse(substitute(listw )),
43
- " is not a listw object" ))
44
- if (! is.numeric(x )) stop(paste(deparse(substitute(x )),
45
- " is not a numeric vector" ))
46
- if (any(is.na(x ))) stop(" NA in X" )
43
+ wname <- deparse(substitute(listw ))
44
+ if (! inherits(listw , " listw" )) stop(wname , " is not a listw object" )
45
+ xname <- deparse(substitute(x ))
46
+ if (! is.numeric(x )) stop(xname ,
47
+ " is not a numeric vector" )
48
+ if (deparse(substitute(na.action )) == " na.pass" )
49
+ stop(" na.pass not permitted" )
50
+ x <- na.action(x )
51
+ na.act <- attr(x , " na.action" )
52
+ if (! is.null(na.act )) {
53
+ subset <- ! (1 : length(listw $ neighbours ) %in% na.act )
54
+ listw <- subset(listw , subset , zero.policy = zero.policy )
55
+ }
47
56
n <- length(listw $ neighbours )
48
57
if (n != length(x )) stop(" objects of different length" )
49
58
if (is.null(spChk )) spChk <- get.spChkOption()
@@ -52,9 +61,9 @@ geary.test <- function(x, listw, randomisation=TRUE, zero.policy=attr(listw, "ze
52
61
53
62
wc <- spweights.constants(listw , zero.policy , adjust.n = adjust.n )
54
63
S02 <- wc $ S0 * wc $ S0
55
- res <- geary(x , listw , wc $ n , wc $ n1 , wc $ S0 , zero.policy )
64
+ res <- geary(as.vector( x ) , listw , wc $ n , wc $ n1 , wc $ S0 , zero.policy )
56
65
C <- res $ C
57
- if (is.na(C )) stop(" NAs generated in geary - check zero.policy" )
66
+ if (is.na(C )) stop(" NAs generated in geary - check zero.policy and na.action " )
58
67
K <- res $ K
59
68
EC <- 1
60
69
if (randomisation ) {
@@ -85,29 +94,43 @@ geary.test <- function(x, listw, randomisation=TRUE, zero.policy=attr(listw, "ze
85
94
names(vec ) <- c(" Geary C statistic" , " Expectation" , " Variance" )
86
95
method <- paste(" Geary C test under" , ifelse(randomisation ,
87
96
" randomisation" , " normality" ))
88
- data.name <- paste(deparse(substitute(x )), " \n weights:" ,
89
- deparse(substitute(listw )), " \n " )
97
+ data.name <- paste(xname , " \n weights:" , wname ,
98
+ ifelse(is.null(na.act ), " " , paste(" \n omitted:" ,
99
+ paste(na.act , collapse = " , " ))),
100
+ ifelse(adjust.n && isTRUE(any(sum(card(listw $ neighbours ) == 0L ))),
101
+ " \n n reduced by no-neighbour observations" , " " ), " \n " )
90
102
res <- list (statistic = statistic , p.value = PrC , estimate = vec ,
91
103
alternative = ifelse(alternative == " two.sided" , alternative ,
92
104
paste(" Expectation" , alternative , " than statistic" )),
93
105
method = method , data.name = data.name )
106
+ if (! is.null(na.act )) attr(res , " na.action" ) <- na.act
94
107
class(res ) <- " htest"
95
108
res
96
109
}
97
110
98
111
geary.mc <- function (x , listw , nsim , zero.policy = attr(listw , " zero.policy" ),
99
- alternative = " greater" , spChk = NULL , adjust.n = TRUE , return_boot = FALSE ) {
112
+ alternative = " greater" , spChk = NULL , adjust.n = TRUE , return_boot = FALSE ,
113
+ na.action = na.fail ) {
100
114
if (is.null(zero.policy ))
101
115
zero.policy <- get(" zeroPolicy" , envir = .spdepOptions )
102
116
stopifnot(is.logical(zero.policy ))
103
117
stopifnot(is.vector(x ))
104
118
alternative <- match.arg(alternative , c(" less" , " greater" , " two.sided" ))
105
- if ( ! inherits( listw , " listw " )) stop(paste( deparse(substitute(listw )),
106
- " is not a listw object" ) )
107
- if ( ! is.numeric( x )) stop(paste( deparse(substitute(x )),
108
- " is not a numeric vector" ) )
119
+ wname <- deparse(substitute(listw ))
120
+ if ( ! inherits( listw , " listw " )) stop( wname , " is not a listw object" )
121
+ xname <- deparse(substitute(x ))
122
+ if ( ! is.numeric( x )) stop( xname , " is not a numeric vector" )
109
123
if (missing(nsim )) stop(" nsim must be given" )
110
- if (any(is.na(x ))) stop(" NA in X" )
124
+ if (deparse(substitute(na.action )) == " na.pass" )
125
+ stop(" na.pass not permitted" )
126
+ x <- na.action(x )
127
+ na.act <- attr(x , " na.action" )
128
+ if (! is.null(na.act )) {
129
+ subset <- ! (1 : length(listw $ neighbours ) %in% na.act )
130
+ listw <- subset(listw , subset , zero.policy = zero.policy )
131
+ if (return_boot )
132
+ message(" NA observations omitted: " , paste(na.act , collapse = " , " ))
133
+ }
111
134
n <- length(listw $ neighbours )
112
135
if (n != length(x )) stop(" objects of different length" )
113
136
if (is.null(spChk )) spChk <- get.spChkOption()
@@ -134,7 +157,7 @@ geary.mc <- function(x, listw, nsim, zero.policy=attr(listw, "zero.policy"),
134
157
res <- numeric (length = nsim + 1 )
135
158
for (i in 1 : nsim ) res [i ] <- geary(sample(x ), listw , n , wc $ n1 , wc $ S0 ,
136
159
zero.policy )$ C
137
- res [nsim + 1 ] <- geary(x , listw , n , wc $ n1 , wc $ S0 , zero.policy )$ C
160
+ res [nsim + 1 ] <- geary(as.vector( x ) , listw , n , wc $ n1 , wc $ S0 , zero.policy )$ C
138
161
rankres <- rank(res )
139
162
xrank <- rankres [length(res )]
140
163
diff <- nsim - xrank
@@ -153,12 +176,13 @@ geary.mc <- function(x, listw, nsim, zero.policy=attr(listw, "zero.policy"),
153
176
parameter <- xrank
154
177
names(parameter ) <- " observed rank"
155
178
method <- " Monte-Carlo simulation of Geary C"
156
- data.name <- paste(deparse(substitute( x )) , " \n weights:" ,
157
- deparse(substitute( listw )) , " \n number of simulations + 1 :" ,
158
- nsim + 1 , " \n " )
179
+ data.name <- paste(xname , " \n weights:" , wname ,
180
+ ifelse(is.null( na.act ) , " " , paste( " \n omitted :" , paste( na.act ,
181
+ collapse = " , " ))), " \n number of simulations + 1: " , nsim + 1 , " \n " )
159
182
lres <- list (statistic = statistic , parameter = parameter ,
160
183
p.value = pval , alternative = alternative , method = method ,
161
184
data.name = data.name , res = res )
185
+ if (! is.null(na.act )) attr(res , " na.action" ) <- na.act
162
186
class(lres ) <- c(" htest" , " mc.sim" )
163
187
lres
164
188
}
0 commit comments