实现 GCV for Kernel Ridge 回归

问：

我正在用 R 语言实现 Kernel Ridge Regression。具体来说，我想使用广义交叉验证 （GCV） 确定 lambda 的最佳值，但我遇到了一些问题。

我编写了以下代码：

# gaussian kernel ---------------------------------------------------------

gKernel = function(xi, xj, nu=1){
  temp = exp(-nu * sum((xi - xj)^2))
  return(temp)
}




# generate data -------------------------------------------------------------------

generateDataSet = function(size, 
                           xlim = c(0, 1), 
                           scale = NULL) {
  x = runif(size, 
            min = xlim[1], 
            max = xlim[2])
  y = -x^2
  if (!is.null(scale)) {
    noise = rnorm(size, 
                  mean = 0, 
                  sd = scale)
    y = y + noise
  }
  df = data.frame(x = x, y = y)
  return(df)
}


data = generateDataSet(200, 
                       xlim = c(-2, 2), 
                       scale = 0.25)

# Plot the data
plot(data$x, data$y)


# Number of data points
N = nrow(data)

# Calculate the Gram matrix
K = matrix(0, nrow = N, ncol = N)
for (i in 1:N) {
  for (j in i:N) {
    K[i, j] = gKernel(data$x[i], data$x[j], 1)
    K[j, i] = K[i, j]
  }
}



kernelPredict = function(X, x, alpha, beta) {
  Y = 0
  for (i in 1:length(X)) {
    Y = Y + alpha[i] * gKernel(X[i], x, beta)
  }
  return(Y)
}


# implement GCV -----------------------------------------------------------

lams = c(1e-13,1e-12, 1e-11, 1e-10, 1e-9, 1e-8, 1e-7, 1e-6,
         1e-5, 1e-4, 1e-3, 1e-2, 1e-1, 0.25, 0.5, 0.75, 1)


gcvs = c()
for(j in 1:length(lams)){
  
  lam = lams[j]
  
  matH = solve(K + lam * diag(N))
  alpha_kr = matH %*% data$y
  
  Y_predict_kr = rep(0, length(data$x))
  for (i in 1:length(data$x)) {
    Y_predict_kr[i] = kernelPredict(data$x, 
                                    data$x[i], 
                                    alpha_kr,
                                    1)
  }
  
  trH = sum(diag(matH))
  
  gcv.den = (1 - trH / N)^2
  gcv.num = mean((data$y - Y_predict_kr)^2)
  gcv = gcv.num / gcv.den
  
  gcvs[j] = gcv
}

plot(lams, gcvs)

lams[which.min(gcvs)]

该代码将最佳 lambda 计算为 1e-13，但当我查看随附的图时，它似乎不是最佳选择。例如，当我将 lambda 设置为 1e-2 时，它似乎更适合数据（请参阅下图）。

在此处输入图像描述

在我的代码中造成这种现象的原因可能是什么？ 如果您能提供一些见解，我将不胜感激。