swh:1:snp:5e11768f61bcb1191e6cfb268a5997ff31b1baee
Tip revision: 92269b50ae23297ededa6a58d77823e1302a6b8e authored by Roger Koenker on 06 June 2021, 16:10:02 UTC
version 5.86
version 5.86
Tip revision: 92269b5
MCV.R
MCV <- function(lambdas, formula, data, tau = 0.5, k = 10){
F <- Munge(formula, lambdas = lambdas)
f <- rqss(F, data, tau = tau)
n <- f$n
m <- length(f$qss)
y <- f$y[1:n]
folds = sample(rep(1:k, length = n))
U = NULL
for(i in 1:k){
s = which(folds != i)
M = rqss(F, data = data[s,], tau = tau)
nd = data[-s,]
G = matrix(0,nrow(nd),m)
for(j in 1:m){ #remove extrapolates, if any
g = f$qss[[j]]$xyz[,1]
G[,j] = (min(g[s]) < g[-s]) & (g[-s] < max(g[s]))
}
h = as.logical(apply(G,1,prod))
u = predict(M, newdata = nd[h,]) - (y[-s])[h]
U = c(U,(u * (tau - (u < 0))))
}
mean(U)
}
set.seed(1729)
n <- 200
x <- sort(runif(n, 0, 20))
g0 <- function(x, tau)
log(x) + 0.2*(log(x))^3 + log(x) * qnorm(tau)/4
y <- g0(x, runif(n))
D <- data.frame(y = y, x = x)
lams <- mcvs <- seq(.02, 5, by = 0.2)
for(i in 1:length(mcvs))
mcvs[i] <- MCV(lams[i], y ~ qss(x, lambda = lambdas[1]), D)
par(mfrow = c(1,2))
plot(lams, mcvs, cex = .5, lwd = 2, type = 'l',
xlab = expression(lambda), ylab = expression(MCV( lambda )))
lambdastar <- lams[which.min(mcvs)]
plot(x, y, cex = .5, col = "grey")
f <- rqss(y ~ qss(x, lambda = lambdastar), data = D)
plot(f, add = TRUE, lwd = 2)
lines(x,g0(x, 0.5), col = "red", lwd = 2)
text(10, 1,bquote(lambda == ~ .(lambdastar)))