SciPy是Python编程语言的一个数学程序库,它为建模和解决科学问题提供了构建块。SciPy包括优化、积分、插值、特征值问题、代数方程、微分方程和许多其他类问题的算法;它还提供特殊的数据结构,如稀疏矩阵和k-维树。SciPy是建立在NUMPY之上的,它提供了array数据结构和相关的快速的数学程序,而SciPy本身就是构建更高级别的科学库的基础,包括scikit-learn和scikit-image。在这篇文章中,我们概述了SciPy 1.0的功能和开发实践,并重点介绍了最近的一些技术进展。
截至2019年2月,SciPy库包含将近600,000行开源代码,这些代码按图1汇总的16个子包进行组织。开发团队和社区目前主要在GitHub上进行交互和操作,这是一个在线版本控制和任务管理平台。超过110,000个GitHub存储库和6,500个软件包依赖SciPy。
图1:SciPy包组织
一、架构和实现
1.1 项目范围
SciPy提供了用于科学计算的基本算法。在发展相对较缓慢的领域,SciPy旨在提供完整的覆盖范围。在其他领域,SciPy的目的是提供基本的构建块,同时与该领域的其他包进行良好的交互。例如,SciPy提供了概率分布,假设检验,频率统计,相关函数等,而Statsmodels提供了更高级的统计估计量和推断方法,scikit- learn3涵盖了机器学习,而PyMC3,emcee和PyStan (http://mc-stan.org) 涵盖了贝叶斯统计和概率建模。scikit-image4提供了比SciPy的ndimage模块更加优秀的的图像处理功能,SymPy提供了用于符号计算的Python交互界面,而sparse. csgraph和spatial与诸如NetworkX之类的专用库相比提供了与图形和网络一起使用的基本工具。
我们使用以下标准来确定是否要在SciPy中加入新特征:
•该算法与多个科学领域相关。
•该算法非常重要。例如,它足够经典,可以包含在教科书中,或者基于经过大量同行引用的同行评审文章。
在软件系统和体系结构方面,SciPy的范围与NumPy的范围相匹配。
1.2 语言选择
根据linguist库 (https://github.com /github/linguist) 的分析,SciPy大约是Python占50%,Fortran占25%,C占20%,Cython占3%和C ++占2%,以及少量的TeX, Matlab,shell 脚本和Make。
尽管Fortran已有很长的历史,但它仍然是一种高性能的科学编程语言。因此,我们包装了以下出色的且经过实测的Fortran库,以在为Python带来便利的同时提高性能:QUADPACK和ODEPACK用于数值积分和初值问题的解决;FITPACK,ODRPACK和MINPACK用于曲线拟合和最小二乘最小化;FFTPACK用于执行傅立叶变换;ARPACK用于解决特征值问题;用于计算贝塞尔函数的ALGORITHM 644;和CDFLIB用于评估累积密度函数。
在SciPy中排名前三位的语言还有一个是C语言,我们在SciPy中包装的C库包括:用于优化的trlib,用于解决稀疏线性系统的SuperLU,用于计算几何图形的Qhull和用于特殊函数的Cephes (http://www.netlib.org /cephes/)。
Cython被描述为一种混合语言,它混合了Python的最佳部分和低级C/C++范例。我们经常使用Cython作为C/C++编写的成熟的低级科学计算库与SciPy提供的Python接口之间的粘合剂。我们还使用Cython在Python代码中实现性能增强,特别是对于大量使用的内部循环受益于具有静态类型的编译代码的情况。
为了实现新功能,Python仍然是首选的语言。如果Python性能存在问题,那么我们更喜欢使用Cython,然后使用C,C++或Fortran(按此顺序)。其主要考虑因素是可维护性:Cython具有最高的抽象级别,大多数Python开发人员都会理解它。C也是很常用的,比起C++尤其是Fortran,当前的核心开发团队更容易管理。
1.3 API和ABI的演变
SciPy的API由大约1500个函数和类组成。我们优化API的策略是添加新功能,而只有在优点超过用户的损失以及向这些用户提供明确的舍弃警告至少一年之后才能删除或更改现有功能。一般而言,我们鼓励修改从而提高库的API的清晰度,但不建议破坏向后兼容性,因为SciPy位置接近科学Python计算堆栈的基础。
除了Python API之外,SciPy还有C和Cython接口。因此,我们还必须考虑应用程序二进制接口(ABI)。这个ABI已经稳定了很长时间,我们的目标是以向后兼容的方式发展它。
二、关键技术改进
在这里,我们描述了过去三年中取得的关键技术改进。
2.1 数据结构
稀疏矩阵。scipy. sparse提供了七个稀疏矩阵数据结构,也称为稀疏格式。最重要的是行和列压缩格式(分别为CSR和CSC)。它们提供快速的轴索引和快速的矩阵向量乘法,在整个SciPy和相关程序包中大量使用。
在过去三年中,我们对稀疏矩阵处理内部进行了重写,提高了性能。现在,CSC和CSR矩阵的迭代和切片速度最高可提高35%,并且COO/DIA转换为CSR/ CSC格式的速度有所提高。SuperLU已更新至版本5.2.1,增强了部分稀疏矩阵功能的低级实现。
从新功能的角度来看,scipy.sparse矩阵和线性运算操作现在支持Python矩阵乘法(@)运算符。我们添加了scipy. sparse. norm和scipy. sparse. random分别用于计算稀疏矩阵范数和从任意分布中得出随机变量。此外,我们尽力使scipy. sparse API与NumPy API保持一致。
cKDTree。scipy. spatial. ckdtree模块实现了组织k维空间中的点的空间分区数据结构,并使用模板化类在C++中进行了重写。增加了对周期性边界条件的支持,该条件通常用于物理过程的仿真中。
2013年,cKDTree.query的k最近邻搜索的时间复杂度大约为loglinear,与其正式描述一致。从那时起,我们增强了cKDTree,通过在C ++中重新实现该查询,消除内存泄漏并允许释放全局解释器锁(GIL),以便可以使用多个线程来进行查询。这通常可以在保留渐近复杂性的同时提高任何给定问题的性能。
在2015年,SciPy添加了sparse_ distance _matrix例程,通过忽略所有超过用户提供的值的距离来生成KDTree对象之间的近似稀疏距离矩阵。该例程不限于常规的L2(欧几里得)范数,而是支持1到无穷大之间的任何Minkowski p-norm。默认情况下,返回的数据结构是基于字典(DOK)的稀疏矩阵,这对于矩阵构造非常有效。这种用于稀疏矩阵汇编的哈希方法可以比使用CSR 格式构造快7倍,并且C++级别的稀疏矩阵构造释放了Python GIL以提高性能。
2015年,对cKDTree算法进行了增强,以支持权重,这在许多科学应用中都是必不可少的,例如,计算星系的相关函数。
2.2 对已编译代码的统一绑定
LowLevelCallable。从SciPy 0.19版开始,用户可以将低级函数包装在scipy. LowLevelCallable对象中,从而减少了直接从Python调用已编译的C函数(例如使用Numba或Cython生成的函数)的开销。支持的低级功能包括PyCapsule对象,ctypes函数指针和cffi函数指针。此外,可以从Cython模块使用scipy. LowLevelCallable. from _ cython自动生成低级回调函数。
2.3 BLAS,LAPACK和special的Cython绑定
SciPy现在包括许多BLAS和LAPACK程序以及scipy. special子包中提供的特殊函数的Cython包装器。这些用于Cython的低级接口也可以在SciPy代码库之外使用,以访问包装的库中的函数,同时避免Python函数调用的开销。对于许多用例,这可以使性能提高一两个数量级。开发人员也可以使用低级Cython接口,而无需链接到所包装的库。这使其他扩展功能避免了查找和使用正确库的复杂性。
包装以Fortran编写的库时,避免这种复杂性尤为重要。这些低级包装器不仅可以在没有Fortran编译器的情况下使用,而且还可以在无需处理所有不同的Fortran编译器ABI和名称修改方案的情况下使用。这些低级Cython包装器大多数都是自动生成的,以保证正确性和易于维护性。
由于SciPy可以使用LAPACK 3.4.0或更高版本构建,因此仅为一致的程序提供Cython包装器所有支持的LAPACK版本之间的接口。各种现有BLAS库提供的标准BLAS接口当前未更改,因此SciPy提供的包装程序通常不需要更改。scipy. special的Cython包装器遵循对该子包的接口的更改进行相应更改。
2.4 数值优化
scipy. optimize子程序包提供了一些用于求根和优化问题函数。在这里,我们重点介绍SciPy1.0的最新功能。
线性优化。与SciPy 1.0一起发布了一种新的针对连续线性编程问题的内点优化器linprog,method = ’interior-point’,实现了商用求解器MOSEK的核心算法,它可以解决所有90多个经过测试的NETLIB LP基准测试问题。
预求解程序解决了一些琐碎的问题,在形式上简化了问题,例如束缚收紧和移除固定变量,并且自动选择了一些用于消除冗余等式约束的程序,以减少由奇异矩阵引起的数值计算困难的可能性。尽管主要的求解器是由Python实现,但由于端到端稀疏矩阵的支持和SciPy编译的线性系统求解器的大量使用,可提供足够的速度来解决有成千上万的变量和约束的问题。
非线性优化:局部最小化。minimize功能提供了一个统一的接口,用于查找非线性优化问题的局部最小值。在最新版本的SciPy中,在minimize中添加了四种用于无约束优化的新方法:dogleg,trust-ncg,trust-exact和trust-krylov。所有方法都是信赖域方法,它们基于一阶和二阶导数信息建立目标函数的局部模型,逼近局部“信赖域”内的最佳点,并进行迭代,直到达到原始目标函数的局部最小值为止,但是每种方法都有其独特的特征,使其适合某些类型的问题。例如,trust-exact通过几乎精确地解决信赖区域子问题来实现快速收敛,但是它要求存储第二阶导数Hessian矩阵并在每次迭代中将其分解,这可能会影响大规模问题(≥1,000个变量)的求解。相反,trust-ncg和trust-krylov非常适合大规模优化问题,因为它们不需要显式存储和分解Hessian矩阵,而可以以更快的近似方式使用二阶导数信息。我们在表1中详细比较了所有最小化方法的特性,这些特性说明了SciPy涵盖数值方法时所要达到的完整性水平。
表1:minimize的优化算法
非线性优化:全局最小化。由于最小化可能返回任何局部最小值,某些问题需要使用全局最小值。新的scipy. optimize. dif - ferentialevolution函数是一种随机全局优化器,它通过演化大量候选解而起作用。
在每个迭代中,通过组合现有组中的候选来生成试验候选组。如果候选试验代表一种优化,那么将更新组。最近,SciPy基准测试套件获得了196个全局优化问题的综合集合,这些问题可以随着时间的推移跟踪现有求解器的性能,并评估新求解器的性能是否值得将其包含在软件包中。
2.5 统计分布
scipy. stats子程序包包含100多个概率分布:96个连续分布和13个离散单变量分布,以及10个多元分布。
该实现依赖于一个一致的框架,该框架提供了对随机变量进行采样,评估累积分布函数(CDF)和概率密度函数(PDF)以及为每个分布拟合参数的方法。通常,这些方法依赖于每种分布的特定实现,例如CDF的闭式表达式或采样算法(如果可用)。否则,将基于通用代码使用默认方法,例如,对PDF进行数值积分以获得CDF。最近添加到scipy. stats中的主要分布包括在scipy. stats中基于直方图的分布,rv _ histogram和scipy中的多项式分布stats. multinomial(用于自然语言处理)
2.6 多项式插值器
从历史上看,SciPy严重依赖于P. Dierckx的FITPACK Fortran库进行数据的单变量插值和逼近,但是SciPy和FITPACK之间交互的原始单片设计和API限制了用户和开发人员。
实现多项式插值器的一种新的模块化设计的过程遍布多个版本。这项工作的目标是要有一组代表分段多项式的基本对象,以实现用于构建各种插值器的算法的集合,并为用户提供用于构建其他插值器的构件。
在新设计的最低层是代表单变量分段多项式的类:PPoly (SciPy 0.13) ,BPoly (SciPy 0.13) 和BSpline (SciPy0.19) ,它们允许进行有效的矢量化评估,微分,积分和求根。PPoly以幂为单位,以断点和每个区间的系数表示分段多项式。BPoly是类似的,并且在Bernstein基础上表示分段多项式(例如,构造Bézier曲线)。BSpline表示样条曲线,即B样条基础元素的线性组合。
在下一层中,这些多项式类用于构造几种常见的数据插值方式:CubicSpline (SciPy 0.18) 构造一个二次可微的分段三次函数,Akima1DInterpolator和PCHIPInterpolator用于构造C连续单调形状保全插值器。
2.7 测试和基准套件
测试套件。对于SciPy的每个组件,我们编写了多个小型可执行测试,以验证其预期的行为。根据持续集成的实践,在运行测试套件之前,所有对SciPy的建议贡献都将与库的主分支暂时集成在一起,并且在永久合并贡献之前必须通过所有测试。持续监视单元测试所涵盖的SciPy中的代码行数是我们可以确定更改和正确实现新功能的一种方法。
SciPy测试套件由一个连续集成矩阵编排而成,该矩阵包括分别由Travis CI和AppVeyor管理的POSIX和Windows (32/64位) 平台。我们的测试涵盖Python版本2.7、3.4、3.5、3.6,并使用pyflakes和pycodestyle进行代码替换。测试套件中有13,000多个单元测试,这些单元测试是为与pytest (https:// docs.pytest.org /en /latest) 框架一起使用而编写的。在图2中,我们显示了使用基于Docker的方法 (https:// github.com/ tylerjereddy/ scipy-cov-track) 生成的历史测试覆盖率数据。根据pytest-cov (https:// pypi.org/ project/ pytest-cov/) ,对于Python代码,SciPy 1.0发行点的测试覆盖率为87%。CircleCI服务会自动生成和发布代码文档,以便于评估文档更改/完整性。
图2:Python和SciPy中的编译代码量
基准套件。除了确保单元测试通过外,重要的是要确认SciPy代码库的性能会随着时间的推移而提高。自2015年2月以来,SciPy的性能已通过Airspeed Velocity (asv https:// github.com/ airspeed-velocity/ asv) 进行监控。SciPy的run.py脚本可以方便地包装asv功能,以便可以通过单个控制台命令随时间生成基准测试结果。例如,在图3中,展示了scipy. spatial的改进。cKDTree. query大约有九年的项目历史。在基准测试中使用的树是在没有应用环形拓扑的情况下生成的 (boxsize = None) ,并且使用Airspeed Velocity 0.4,Python 2.7,NumPy 1.8.2和Cython 0.27.3、0.21.、0.18版本进行了测试,向后兼容。将cKDTree完全进行Cython化后,再用C++对其进行重写时,实现了显着的性能改进。
图3:从cKDTree的引入到SciPy 1.0发行的基准测试结果
三、讨论
SciPy具有强大的开发者社区和庞大的用户群。尽管如此,SciPy一直在努力改进。许多开源项目面临的问题是吸引和留住开发人员。因为代码交接过多会导致记忆丧失,从而导致过去的错误重复出现。幸运的是,SciPy项目继续吸引着热情而有才干的新开发人员,同时保持了规模虽小但专职的前期贡献者的参与。
要解决的最后一个重要挑战是如何容纳GPU和分布式计算,而又不会破坏我们传统的架构。尽管我们仍不清楚如何在整个库中采用这些新兴技术的确切方法,但是我们现在有一个具体的子包实现,该子包允许实验性地使用多个后端,这将在以后的报告中详细描述。