我有一个最小二乘问题,无需任何已知的参数估计即可解决。我施加了一个约束,即我所需的解是平滑的(模型参数缓慢变化),因此我将相邻参数之间的差异最小化(用于此地质问题的传统疗法)。
通过将约束方程式排列为原始数据方程式d = Gm中的行来实现约束。辅助参数w是通过反复试验选择的(某些教科书将w称为拉格朗日乘数)。
我有以下内容:
G = np.array([[1,0,1,0,0,6],
[1,0,0,1,0,6.708],
[1,0,0,0,1,8.485],
[0,1,1,0,0,7.616],
[0,1,0,1,0,7],
[0,1,0,0,1,7.616]])
d = np.array([[2.323],
[2.543],
[2.857],
[2.64],
[2.529],
[2.553]])
现在添加任意w加权平滑度的约束(w = 0.01):
w = 0.01
G = np.array([[1,0,1,0,0,6],
[1,0,0,1,0,6.708],
[1,0,0,0,1,8.485],
[0,1,1,0,0,7.616],
[0,1,0,1,0,7],
[0,1,0,0,1,7.616],
[w,-w,0,0,0,0],
[0,w,-w,0,0,0],
[0,0,w,-w,0,0],
[0,0,0,w,-w,0],
[0,0,0,0,w,-w]])
d = np.array([[2.323],
[2.543],
[2.857],
[2.64],
[2.529],
[2.553],
[0],
[0],
[0],
[0],
[0]])
但是,为w选择合适的值似乎是限制模型参数良好解决方案的关键步骤。
所以我的问题是:使用Python,有没有办法让我循环遍历许多具有不同w值的计算出的解决方案,并选择用于获得最佳质量的解决方案的值?
在提出的解决方案中,我将G_0称为G,没有附加约束,类似地,d_0是d,没有附加零。我还假设您正在从某处读取G_0和d_0,并且我称它们为已知。
import numpy as np
def create_W(n_rows, w):
W = -np.diagflat(np.ones(n_rows), 1)
np.fill_diagonal(W, 1)
return W
def solution_quality_metric(m):
# this need to be implemented to determine what you mean by "best"
n_rows = 5
d_w = np.zeros(n_rows)
# choose range for w values for example w_min = 0, w_max = 1, dw = 0.01
best_m = -np.inf
best_w = w_min
for w in np.arange(w_min, w_max, dw):
W = create_W(n_rows, w)
G = np.concatenate([G_0, W], axis=0)
d = np.concatenate([d_0, d_w])
m = np.lstsq(G, d)
if solution_quality_metric(m) > best_m:
best_m = solution_quality_metric(m)
best_w = w
此代码显然将无法按原样工作,因为您未指定“具有最佳质量的解决方案”的含义。为此,您需要实现solution_quality_metric
函数