CVXPY下2×2矩阵混合整数SDP可否转换为MISOCP求解

可以的，利用2×2实对称矩阵半正定的充要条件，完全可以将原问题转换为混合整数二阶锥规划（MISOCP）求解，兼容性更好、求解精度也更高。

转换原理 #

对于2×2实对称矩阵A = [[a, c], [c, b]]，半正定（PSD）的充要条件可以表示为3个纯凸约束，无需调用半正定锥求解器：

• 对角线元素非负：a ≥ 0、b ≥ 0
• 行列式非负：a*b ≥ c²（该约束属于旋转二阶锥约束，所有支持SOCP的求解器都可原生处理）

目标函数max log_det(A)也可以对应转换为SOCP兼容的形式：对于2阶矩阵，det(A) = a*b - c²，最大化log(det(A))等价于最大化行列式的对数，CVXPY可自动将该目标转换为对应锥约束，不需要额外手动处理。

另外原示例代码存在一处笔误：off是控制每个样本是否为异常点的布尔变量，应该定义为长度为n的向量，而非标量，否则sum(off) ≤20的约束完全无效。

转换后可运行代码 #

import cvxpy as cp
import matplotlib.pyplot as plt
import numpy as np

np.random.seed(271828) 
m = 2; n = 50
x = np.random.randn(m,n)

# 修正off变量定义：每个样本对应一个布尔开关
off = cp.Variable(n, boolean=True)
# 定义2x2对称矩阵变量，不用直接标记PSD，后续用约束实现
a = cp.Variable(nonneg=True)
b = cp.Variable(nonneg=True)
c = cp.Variable()
A = cp.vstack([cp.hstack([a, c]), cp.hstack([c, b])])
b_shift = cp.Variable(2)

# 目标：最大化log_det(A)
obj = cp.Maximize(cp.log_det(A))
# 约束1：2x2矩阵半正定的充要条件
constraints = [a * b >= c ** 2]
# 约束2：椭球包含正常样本，异常样本不受约束（大M设置为20）
constraints += [cp.norm(A @ x[:,i] + b_shift) <= 1 + 20 * off[i] for i in range(n)]
# 约束3：异常点数量不超过20
constraints += [cp.sum(off) <= 20]

prob = cp.Problem(obj, constraints)
# 可替换为任意支持MISOCP的求解器，如MOSEK、GUROBI、SCIP等，精度远高于XPRESS
optval = prob.solve(solver='MOSEK', verbose=False)

print(f"最优值: {optval}")
print(f"A的最优解:\n{A.value}")
print(f"识别出的异常点数量: {int(np.sum(off.value))}")

# 绘制椭球和样本点
angles = np.linspace(0, 2*np.pi, 200)
rhs = np.row_stack((np.cos(angles) - b_shift.value[0], np.sin(angles) - b_shift.value[1]))
ellipse = np.linalg.solve(A.value, rhs)

plt.scatter(x[0,:], x[1,:], c=off.value, cmap='coolwarm', label='正常/异常样本')
plt.plot(ellipse[0,:].T, ellipse[1,:].T, color='black', label='最小体积椭球')
plt.xlabel('维度1')
plt.ylabel('维度2')
plt.title('带异常点的最小体积椭球拟合')
plt.legend()
plt.show()

效果说明 #

转换后的MISOCP问题可以用所有支持混合整数二阶锥规划的求解器求解，相比原MISDP形式的求解速度更快、精度更高，且开源求解器也可直接处理，不需要商用半正定规划求解器。

内容的提问来源于stack exchange，提问作者Richard