[½»Á÷] ¡¾ÇóÖú¡¿ÖðÒ»·¨½»²æÑéÖ¤Ïà¹ØϵÊý ÒÑÓÐ1È˲ÎÓë

3Â¥: Originally posted by snoopyzhao at 2010-11-05 20:09:14:
ÌùÒ»¸öÎÒÓà R ÓïÑÔдµÄ LOO-LGO CV µÄ³ÌÐò£¬Ó¦¸ÃÄܹ»¿´µÃ¶®£¬ÖÁÓÚÆäËüÓïÑÔÈçºÎʵÏÖ£¬ÎÒ¾ÍÎÞÄÜΪÁ¦ÁË¡­¡­

[code]
lmcv <- function(obj, ng = NULL, random = FALSE) {
#
# The Leave-One-Out (LOO) and/ ...

lmcv <- function(obj, ng = NULL, random = FALSE) {
#
# The Leave-One-Out (LOO) and/or Leave-Group-Out (LGO) Cross-Validation in R for (Multiple) Linear Regression.
#
# Input:
#   obj: the model of MLR
#   ng: number of group, if missing, do LOO
#   random: logical, if TRUE, do random CV
#
# Output:
#   q.squared: cross-validation relation coefficient.
#   SDEP: Standard Deviation of Error of Predictions
#   newsq: variance in Y explained only for LOO CV
#
# Usage:
#   loo <- cv.lm(obj)
#   lgo <- cv.lm(obj, ng = 5)
#
    data <- data.frame(y = obj$model[,1], x = obj$model[, -1])
    col.names <- colnames(data)
    N <- nrow(data)
    if (random == TRUE) data <- data[sample(1:N),]
    if (missing(ng)) ng <- N # LOO CV
    ytest <- numeric(N)
    ypred <- numeric(N)
    newrsq <- numeric(ng)
    g <- N %/% ng
    for (i in 1:ng) {
        if (g == 1) {
            index <- i
        }
        else {
            index <- c(i, ng * seq(1, (g - 1)) + i)
            if (N %% ng != 0 & i <= N %% ng) index <- c(index, (g * ng + i))
        }
        ytest[index] <- data[index,1]
        newtrain <- data[-index, ]
        xtest <- data.frame(x = data[index, -1])
        colnames(xtest) <- col.names[-1]
        newfm <- lm(y ~., data = newtrain)
        newrsq[i] <- summary(newfm)$r.squared
        ypred[index] <- predict(newfm, xtest)
    }
    q.squared <- 1 - sum((ytest - ypred)^2) / sum((ytest - mean(ytest))^2)
    SDEP <- sqrt(sum((ytest - ypred)^2) / N)
    return(list(q.squared = q.squared, SDEP = SDEP, newrsq = newrsq))
}