機(jī)器之心報(bào)道

機(jī)器之心編輯部

每年千萬(wàn)下載量，科學(xué)計(jì)算開源庫(kù) SciPy，你已經(jīng)是個(gè)成熟的小伙伴了。

作為科學(xué)計(jì)算中的中流砥柱，SciPy 從 2001 年到現(xiàn)在已經(jīng)走過(guò)了十九個(gè)年頭，它為最優(yōu)化、積分、微分方程等各種數(shù)值計(jì)算提供了完整的流程，也為科研分析人員提供了最好用與高效的開源庫(kù)。

前天 2 月 3 日，SciPy 的維護(hù)者在 Nature Methods 上發(fā)表了一篇論文，其回顧了 SciPy 發(fā)展的里程碑與關(guān)鍵技術(shù)。并借助 SciPy 1.0 這個(gè)成熟的象征，展現(xiàn)了當(dāng)前科學(xué)計(jì)算以及未來(lái)發(fā)展方向都是什么樣的。

• 論文地址：https://www.nature.com/articles/s41592-019-0686-2
• SciPy 項(xiàng)目地址：https://github.com/scipy/scipy

為什么 SciPy 那么重要？

SciPy 是一個(gè)面向 Python 的開源科學(xué)計(jì)算庫(kù)。自 2001 年首次發(fā)布以來(lái)，SciPy 已經(jīng)成為 Python 語(yǔ)言中科學(xué)算法的行業(yè)標(biāo)準(zhǔn)。該項(xiàng)目擁有超過(guò) 800 個(gè)獨(dú)特的代碼貢獻(xiàn)者，數(shù)以千計(jì)的相關(guān)開發(fā)包，和超過(guò) 150,000 個(gè)依賴存儲(chǔ)庫(kù)以及每年數(shù)以百萬(wàn)計(jì)的下載量。在下述簡(jiǎn)介中，會(huì)概述 SciPy 1.0 的功能和開發(fā)實(shí)踐，并著重闡述一些最新的技術(shù)發(fā)展與更新。（注：數(shù)據(jù)更新來(lái)自 Github 最新反饋）

SciPy 是一個(gè)基于 Python 和 Numeric 的開源代碼包，當(dāng)前模塊集包括了上圖中的內(nèi)容。

項(xiàng)目覆蓋范圍

SciPy 提供了科學(xué)計(jì)算的基本算法，這些算法涵蓋了現(xiàn)有數(shù)學(xué)軟件分類系統(tǒng)中的算法。在迭代相對(duì)緩慢的領(lǐng)域（如：線性代數(shù)），SciPy 旨在提供完整的算法覆蓋。

而在其他領(lǐng)域，它提供基本的構(gòu)件，并與該領(lǐng)域的其他軟件包進(jìn)行良好的互動(dòng)于兼容。例如，SciPy 提供了人們期望在統(tǒng)計(jì)學(xué)教科書中能找到的基本算法（概率分布、假設(shè)檢驗(yàn)、頻率統(tǒng)計(jì)、相關(guān)函數(shù)等），但 Statsmodels 提供了更先進(jìn)的統(tǒng)計(jì)預(yù)估及推斷方法。雖然 scikit-learn 涵蓋了機(jī)器學(xué)習(xí)，但 PyMC、emcee 和 PyStan 涵蓋了貝葉斯統(tǒng)計(jì)和概率建模等。

SciPy 能干什么

我們知道 SciPy 是一個(gè)用于數(shù)學(xué)、科學(xué)、工程領(lǐng)域的常用庫(kù)，可以處理插值、積分、最優(yōu)化、常微分方程數(shù)值解的求解、信號(hào)處理等問(wèn)題。因此自然科學(xué)領(lǐng)域絕大多數(shù)涉及計(jì)算的工作都能用它來(lái)完成，例如我們熟知的統(tǒng)計(jì)學(xué)習(xí)，擬合個(gè)分布、做了 K 最近鄰算法都是非常便捷的。

當(dāng)然目前新冠肺炎疫情廣受關(guān)注，研究者也可以用它模擬各種關(guān)鍵信息。例如之前中國(guó)疾病預(yù)防控制中心、國(guó)內(nèi)各省市疾控中心等機(jī)構(gòu)在新英格蘭醫(yī)學(xué)雜志發(fā)表的《新型冠狀病毒肺炎在中國(guó)武漢的早期傳播》論文。

在獲取數(shù)據(jù)之后，進(jìn)行各種統(tǒng)計(jì)學(xué)分析很多都可以用 Scipy 完成，具體而言：

• 研究者根據(jù)發(fā)病日期構(gòu)建傳染曲線；
• 使用對(duì)數(shù)高斯分布擬合暴露歷史和發(fā)病日期數(shù)據(jù)，估計(jì)潛伏期分布；
• 使用韋伯分布擬合發(fā)病日期、首次就診日期和住院日期，并估計(jì)發(fā)病離就診的時(shí)間間隔分布、發(fā)病離住院的時(shí)間間隔分布；
• 使用伽瑪分布擬合病例集群數(shù)據(jù)，從而估計(jì)人際傳播的時(shí)間間隔（serial interval）分布。

這些分析任務(wù)主要在于利用統(tǒng)計(jì)分布擬合對(duì)應(yīng)的數(shù)據(jù)，該肺炎論文的研究者采用 MATLAB 做的擬合。但實(shí)際上，scipy.status 包含了 100 多個(gè)概率分布，這些統(tǒng)計(jì)分析也能通過(guò) SciPy 完成。

面對(duì)洶涌的疫情，不論我們是有第一手?jǐn)?shù)據(jù)，還是從各網(wǎng)站爬取疫情信息，利用 SciPy 建模與分析都是非常好的選擇。

SciPy 發(fā)展里程碑

20 世紀(jì) 90 年代末期，美國(guó)梅奧醫(yī)學(xué)中心的博士生 Travis Oliphant 發(fā)布了一系列構(gòu)建于數(shù)值數(shù)組之上的包，并提供了用于信號(hào)處理、特殊函數(shù)、稀疏矩陣、正交、最優(yōu)化和快速傅里葉變換等的算法。

這些包中的 Multipack 是一組包裝了 Fortran 和 C 語(yǔ)言的擴(kuò)展模塊，用于解決非線性方程和最小二乘問(wèn)題、求微分方程的積分以及擬合曲線。

隨后，編程環(huán)境愈加豐富，也具備了適當(dāng)?shù)臄?shù)值數(shù)組對(duì)象，開發(fā)全?？茖W(xué)軟件的時(shí)機(jī)成熟了。2001 年，Eric Jones 和 Travis Vaught 創(chuàng)建了 Enthought 科學(xué)計(jì)算解決方案。

之后，為了簡(jiǎn)化工具堆棧，他們創(chuàng)建了以 SciPy 庫(kù)為中心的 SciPy 項(xiàng)目。該項(xiàng)目的發(fā)展勢(shì)頭很猛，2001 年 2 月推出網(wǎng)站和代碼庫(kù)，6 月宣布郵件列表，8 月推出了 SciPy 0.1 版。

SciPy 早期版本的文檔較少，但隨著 2006 年發(fā)布 Numpy 指南（Guide to Numpy），這種情況開始改變。2007 年，Sphinx 文檔生成器使得 SciPy 能夠從包含 Python 代碼的純文本中自動(dòng)呈現(xiàn)超文本和 PDF 文檔。2008 年，Pydocweb 工具使得 SciPy 又能夠以維基百科的方式進(jìn)行協(xié)作文檔開發(fā)。

在早期的 SciPy workshop 中，反復(fù)出現(xiàn)的一些主題反映了 SciPy 的開發(fā)狀態(tài)，它將重心放在了底層數(shù)組包、繪圖、并行處理、加速/包裝和用戶界面上。到了 2004 年，關(guān)于 SciPy 應(yīng)用于科學(xué)計(jì)算問(wèn)題上的內(nèi)容開始出現(xiàn)。

圖 1：自 2001 年發(fā)布 0.1 版到 2017 年推出 1.0 版本，SciPy 發(fā)展過(guò)程中的一些里程碑式事件。

SciPy 近三年的關(guān)鍵技術(shù)

從 2001 年發(fā)布的 0.1，到近幾年發(fā)布成熟版的 SciPy 1.0，最近三年 SciPy 在科學(xué)計(jì)算上有了更多的技術(shù)累積。我們可以用更少的算力運(yùn)行更大的矩陣計(jì)算，用更精簡(jiǎn)的方式擬合更復(fù)雜與多樣的概率分布，也可以跑一跑最新的最優(yōu)化方法。研究者在這篇論文中著重介紹了 SciPy 一路走來(lái)的關(guān)鍵技術(shù)。

數(shù)據(jù)結(jié)構(gòu)：稀疏矩陣

scipy.sparse 提供了 7 種稀疏矩陣數(shù)據(jù)結(jié)構(gòu)，或者稱之為稀疏格式。其中最重要的一種是壓縮行/壓縮列的稀疏格式，它們分別為 CSR 與 CSC。這兩種方法都提供了快速的主軸索引與快速的矩陣-向量乘法，這兩種稀疏格式在 SciPy 及依賴的庫(kù)中得到了廣泛的應(yīng)用。

從新特性的角度來(lái)看，scipy.sparse 矩陣與線性運(yùn)算子現(xiàn)在都已經(jīng)支持 Python 矩陣乘法（@）。

cKDTree

scipy.spatial.ckdtree 模塊實(shí)現(xiàn)了空間分割的數(shù)據(jù)結(jié)構(gòu)，該結(jié)構(gòu)會(huì)在 K 維空間中組織數(shù)據(jù)點(diǎn)。整個(gè) cKDTree 模塊通過(guò)模板化類用 C 重寫了，并新增對(duì)周期性邊界條件的支持，它經(jīng)常用于物理過(guò)程的模擬。

2013 年，基于 cKDTree.query 的 K 最近鄰算法時(shí)間復(fù)雜度逼近了對(duì)數(shù)線性。2015 年，cKDTree 二元樹計(jì)數(shù)算法通過(guò)加強(qiáng)以支持加權(quán)，這對(duì)于很多科學(xué)應(yīng)用來(lái)說(shuō)都是非常重要的，例如計(jì)算星系的相關(guān)性函數(shù)。

統(tǒng)一捆綁到已編譯代碼：LowLevelCallable

到了 SciPy 0.19，用戶就可以直接使用 scipy.LowLevelCallable 對(duì)象包裝底層函數(shù)，從而減少直接從 Python 調(diào)用已編譯 C 函數(shù)的開銷，這種編譯的 C 函數(shù)可能是由 Numba 或 Cython 生成的。

數(shù)學(xué)優(yōu)化

scipy.optimize 子包提供了數(shù)學(xué)解決方案，用于解決多種類型的「root finding」和優(yōu)化問(wèn)題。

研究者在表 1 中詳細(xì)比較了所有最小化方法的特征，這些特征說(shuō)明了 SciPy 如果要達(dá)到比較完整的水平，它需要涵蓋的數(shù)值方法或主題。

統(tǒng)計(jì)分布

scipy.status 包含了 100 多個(gè)概率分布：96 個(gè)連續(xù)分布和 13 個(gè)離散單變量分布，以及 10 個(gè)多變量分布。該實(shí)現(xiàn)依賴于一個(gè)一致的框架，該框架提供了抽樣隨機(jī)變量的方法，用以評(píng)估累積分布函數(shù)指數(shù)（CDF）和概率密度函數(shù)指數(shù)（PDF），并適合每一個(gè)分布的參數(shù)。

測(cè)試套件

在更改代碼時(shí)，測(cè)試驅(qū)動(dòng)的開發(fā)被認(rèn)為是一種管理不確定性的方法。對(duì)于 SciPy 的每個(gè)組件，研究者都編寫了多種能夠驗(yàn)證它們預(yù)期行為的可執(zhí)行小測(cè)試。這些可執(zhí)行小測(cè)試的集合被稱為『測(cè)試套件』，增強(qiáng)了對(duì)正確性和準(zhǔn)確度的置信度，并允許用戶在修改已有代碼時(shí)不會(huì)改變預(yù)期行為。

圖 2：不同版本 SciPy 中的 Python 和編譯代碼量。

除了保證單元測(cè)試通過(guò)之外，重要的是確保 SciPy 代碼庫(kù)的性能能夠隨時(shí)間推移而提升。比如，下圖 3 展示了在大約 9 年的項(xiàng)目發(fā)展歷程中，scipy.spatial. cKDTree.query 的性能提升情況。

圖 3：從 cKDTree 的提出到 SciPy 1.0 的發(fā)布，scipy.spatial.cKDTree.query 的基準(zhǔn)測(cè)試的結(jié)果。圖中的每個(gè)標(biāo)記表示 SciPy 主分支中提交的基準(zhǔn)測(cè)試的執(zhí)行時(shí)間。

SciPy 仍在路上

SciPy 項(xiàng)目每 6 個(gè)月進(jìn)行一次更新。任何感興趣、有技術(shù)能力的開發(fā)者都可以參與貢獻(xiàn)代碼。雖然 SciPy 的研發(fā)成本已經(jīng)超過(guò)了 1 千萬(wàn)美元，但是項(xiàng)目依然是沒(méi)有資金支持的。這些代碼都是由大學(xué)研究生、學(xué)術(shù)界和工業(yè)界的人們?cè)陂e暇時(shí)間完成的。

SciPy 有著一個(gè)很大的開發(fā)社區(qū)，用戶基數(shù)也很龐大。在 2017 年，通過(guò) PyPI 下載的次數(shù)是 13,096,468 次，而通過(guò) conda 下載的次數(shù)則有 5,776,017 次。

盡管如此，SciPy 依然在繼續(xù)進(jìn)步。下圖的表格是一個(gè)持續(xù)更新的文檔，描述了團(tuán)隊(duì)正在項(xiàng)目中進(jìn)行改進(jìn)和提升的工作。這份文檔也提到了一些需要改進(jìn)的地方。