intro_to_ml/15_loocv.R

40 lines
1.1 KiB
R

#loocv
ridge.svd <- function(data, degre) {
X <- data$X
n <- nrow(X)
A <- scale(outer(c(X), 1:degre, "^"))
Ym <- mean(data$Y)
Y <- data$Y - Ym
As <- svd(A)
d <- As$d
function(lambda) {
coef <- c(As$v %*% ((d / (d^2 + lambda)) * (t(As$u) %*% Y)))
coef <- coef / attr(A,"scaled:scale")
inter <- Ym - coef %*% attr(A,"scaled:center")
coef <- c(inter, coef)
trace.H <- sum(d^2 / (d^2 + lambda))
pred <- polyeval(coef, X)
gcv <- sum( ((Y - pred) / (1 - (trace.H / n))) ^ 2 ) / n
list(coef = coef, gcv = gcv)
}
}
ridge.loocv <- function(data, deg, lambdas) {
errs <- double(length(lambdas))
coefs <- matrix(data = NA, nrow = nrow(data$X), ncol = 1+deg)
ridge <- ridge.svd(data, deg)
idx <- 1
for(lambda in lambdas) {
res <- ridge(lambda)
coefs[idx,] <- res$coef
errs[idx] <- res$gcv
idx <- idx + 1
}
err.min <- min(errs)
lambda.best <- lambdas[which(errs == err.min)]
coef.best <- coefs[which(errs == err.min),]
pred <- polyeval(coef.best, data$X)
mae <- mean(abs(pred - data$Y))
list(coef = coef.best, lambda = lambda.best, mae = mae)
}