python用于数学计算的工具介绍：scipy和numpy

夜勿眠

SciPy入门指南


这个教程是为了帮助初学者掌握scipy并且肯能快地实际使用。

• 什么是scipy、numpy、matplotlib？
  
• 它们是用来干什么的？
  
• 如何使用scipy工作
  
• 学习使用scipy
  
• 
  示例会话
  
  
  
  • 交互工作
    
  • 运行脚本
    
  • 一些关于导入(import)的笔记
    
  
  

什么是scipy、numpy、matplotlib？


Python是一种通用语言。它被解释运行，是动态类型语言，并且非常适合交互工作和快速实现原型，然而又足够强大用来写大型应用。

NumPy是一个定义了数值数组和矩阵类型和它们的基本运算的语言扩展。

SciPy是另一种使用NumPy来做高等数学、信号处理、优化、统计和许多其它科学任务的语言扩展。

Matplotlib是一个帮助绘图的语言扩展。


它们是用来干什么的？


SciPy和其它这些能用来完成许多任务：

• 
  首先它们很适于用来进行严重依赖数学和数值运算的计算,它们可以原生地使用数组、矩阵和对数组和矩阵的运算，求出特征根，计算积分，解决微分方程。
  NumPy的数组类(被用来实现矩阵类的基础)是考虑速度的实现，所以存取NumPy数组比存取Python列表速度更快。此外，NumPy实现一种数组语言，以致于不需要大多数循环。例如，普通Python(C等等也是相似的):
  
    a = range(10000000)
  b = range(10000000)
  c = []
  for i in range(len(a)):
  c.append(a + b)
  
  这个循环将在几GHz的处理器上耗时5到10秒，而NumPy：
  
    import numpy as np
  a = np.arange(10000000)
  b = np.arange(10000000)
  c = a + b
  
  不仅这个更加简洁和易读，而且相比几乎是瞬间的，甚至NumPy的导入都比普通Python中的循环更快。为什么？Python是一种动态类型的解释语言，这意味着在每次循环迭代时它都必须检查运算对象a和b的类型来选则+正确的意义。(在python中+被用到许多地方，比如连接字符串、可以有不同元素类型的列表)当’+’的操作对象之一是一个NumPy数组时，NumPy的add函数将被Python自动选择，仅仅检测一次类型。然后它通过编译C函数执行”真正的”加法循环。这和普通python中的解释循环相比是非常快的。
  
  
• 
  有许多通用的或特定用途的数值代码使用了numpy和scipy。参考局部软件索引来查看部分列表。Python有许多用来创建交互应用的高级的模块(例如TraitsUI和wxPython)。和scipy一同使用这些是最快的创建科学应用的方法。
  
  
• 使用ipython使得交互工作更加简单。数据处理，数值模型探索，在运行中尝试运算可以很快地从一个想法得到结果(参考artical
  on ipython)。
  
• 
  matplotlib模块制造高质量的绘图，通过它你可以把你的数据和模型转化为展示或文章用的图像。不必在一个程序中做数值计算，保存数据，然后用另一个绘图。
  

如何使用scipy工作


Python是一门语言，它有几个用户界面。没有一个单独的程序可以开始并且给一个集成的用户体验。取而代之的是各种使用python的方法。1

最普遍的方法是使用高级python交互shell来输入命令和运行脚本。脚本可以用任何编辑器来写，例如SPE,PyScripter,甚至notepad，emacs或者vi/vim。

scipy和numpy默认都没有提供绘图函数。它们仅仅是数值工具。推荐的绘图工具包是matplotlib。

在Windows、Mac OS X和Linux下，所有这些工具都被Enthought Python发行版提供，获取更多关于安装这些的指导参考此站点安装scipy部分。


学习使用scipy


最快的使用scipy工作的方法可能就是这个交互数据分析教程

去学习Python语言，python教程可以让你迅速熟悉python语法和对象。你可以从http://docs.python.org/download.html下载这个教程。

Dave Kuhlman的numpy和scipy教程是另一个很好的介绍：http://www.rexx.com/~dkuhlman/scipy_course_01.html.

本站点上的文档和Cookbook部分提供更多学习材料。


示例会话


交互工作


让我们看看在矩形窗函数的傅利叶变换。我们将使用一个交互python shell——ipython来做这个。因为我们想要通过交互绘图呈现结果，我们将启动使用--pylab参数启动ipython，这个参数允许使用matplotlib交互。

$ ipython --pylab
Python 2.5.1 (r251:54863, May  2 2007, 16:27:44)
Type "copyright", "credits" or "license" for more information.
IPython 0.7.3 -- An enhanced Interactive Python.
?       -> Introduction to IPython's features.
%magic  -> Information about IPython's 'magic' % functions.
help    -> Python's own help system.
object? -> Details about 'object'. ?object also works, ?? prints more.
Welcome to pylab, a matplotlib-based Python environment.
For more information, type 'help(pylab)'.


ipython提供了许多方便的特性，想tab补全python函数和一个优秀的帮助系统。

In [1]: %logstart
Activating auto-logging. Current session state plus future input saved.
Filename       : ipython_log.py
Mode           : rotate
Output logging : False
Raw input log  : False
Timestamping   : False
State          : active


这个激活一个登录会话到一个文件。登录文件格式允许在以后就像一个python脚本一样被简单的执行，或编辑进一个程序。ipython也记录所有的输入输出(并把它们保存在叫In和Out的列表中)，因此你可以启动有追溯的登录。

In [2]: from scipy import *


因为numpy和scipy不是构建在python中的，你必须显示地告诉python加载它们的特性。Scipy提供numpy所以当导入scipy时导入它是不必要的。

现在开始实际的数学：

In [3]: a = zeros(1000)
In [4]: a[:100]=1


第一行正如你所期望的那样简单地创建了一个有1000个0的数组;numpy默认使这些0是双精度浮点数，但是如果我想要单精度或复数，我可以指定zeors的额外参数。第二行把一百个元素设置成-1.

然后我想要对这个数组进行傅利叶变换，scipy提供fft函数来完成这些：

b = fft(a)


为了看看b是什么样的，我将使用matplotlib库。如果你使用”–pylab”启动ipython将不需要导入matplotlib。否则你应导入它:from pylab import *但是你将没有交互功能(当你创建时自动绘图)。

In [6]: plot(abs(b))
Out[6]: [<matplotlib.lines.Line2D instance at 0xb7b9144c>]
In [7]: show()




这将出现一个显示b的图像的窗口，如果你启动ipython时使用--pylab的话show命令是不必要的。

我注意到如果我把b的0频移动到中间看起来更好。我可以通过连接b的后半部分和前半部分来实现，但是我记不清concatenate的语法了:

In [8]: concatenate?
Type:           builtin_function_or_method
Base Class:     <type 'builtin_function_or_method'>
String Form:    <built-in function concatenate>
Namespace:      Interactive
Docstring:
concatenate((a1, a2, ...), axis=0)
Join arrays together.
The tuple of sequences (a1, a2, ...) are joined along the given axis
(default is the first one) into a single numpy array.
Example:
>>> concatenate( ([0,1,2], [5,6,7]) )
array([0, 1, 2, 5, 6, 7])
In [9]: f=arange(-500,500,1)
In [10]: grid(True)
In [11]: plot(f,abs(concatenate((b[500:],b[:500]))))
Out[11]: [<matplotlib.lines.Line2D instance at 0xb360ca4c>]
In [12]: show()




这得到了我想要的图像。我可以使用交互控制上下移动和缩放图片，并且为包含在出版物中产生postscript输出(如果你想要学习更多绘图知识，建议你阅读matplotlib教程)。


运行脚本


当你一遍又一遍重复做着相同的工作，把一些命令保存在文件中并把它们作为脚本在ipython中运行将很有用。你可以使用”Ctrl-D”退出当前的ipython会话并且编辑ipython_log.py文件。当你想要执行这个文件中的指令时，你可以打开一个新的ipython会话输入命令%run
-i ipython_log.py。

当编辑一个脚本文件时，在ipython中尝试一些命令也很方便。这将允许你在保存和运行之前，对一些简单的情况逐行尝试你的脚本。


一些关于导入(import)的笔记


如果你仅仅初学scipy和与其伴侣，以下的东西对你并不那么重要，不要太操心。但是当你开发一些大型的应用时最好记住它。

对交互工作(在ipython中)和一些小型的脚本使用from scipy import *没什么。这样将会有个优点就是所有的功能在当前命名空间都是立即可用的。然而，对大型的程序和软件包来说，建议只导入你真正需要的函数和模块。让我们考虑这种情况：你(为了无论什么理由)想要比较numpy和scipy的fft函数。在你的脚本中你应这样写：

# import from module numpy.fft
from numpy.fft import fft
# import scipy's fft implementation and rename it;
# Note: `from scipy import fft` actually imports numpy.fft.fft (check with
# `scipy.fft?` in Ipython or look at .../site-packages/scipy/__init__.py)
from scipy.fftpack import fft as scipy_fft


这个优势就是，当你查看代码时，你可以显式的知道你在导入什么，代码便因此变得清晰和可读。而且，这通常比通过import *导入所有东西更快，特别是如果你是从一个像scipy一样特别大的库中导入。

然而，如果你使用许多不同的numpy函数，如果你显式导入每一个函数导入声明将变得非常长。但是你可以导入整个包来代替使用import *。

from numpy import *  # bad
from numpy import abs, concatenate, sin, pi, dot, amin, amax, asarray, cov, diag, zeros, empty, exp, eye, kaiser # very long
import numpy         # good
# use numpy.fft.fft() on array 'a'
b = numpy.fft.fft(a)


没关系，因为通常import numpy非常快。另一方面，scipy相当大(有很多子包)。因此from
scipy import *第一次导入可能非常慢(所有接下来的导入声明将会更迅速地执行，因为实际上没有再次导入)。这是为什么当你导入scipy时子包的导入默认被禁止(像scipy.fft)，这样它才能像import
numpy一样快。如果你想使用比如说scipy.fft，你不得不显式的导入它(这无论如何是个好的想法)。如果你想要一次加载所有子包，你将不得不import scipy; scipy.pkgbuild()。使用ipython的交互会话，你可以通过scipy
profile调用它(ipython -p scipy),为你读取scipy的配置文件(通常在~/.ipython/ipythonrc-scipy)和加载所有scipy。对在即时交互环境使用scipy和matplotlib绘图，你可以用像这样的命令ipython
--pylab -p scipy。




文章摘自：http://reverland.org/python/2012/08/24/scipy/

获得更多有关全面概述包结构和”pythonic”导入惯例的信息，看一看Python教程的这一部分