python – NumPy中最小二乘算法的高效计算
作者:互联网
我需要在最小二乘意义上解决大量的线性系统.我在理解numpy.linalg.lstsq(a,b),np.dot(np.linalg.pinv(a),b)和数学实现的计算效率差异方面遇到了麻烦.
我使用以下矩阵:
h=np.random.random((50000,100))
a=h[:,:-1].copy()
b=-h[:,-1].copy()
并且算法的结果是:
# mathematical implementation
%%timeit
np.dot(np.dot(np.linalg.inv(np.dot(a.T,a)),a.T),b)
10个循环,最佳3:36.3 ms每个循环
# numpy.linalg.lstsq implementation
%%timeit
np.linalg.lstsq(a, b)[0]
10个循环,最佳3:每循环103毫秒
%%timeit
np.dot(np.linalg.pinv(a), b)
1个循环,最佳3:216 ms每个循环
为什么会有区别?
解决方法:
例程lstsq
处理任何系统:过度确定,未确定或未确定.它的输出是你从pinv(a)* b得到的,但它比计算伪逆更快.原因如下:
一般建议:除非您确实需要,否则不要计算逆矩阵.为特定右手侧求解系统比反转其矩阵更快.
然而,即使您正在反转矩阵,通过求解aTa = aTb的方法也会更快.是什么赋予了?反之,反转aTa仅在具有完整列等级时才有效.所以你已经把这个问题局限于这种特殊情况,并且作为一种权衡获得了不那么普遍的权利,并且如下所示,安全性较低.
但是反转矩阵仍然是低效的.如果您知道a具有完整的列排名,则以下内容比您的三次尝试中的任何一次更快:
np.linalg.solve(np.dot(a.T, a), np.dot(a.T, b))
也就是说,当处理条件差的矩阵时,lstsq仍然优于上述.形成产品aTa基本上对条件数进行平方,因此您更有可能获得无意义的结果.这是一个警示的例子,使用SciPy的linalg模块(它基本上等同于NumPy但有更多方法):
import numpy as np
import scipy.linalg as sl
a = sl.hilbert(10) # a poorly conditioned square matrix of size 10
b = np.arange(10) # right hand side
sol1 = sl.solve(a, b)
sol2 = sl.lstsq(a, b)[0]
sol3 = sl.solve(np.dot(a.T, a), np.dot(a.T, b))
这里lstsq提供与solve(该系统的唯一解决方案)几乎相同的输出.然而,sol3完全错误是因为数字问题(你甚至不会被警告).
SOL1:
[ -9.89821788e+02, 9.70047434e+04, -2.30439738e+06,
2.30601241e+07, -1.19805858e+08, 3.55637424e+08,
-6.25523002e+08, 6.44058066e+08, -3.58346765e+08,
8.31333426e+07]
SOL2:
[ -9.89864366e+02, 9.70082635e+04, -2.30446978e+06,
2.30607638e+07, -1.19808838e+08, 3.55645452e+08,
-6.25535946e+08, 6.44070387e+08, -3.58353147e+08,
8.31347297e+07]
SOL3:
[ 1.06913852e+03, -4.61691763e+04, 4.83968833e+05,
-2.08929571e+06, 4.55280530e+06, -5.88433943e+06,
5.92025910e+06, -5.56507455e+06, 3.62262620e+06,
-9.94523917e+05]
标签:python,numpy,least-squares 来源: https://codeday.me/bug/20190519/1133890.html