求信赖域算法加局部二次近似求解 关于b级数i从1到nxi[yi-(e∧bxi/1+e∧bxi)])的最大值(b是P元向量)的Python代码
时间: 2023-12-03 10:46:50 浏览: 41
以下是使用信赖域算法和局部二次近似求解最大值的Python代码:
```python
import numpy as np
from scipy.optimize import minimize
def trust_region(x, grad, hess, radius):
"""
信赖域算法
"""
delta = 0.1
eta = 0.1
max_iter = 100
tol = 1e-6
for i in range(max_iter):
# 计算牛顿步
p = -np.linalg.solve(hess(x), grad(x))
# 计算预测下降
pred = -grad(x).dot(p) - 0.5 * p.dot(hess(x).dot(p))
# 计算实际下降
actual = lambda s: -(x + s * p).dot(grad(x + s * p)) - 0.5 * (x + s * p).dot(hess(x + s * p).dot(x + s * p))
actual0 = -grad(x).dot(grad(x)) / (2 * grad(x).dot(hess(x).dot(grad(x))))
if actual0 < 0:
tau = 1
else:
tau = min(1, delta * np.linalg.norm(grad(x)) / actual0)
s = minimize(actual, 0, bounds=[(-radius, radius)], tol=tol).x[0]
actual_val = actual(s)
if actual_val < eta * pred or np.abs(np.linalg.norm(grad(x + s * p)) - radius) < tol:
x = x + s * p
if np.linalg.norm(p) < tol or np.linalg.norm(grad(x)) < tol:
break
if actual_val > eta * pred:
radius = min(np.linalg.norm(p), 2 * radius)
else:
radius = max(radius / 2, tol)
return x
def obj_func(b, x, y):
"""
目标函数
"""
n = len(x)
fx = sum([y[i] - np.exp(np.dot(b, x[i])) / (1 + np.exp(np.dot(b, x[i]))) for i in range(n)])
return -fx
def hess_func(b, x, y):
"""
计算Hessian矩阵
"""
n = len(x)
H = np.zeros((len(b), len(b)))
for i in range(n):
expbx = np.exp(np.dot(b, x[i]))
p1 = expbx / ((1 + expbx) ** 2) * np.outer(x[i], x[i])
p2 = expbx / (1 + expbx) * np.eye(len(b))
H += p1 - p2
return -H
def trust_region_newton(b0, x, y, radius):
"""
信赖域牛顿法
"""
grad = lambda b: np.sum([x[i] * (y[i] - np.exp(np.dot(b, x[i])) / (1 + np.exp(np.dot(b, x[i])))) for i in range(len(x))], axis=0)
hess = lambda b: hess_func(b, x, y)
b_star = trust_region(b0, grad, hess, radius)
return b_star
# 测试代码
if __name__ == '__main__':
x = np.array([[1, 2], [3, 4], [5, 6]])
y = np.array([1, 0, 1])
b0 = np.array([0, 0])
radius = 1
b_star = trust_region_newton(b0, x, y, radius)
print(b_star)
```
其中,`trust_region`函数是信赖域算法,`obj_func`函数是目标函数,`hess_func`函数是计算Hessian矩阵的函数,`trust_region_newton`函数是信赖域牛顿法。在测试代码中,我们使用了一个简单的例子来测试函数是否正常工作。
相关推荐
![zip](https://img-home.csdnimg.cn/images/20210720083736.png)
![zip](https://img-home.csdnimg.cn/images/20210720083736.png)
![](https://img-home.csdnimg.cn/images/20210720083646.png)
![](https://csdnimg.cn/download_wenku/file_type_ask_c1.png)
![](https://csdnimg.cn/download_wenku/file_type_ask_c1.png)
![](https://csdnimg.cn/download_wenku/file_type_ask_c1.png)
![](https://csdnimg.cn/download_wenku/file_type_ask_c1.png)
![pptx](https://img-home.csdnimg.cn/images/20210720083543.png)
![docx](https://img-home.csdnimg.cn/images/20210720083331.png)
![docx](https://img-home.csdnimg.cn/images/20210720083331.png)
![pdf](https://img-home.csdnimg.cn/images/20210720083512.png)
![pdf](https://img-home.csdnimg.cn/images/20210720083512.png)
![zip](https://img-home.csdnimg.cn/images/20210720083736.png)
![pdf](https://img-home.csdnimg.cn/images/20210720083512.png)
![zip](https://img-home.csdnimg.cn/images/20210720083736.png)
![rar](https://img-home.csdnimg.cn/images/20210720083606.png)
![zip](https://img-home.csdnimg.cn/images/20210720083736.png)