如何在Python中替代Excel Solver调整两个变量最小化误差项

需求说明 #

我正在将旧的VBA代码转换为Python脚本，核心需求是通过调整asset_val（资产价值）和asset_vol（资产波动率）两个列的取值，最小化误差项。

原VBA核心代码 #

原VBA调用Excel的Solver功能逐行完成优化，其中第28列为误差值，第20、21列分别为初始资产价值、初始资产波动率：

For k = startobs To (endobs - 1) Step 1

SolverReset
SolverOK SetCell:=Worksheets("Select").Cells(k, 28).address, MaxMinVal:=2, ByChange:=Range(Cells(k, 20), Cells(k, 21)).address
SolverReset
SolverOK SetCell:=Worksheets("Select").Cells(k, 28).address, MaxMinVal:=2, ByChange:=Range(Cells(k, 20), Cells(k, 21)).address
    
SolverOptions precision:=0.00001, Iterations:=300, Scaling:=True, Convergence:=0.00001

现有Python基础实现 #

误差项计算 #

误差项综合考虑期权定价偏差和波动率偏差，计算公式如下：

def Error(S, K, T, r, sigma, E, ewma_vol):
    return ((BnSCall(S, K, T, r, sigma) - E)**2 + (SigmaByEquity(S, K, T, r, sigma) - ewma_vol * E)**2)

EWMA（指数加权移动平均）波动率计算 #

def EstVolEWMA(log_ret, lamda): 
    est_vol_ewma = abs(log_ret*0)
    est_vol_ewma_prev = 0           
    for date in est_vol_ewma.index[50:]:
        est_vol_ewma_t = (lamda*(est_vol_ewma_prev**2)+(1-lamda)*(log_ret[date]**2)*52)**0.5 
        est_vol_ewma[date] = est_vol_ewma_t
        est_vol_ewma_prev = est_vol_ewma_t
    return est_vol_ewma

布莱克-斯科尔斯看涨期权定价 #

def BnSCall(S,K,T,r,sigma):
    return S*norm.cdf(BnSd1(S,K,T,r,sigma), 0, 1)-K*np.exp(-np.log(1+r)*T)*norm.cdf(BnSd2(S,K,T,r,sigma))

初始值计算 #

初始资产价值

def InitAssetVal(debt, avv_mcap): 
    return debt + avv_mcap
data['asset_val'] = InitAssetVal(data['debt interpolated'], data['av_m_cap'])

初始资产波动率

def InitAssetVol(equity_vol): 
    return equity_vol*0.5
data['asset_vol'] = InitAssetVol(data['est_vol_ewma'])

Python优化实现（替代VBA Solver功能） #

使用scipy.optimize.minimize实现和VBA Solver参数完全对齐的逐行优化逻辑，代码如下：

import scipy.optimize as opt
import numpy as np

# 提前计算全量EWMA波动率，避免优化时重复计算
data['ewma_vol'] = EstVolEWMA(data['log_ret'], lamda=0.94)
# 替换为实际需要优化的观测范围
startobs = 50
endobs = len(data)
optimize_range = data.index[startobs: endobs-1]

for idx in optimize_range:
    # 提取当前行固定参数
    row = data.loc[idx]
    K = row['debt interpolated']
    T = 1
    r = row['12m rate'] / 100
    E = row['av_m_cap']
    ewma_vol = row['ewma_vol']
    
    # 构造当前行专属误差函数，入参为[资产价值, 资产波动率]
    def row_error(x):
        S = x[0]
        sigma = x[1]
        # 资产价值、波动率不能为负，异常值返回无穷大
        if S <= 0 or sigma <= 0:
            return np.inf
        return Error(S, K, T, r, sigma, E, ewma_vol)
    
    # 取预计算的初始值
    x0 = [row['asset_val'], row['asset_vol']]
    # 求解参数完全对齐原VBA Solver配置
    res = opt.minimize(row_error, x0, method='L-BFGS-B', 
                      bounds=((1e-6, None), (1e-6, None)),
                      tol=1e-5,
                      options={'maxiter': 300, 'disp': False})
    # 写入优化结果
    data.loc[idx, 'asset_val'] = res.x[0]
    data.loc[idx, 'asset_vol'] = res.x[1]
    data.loc[idx, 'Error'] = res.fun

注意事项 #

• 需提前实现BnSd1、BnSd2、SigmaByEquity三个依赖函数，保证逻辑和原实现一致
• 若优化收敛效果不佳，可将求解方法替换为Nelder-Mead适配非平滑误差场景
• 批量优化速度较慢时，可通过joblib等并行库对逐行逻辑做并行加速

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