SciPy入门指南(译)

翻译scipy之getting started

SciPy入门指南

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

  • toc {:toc }
p_large_LzVl_40d700008a811261.jpg
Figure 1: optimization

什么是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[i] + b[i])

这个循环将在几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有许多用来创建交互应用的高级的模块(例如TraitsUIwxPython)。和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()
p_large_34di_40d700008a7f1261.jpg
Figure 2: fft

这将出现一个显示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()
p_large_fk20_40d700008a801261.jpg
Figure 3: scipy-2

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

运行脚本

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

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


Footnotes:

1

是这样?还是我理解错了?