如何用scipy.optimize.nnls获取非负最小二乘的多组解？

首先得明确两个核心关键点：

1. 你当前用的scipy.optimize.nnls本身就是返回唯一的最小残差非负解——这是算法的固有特性，它会找到使||C.T x - line||_2最小的非负x，这个最优解是唯一的，所以无论跑多少次结果都一样，完全属于正常现象。
2. 你提到需要“系数和为1”的约束，但当前代码里根本没加这个限制！现在得到的解只是碰巧和接近1，并非严格满足sum(x)=1的要求。

怎么生成多组符合要求的解？ #

我们的目标是找到多个满足以下条件的x：

• x ≥ 0（非负）
• sum(x) = 1（系数和为1）
• ||C.T x - line||_2 ≤ ε（ε是你设定的残差容差，要大于最小残差）

下面给你两种实用的实现方法：

方法一：随机初始化可行点+带约束优化 #

我们用scipy.optimize.minimize（SLSQP方法）求解带约束的问题，通过随机初始化满足sum(x)=1的非负向量，每次得到不同的可行解：

import numpy as np
from scipy.optimize import minimize, nnls

def residual_obj(x, C_T, target):
    # 定义残差平方和的目标函数
    return np.linalg.norm(C_T @ x - target)**2

def generate_multiple_valid_solutions(C_T, target, num_sols=5, tol=1e-2):
    # 先求无约束（仅非负）的最小残差，作为容差基准
    _, min_residual = nnls(C_T, target)
    # 允许残差平方比最小值大tol
    max_allowed_residual_sq = min_residual + tol

    valid_solutions = []
    for _ in range(num_sols):
        # 用Dirichlet分布随机生成初始点：天然满足x≥0且sum(x)=1
        x_init = np.random.dirichlet(np.ones(C_T.shape[1]))
        
        # 设置约束条件
        constraints = [
            # 系数和为1的等式约束
            {'type': 'eq', 'fun': lambda x: np.sum(x) - 1},
            # 残差平方不超过允许值的不等式约束
            {'type': 'ineq', 'fun': lambda x: max_allowed_residual_sq - residual_obj(x, C_T, target)}
        ]
        # 非负边界约束
        bounds = [(0, None) for _ in range(C_T.shape[1])]
        
        # 求解优化问题
        result = minimize(
            residual_obj,
            x0=x_init,
            args=(C_T, target),
            method='SLSQP',
            bounds=bounds,
            constraints=constraints
        )
        
        if result.success:
            valid_solutions.append(result.x)
    
    return valid_solutions

# 测试你的数据
if __name__ == '__main__':
    C = np.array([[112.771820, 174.429720, 312.175750, 97.348620],
                  [112.857010, 174.208300, 312.185270, 93.467580],
                  [114.897210, 175.661850, 314.275100, 99.015480]
                ])
    d = np.array([[112.7718, 174.4297, 312.1758, 97.3486],
                  [112.7718, 174.4297, 312.1758, 97.3486]])
    
    for idx, line in enumerate(d):
        print(f"=== 处理第{idx+1}组目标向量 ===")
        # 原nnls结果（无sum约束）
        orig_sol, _ = nnls(C.T, line)
        print(f"原nnls解（无sum约束）：{orig_sol}，和为{np.sum(orig_sol):.6f}")
        
        # 生成多组带sum=1的解
        solutions = generate_multiple_valid_solutions(C.T, line, num_sols=3, tol=1e-2)
        for sol_idx, sol in enumerate(solutions):
            residual = np.linalg.norm(C.T @ sol - line)
            print(f"可行解{sol_idx+1}：{sol}，和为{np.sum(sol):.6f}，残差{residual:.6f}")

方法二：扰动目标向量生成近似解 #

如果你的容差要求不那么严格，可以给目标向量line加入小的随机噪声，然后求解带sum=1约束的非负最小二乘，这样每次扰动不同，得到的解也不同：

def perturbed_target_solutions(C_T, target, num_sols=5, noise_scale=1e-4):
    solutions = []
    for _ in range(num_sols):
        # 给目标向量加小噪声
        noisy_target = target + np.random.normal(0, noise_scale, size=target.shape)
        
        # 定义带sum约束的目标函数
        def obj(x):
            return np.linalg.norm(C_T @ x - noisy_target)**2
        
        constraints = [{'type': 'eq', 'fun': lambda x: np.sum(x)-1}]
        bounds = [(0, None)]*C_T.shape[1]
        x_init = np.random.dirichlet(np.ones(C_T.shape[1]))
        
        result = minimize(obj, x0=x_init, method='SLSQP', bounds=bounds, constraints=constraints)
        if result.success:
            solutions.append(result.x)
    return solutions

注意事项 #

• 调整tol（容差）或noise_scale（噪声幅度）可以控制解的多样性：值越大，解的差异越明显，但残差也会越大。
• 如果优化失败（result.success为False），可以尝试调整初始化方式、容差或者换用其他优化方法（比如trust-constr）。

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