问题描述
我正在使用 CVXR
建模包来解决凸优化问题。我确信问题是凸的并且它遵循 DCP 规则,但是如果我使用 CVXR
检查 DCP 规则,它返回 False
。但是,如果我处理完全相同的问题并使用 CVXPY
检查它,它会返回 True
(如预期)
这里发生了什么?我附上了 R 和 Python 中这种行为的最小可重现示例:
使用 CVXR
的 R 代码
library(splines2)
library(CVXR)
deriv_basis = splines2::dbs(seq(0,1,length.out=100),degree=3,intercept=T,df=30,derivs=2)
R = t(deriv_basis) %*% deriv_basis
beta_var = CVXR::::Variable(nrow(R))
q = CVXR::quad_form(beta_var,R)
CVXR::is_dcp(q)
[1] FALSE
write.table(x=R,file='R.csv'),row.names=F,sep=';')
使用 CVXPY
的 Python 代码
import cvxpy
import pandas as pd
R = pd.read_csv('R.csv',sep=';').values
beta_var = cvxpy.Variable(R.shape[1])
q = cvxpy.quad_form(beta_var,R)
q.is_dcp()
Out[1]: True
有人可以解释这里发生了什么以及如何解决它以便我可以使用 CVXR?
解决方法
问题是 R 矩阵中的负特征值。如果您通过将其设置为零来解决该问题,例如,它满足 dcp 条件。我还修复了问题代码中的语法错误并删除了多余的 :: 。另一种可能性(未显示)是在 pracma 包中使用 nearest_spd
来调整 R 矩阵。
library(splines2)
library(CVXR)
deriv_basis <- dbs(seq(0,1,length.out=100),degree = 3,intercept = TRUE,df = 30,derivs = 2)
R <- t(deriv_basis) %*% deriv_basis
e <- eigen(R)
# check decomposition
all.equal(R,e$vectors %*% diag(e$values) %*% t(e$vectors),check.attributes = FALSE)
## [1] TRUE
e$values # note negative value
## [1] 1.095213e+08 1.095213e+08 1.056490e+07 1.055430e+07 1.052481e+07
## [6] 1.046063e+07 1.034247e+07 1.015017e+07 9.866358e+06 9.485145e+06
## [11] 8.643220e+06 8.280963e+06 7.549803e+06 6.731472e+06 5.853402e+06
## [16] 4.949804e+06 4.056714e+06 3.209045e+06 2.437320e+06 1.759963e+06
## [21] 1.214976e+06 7.785251e+05 4.590441e+05 2.428199e+05 1.107300e+05
## [26] 4.060476e+04 1.040537e+04 1.320942e+03 7.239578e-09 -5.019224e-09
# zap negative eigenvalues making them zero
R <- with(e,vectors %*% diag(pmax(values,0)) %*% t(vectors))
beta_var <- Variable(nrow(R))
q <- quad_form(beta_var,R)
is_dcp(q)
## [1] TRUE
,
library(data.table)
library(magrittr)
station_photos <- "
year_unit_station Photo_Number Creation_Datetime bin_name
1: 2016_275_02 275_02_0017.JPG 2016-09-23 11:51:03
2: 2016_275_02 275_02_0035.JPG 2016-09-27 15:58:21
3: 2016_275_02 275_02_0036.JPG 2016-09-27 15:58:49
4: 2016_275_02 275_02_0037.JPG 2016-09-27 16:00:04
5: 2016_275_02 275_02_0038.JPG 2016-09-27 16:00:59
6: 2016_275_02 275_02_0039.JPG 2016-09-27 16:01:27
7: 2016_275_02 275_02_0062.JPG 2016-10-02 12:22:35
8: 2016_275_02 275_02_0075.JPG 2016-10-31 03:09:43" %>%
readr::read_fwf(col_types = "-ccc") %>%
setDT() %>%
setnames(.[1,unlist(.SD)]) %>%
.[-1] %>%
.[,Creation_Datetime := anytime::anytime(Creation_Datetime)]
station_bins <- "
year_unit_station service_end_dttm bin_name
1: 2016_275_02 2016-09-23 11:21:00 2016_275_02_1
2: 2016_275_02 2016-09-30 10:45:00 2016_275_02_2
3: 2016_275_02 2016-10-07 08:31:00 2016_275_02_3" %>%
readr::read_fwf(col_types = "-ccc") %>%
setDT() %>%
setnames(.[1,service_end_dttm := anytime::anytime(service_end_dttm)]
应该给出与 G. Grothendieck 暗示的 cvxpy
相同的结果,因为 DCP 规则是相同的。最近版本的 CVXR
似乎有问题。我在 cvxpy github 上打开了一个问题。