在这种情况下,为什么 L_BFGS_B 优化会跳到可行解决方案的极端范围?
Why does L_BFGS_B optimization skip to extreme range of viable solutions in this instance?
我很感激这是一个非常具体的问题!
为了帮助解释:我正在探索使用线性优化器来演示函数曲面中的尖锐 'cliffs' 会导致非最佳解决方案。 R 中的可重现代码如下:
library(glmnet)
library(mice)
# Load data
df <- read.csv(paste0('https://raw.githubusercontent.com/jbrownlee/Datasets',
'/master/pima-indians-diabetes.data.csv'), header = F)
colnames(df) <- c('Pregnancies', 'Glucose', 'BloodPressure', 'SkinThickness',
'Insulin', 'BMI', 'DiabetesPedigreeFunction', 'Age', 'Outcome')
set.seed(40)
# Impute 0 (missing) values for columns 2 through 8 (Glucose - Age)
df[2:8] <- lapply(df[2:8], function(x) replace(x, x %in% 0, NA))
micedf <- mice(df)
df <- complete(micedf)
# Create train/test split
sample_size <- floor(0.75 * nrow(df))
train_index <- sample(seq_len(nrow(df)), size = sample_size)
train <- df[train_index,]
test <- df[-train_index,]
# Generate model matrix format for glmnet
x <- as.matrix(train[,1:8])
y <- train$Outcome
# Fitting function
GLM_tune <- function(alpha) {
set.seed(40)
cvglmnet <- glmnet::cv.glmnet(x, y, nfolds = 5, family = "binomial",
alpha = alpha, type.measure = "auc",
parallel = F)
return (cvglmnet$cvm[cvglmnet$lambda == cvglmnet$lambda.1se]) }
现在,如果我输入一个介于 0 和 1 之间的值,如下所示:
optim(par = 0.9, fn = GLM_tune, lower = 0, upper = 1,
control = list(fnscale = -1, trace=3), method = c("L-BFGS-B"))
# >> $par = 0.86
优化器爬升到局部最大值 - 我已经通过使用以下方法探索整个表面区域进行了测试:
surf <- data.frame(alpha = 0, auc = 0)
for (a in seq(from=0, to=1000)) {
surf[a+1,1] <- a/1000
surf[a+1,2] <- GLM_tune(a/1000)
}
library(ggplot2)
ggplot() +
geom_point(data=surf, size = 1.2, color = "black", aes(alpha, auc))
但是,当我将 alpha = 1 设置为起点时,算法在第二次迭代时转换为 alpha = 0,然后作为 'final' 解决方案退出:
optim(par = 1, fn = GLM_tune, lower = 0, upper = 1,
control = list(fnscale = -1, trace=3), method = c("L-BFGS-B"))
# >> $par = 0
为什么会这样?显然我不完全理解该算法,但我假设 optim 函数中的默认步进为 0.001(请参阅 ndeps)——那么为什么下一步会走向相反的极端?我是否缺少应针对这些问题设置的关键参数?
从 objective 函数的路径可以清楚地看出它有很多局部最大值,因此,像 "L-BFGS-B" 这样的基于梯度的优化算法 不是 适合寻找全局最大值。
此外,我的R(3.6),
optim(par = 1, fn = GLM_tune, lower = 0, upper = 1,
control = list(fnscale = -1, trace=3), method = c("L-BFGS-B"))$par
## [1] 1
returns 1
而不是您指出的 0
。
要理解为什么收敛到1
,我们可以看一下"L-BFGS-B"算法的优化路径。为此,我更喜欢使用 R 包 optimParallel https://CRAN.R-project.org/package=optimParallel。我是包的作者:
library("optimParallel")
cl <- makeCluster(2); setDefaultCluster(cl=cl)
clusterExport(cl, c("x", "y")) # export implicitly used values
optimParallel(par = 1, fn = GLM_tune,
lower = 0, upper = 1,
control = list(fnscale = -1),
parallel = list(optimParallel.loginfo=TRUE))$loginfo
## step par1 fn gr1
## 1.0000000 1.0000000 -0.8215854 0.0000000
我们看到 1
处的梯度是 0
。因此,算法停在 1
.
也就不足为奇了
我们可以用
检查近似梯度的计算
ndeps <- 0.001 # the default value
(GLM_tune(1) - GLM_tune(1-ndeps))/ndeps
## [1] 0
请注意,如果 1
不是上限,则 optim()
使用中心差分梯度近似。像
(GLM_tune(1+ndeps) - GLM_tune(1-ndeps))/(2*ndeps)
我很感激这是一个非常具体的问题!
为了帮助解释:我正在探索使用线性优化器来演示函数曲面中的尖锐 'cliffs' 会导致非最佳解决方案。 R 中的可重现代码如下:
library(glmnet)
library(mice)
# Load data
df <- read.csv(paste0('https://raw.githubusercontent.com/jbrownlee/Datasets',
'/master/pima-indians-diabetes.data.csv'), header = F)
colnames(df) <- c('Pregnancies', 'Glucose', 'BloodPressure', 'SkinThickness',
'Insulin', 'BMI', 'DiabetesPedigreeFunction', 'Age', 'Outcome')
set.seed(40)
# Impute 0 (missing) values for columns 2 through 8 (Glucose - Age)
df[2:8] <- lapply(df[2:8], function(x) replace(x, x %in% 0, NA))
micedf <- mice(df)
df <- complete(micedf)
# Create train/test split
sample_size <- floor(0.75 * nrow(df))
train_index <- sample(seq_len(nrow(df)), size = sample_size)
train <- df[train_index,]
test <- df[-train_index,]
# Generate model matrix format for glmnet
x <- as.matrix(train[,1:8])
y <- train$Outcome
# Fitting function
GLM_tune <- function(alpha) {
set.seed(40)
cvglmnet <- glmnet::cv.glmnet(x, y, nfolds = 5, family = "binomial",
alpha = alpha, type.measure = "auc",
parallel = F)
return (cvglmnet$cvm[cvglmnet$lambda == cvglmnet$lambda.1se]) }
现在,如果我输入一个介于 0 和 1 之间的值,如下所示:
optim(par = 0.9, fn = GLM_tune, lower = 0, upper = 1,
control = list(fnscale = -1, trace=3), method = c("L-BFGS-B"))
# >> $par = 0.86
优化器爬升到局部最大值 - 我已经通过使用以下方法探索整个表面区域进行了测试:
surf <- data.frame(alpha = 0, auc = 0)
for (a in seq(from=0, to=1000)) {
surf[a+1,1] <- a/1000
surf[a+1,2] <- GLM_tune(a/1000)
}
library(ggplot2)
ggplot() +
geom_point(data=surf, size = 1.2, color = "black", aes(alpha, auc))
但是,当我将 alpha = 1 设置为起点时,算法在第二次迭代时转换为 alpha = 0,然后作为 'final' 解决方案退出:
optim(par = 1, fn = GLM_tune, lower = 0, upper = 1,
control = list(fnscale = -1, trace=3), method = c("L-BFGS-B"))
# >> $par = 0
为什么会这样?显然我不完全理解该算法,但我假设 optim 函数中的默认步进为 0.001(请参阅 ndeps)——那么为什么下一步会走向相反的极端?我是否缺少应针对这些问题设置的关键参数?
从 objective 函数的路径可以清楚地看出它有很多局部最大值,因此,像 "L-BFGS-B" 这样的基于梯度的优化算法 不是 适合寻找全局最大值。
此外,我的R(3.6),
optim(par = 1, fn = GLM_tune, lower = 0, upper = 1,
control = list(fnscale = -1, trace=3), method = c("L-BFGS-B"))$par
## [1] 1
returns 1
而不是您指出的 0
。
要理解为什么收敛到1
,我们可以看一下"L-BFGS-B"算法的优化路径。为此,我更喜欢使用 R 包 optimParallel https://CRAN.R-project.org/package=optimParallel。我是包的作者:
library("optimParallel")
cl <- makeCluster(2); setDefaultCluster(cl=cl)
clusterExport(cl, c("x", "y")) # export implicitly used values
optimParallel(par = 1, fn = GLM_tune,
lower = 0, upper = 1,
control = list(fnscale = -1),
parallel = list(optimParallel.loginfo=TRUE))$loginfo
## step par1 fn gr1
## 1.0000000 1.0000000 -0.8215854 0.0000000
我们看到 1
处的梯度是 0
。因此,算法停在 1
.
我们可以用
检查近似梯度的计算ndeps <- 0.001 # the default value
(GLM_tune(1) - GLM_tune(1-ndeps))/ndeps
## [1] 0
请注意,如果 1
不是上限,则 optim()
使用中心差分梯度近似。像
(GLM_tune(1+ndeps) - GLM_tune(1-ndeps))/(2*ndeps)