Description Usage Arguments Author(s) Examples
Get the sum of squared residuals in for a geographically-weighted GLM
1 | gwglmnet.adaptive.ssr(bw, x, y, colocated, dist, s, verbose, family, prior.weights, gweight, type, target)
|
bw |
|
x |
|
y |
|
colocated |
|
dist |
|
s |
|
verbose |
|
family |
|
prior.weights |
|
gweight |
|
type |
|
target |
Wesley Brooks
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 | ##---- Should be DIRECTLY executable !! ----
##-- ==> Define data, use random,
##-- or do help(data=index) for the standard data sets.
## The function is currently defined as
function (bw, x, y, colocated, dist, s, verbose, family, prior.weights,
gweight, type, target)
{
reps = length(colocated)
loow = gweight(dist, bw)[-colocated]
w <- prior.weights[-colocated] * loow
xx = as.matrix(x[-colocated, ])
yy = as.matrix(y[-colocated])
m <- ncol(xx)
n <- nrow(xx)
one <- rep(1, n)
meanx <- drop(one %*% xx)/n
x.centered <- scale(xx, meanx, FALSE)
normx <- sqrt(drop(one %*% (x.centered^2)))
names(normx) <- NULL
xs = x.centered
for (k in 1:dim(x.centered)[2]) {
if (normx[k] != 0) {
xs[, k] = xs[, k]/normx[k]
}
else {
xs[, k] = rep(0, dim(xs)[1])
normx[k] = Inf
}
}
glm.step = try(glm(yy ~ xs, family = family, weights = w))
if (class(glm.step) == "try-error") {
cat(paste("Couldn't make a model for finding the SSR at bandwidth ",
bw, "\n", sep = ""))
return(Inf)
}
beta.glm = glm.step$coeff[2:(m + 1)]
adapt.weight = abs(beta.glm)
for (k in 1:dim(x.centered)[2]) {
if (!is.na(adapt.weight[k])) {
xs[, k] = xs[, k] * adapt.weight[k]
}
else {
xs[, k] = rep(0, dim(xs)[1])
adapt.weight[k] = 0
}
}
if (family == "binomial") {
model = glmnet(x = xs, y = cbind(1 - yy, yy), weights = w,
family = family, lambda = s)
}
else {
model = glmnet(x = xs, y = yy, weights = w, family = family,
lambda = s)
}
ll = model$lambda
xs.colocated = (x[colocated, ] - meanx) * adapt.weight/normx
predictions = predict(model, newx = matrix(xs.colocated,
nrow = reps, ncol = dim(xs)[2]), s = ll, type = "response",
)
cv.error = colSums(abs(matrix(predictions - matrix(y[colocated],
nrow = reps, ncol = length(ll)), nrow = reps, ncol = length(ll))))
s.optimal = ll[which.min(cv.error)]
fitted = predict(model, newx = xs, s = s.optimal, type = "response")
if (family == "poisson")
pearson.resid = sum(w * (yy - fitted)^2/fitted)
if (family == "binomial")
pearson.resid = sum(w * (yy - fitted)^2/(fitted * (1 -
fitted)))
(abs(pearson.resid - target))^2
}
|
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.