NumPy-快速处理数据--ufunc运算--广播--ufunc方法

本文摘自《用Python做科学计算》，版权归原作者所有。

1. NumPy-快速处理数据--ndarray对象--数组的创建和存取

2. NumPy-快速处理数据--ndarray对象--多维数组的存取、结构体数组存取、内存对齐、Numpy内存结构

接下来介绍ufunc运算、广播、ufunc方法

ufunc是universal function的缩写，它是一种能对数组的每个元素进行操作的函数。NumPy内置的许多ufunc函数都是在C语言级别实现的，因此它们的计算速度非常快。

>>> import numpy as np
>>> x = np.linspace(0, 2*np.pi, 10)#从0到2*PI的等距离的10个数
>>> x
array([ 0.        ,  0.6981317 ,  1.3962634 ,  2.0943951 ,  2.7925268 ,
        3.4906585 ,  4.1887902 ,  4.88692191,  5.58505361,  6.28318531])
>>> y = np.sin(x)#np.sin是一个ufunc函数，对x中的每个元素求正弦值，然后将结果返回，并且赋值给y
#注意：计算之后x中的值并没有改变，而是新创建了一个数组y保存结果。
>>> y
array([  0.00000000e+00,   6.42787610e-01,   9.84807753e-01,
         8.66025404e-01,   3.42020143e-01,  -3.42020143e-01,
        -8.66025404e-01,  -9.84807753e-01,  -6.42787610e-01,
        -2.44929360e-16])
>>> t = np.sin(x, x)
#sin(x, x)第二个参数表示计算结果要覆盖的数组，即对x中每个元素求sin，然后结果存放在x中。
>>> x
array([  0.00000000e+00,   6.42787610e-01,   9.84807753e-01,
         8.66025404e-01,   3.42020143e-01,  -3.42020143e-01,
        -8.66025404e-01,  -9.84807753e-01,  -6.42787610e-01,
        -2.44929360e-16])
>>> id(t) == id(x)#t和x id相同，表示共享一块内存
True

计算sin有两种途径：numpy.sin和math.sin。应当注意的是，math.sin只支持对单个数值求正弦，而numpy.sin既支持单个数值也支持数组。在计算单个数值时，math.sin速度快，而在计算数组时numpy.sin运算速度快。因此在导入numpy时采用import numpy as np的方法，不建议使用*号全部载入，否则后载入的sin会覆盖掉先载入的sin。在使用时根据情况选择使用numpy.sin还是math.sin

此外，numpy.sin返回的数的类型和math.sin返回的类型有所不同，math.sin返回的是Python的标准float类型，而numpy.sin则返回一个numpy.float64类型

>>> type(math.sin(0.5))# math.sin返回的是Python的标准float类型
<type 'float'>
>>> type(np.sin(0.5))  # numpy.sin则返回一个numpy.float64类型
<type 'numpy.float64'>

add()函数：计算两个数组相加

>>> a = np.arange(0, 4)
>>> a
array([0, 1, 2, 3])
>>> b = np.arange(1, 5)
>>> b
array([1, 2, 3, 4])
>>> np.add(a, b)   #数组a + b带返回结果
array([1, 3, 5, 7])
>>> np.add(a, b, a)#第三个参数用来存放计算结果
array([1, 3, 5, 7])
>>> a
array([1, 3, 5, 7])

由于Python的操作符重载功能，计算两个数组相加可以简单地写为a+b，而np.add(a,b,a)则可以用a+=b来表示。下面是数组的运算符和其对应的ufunc函数的一个列表，注意除号"/"的意义根据是否激活__future__.division有所不同。

y = x1 + x2	add(x1, x2 [, y])
y = x1 - x2	subtract(x1, x2 [, y])
y = x1 * x2	multiply (x1, x2 [, y])
y = x1 / x2	divide (x1, x2 [, y]), 如果两个数组的元素为整数，那么用整数除法
y = x1 / x2	true divide (x1, x2 [, y]), 总是返回精确的商
y = x1 // x2	floor divide (x1, x2 [, y]), 总是对返回值取整
y = -x	negative(x [,y])
y = x1 ** x2	power(x1, x2 [, y])
y = x1 % x2	remainder(x1, x2 [, y]), mod(x1, x2, [, y])

数组对象支持这些操作符，极大地简化了算式的编写。

广播

当我们使用ufunc函数对两个数组进行计算时，ufunc函数会对这两个数组的对应元素进行计算，因此它要求这两个数组有相同的大小(shape相同)。如果两个数组的shape不同的话，会进行如下的广播(broadcasting)处理：

1. 让所有输入数组都向其中shape最长的数组看齐，shape中不足的部分都通过在前面加1补齐
2. 输出数组的shape是输入数组shape的各个轴上的最大值
3. 如果输入数组的某个轴和输出数组的对应轴的长度相同或者其长度为1时，这个数组能够用来计算，否则出错
4. 当输入数组的某个轴的长度为1时，沿着此轴运算时都用此轴上的第一组值 

>>> a = np.arange(0, 60, 10).reshape(-1, 1)#产生0~60，间距为10的数组，行(第0轴)为-1则根据列自动匹配
>>> a
array([[ 0],
       [10],
       [20],
       [30],
       [40],
       [50]])
>>> a.shape#a是二维数组
(6, 1)
>>> b = np.arange(0, 5)
>>> b
array([0, 1, 2, 3, 4])
>>> b.shape#b是一维数组
(5,)
>>> c = a + b
#由于a和b的shape长度(也就是ndim属性)不同，根据规则1，需要让b的shape向a对齐，于是将b的shape补齐为(1,5)。
#两个输入数组的shape分别为(6,1)和(1,5)，根据规则2，可知输出数组的shape为(6,5)。
#规则3，当输入数组的某个轴的长度为1时，沿着此轴运算时都用此轴上的第一组值。
>>> c
array([[ 0,  1,  2,  3,  4],
       [10, 11, 12, 13, 14],
       [20, 21, 22, 23, 24],
       [30, 31, 32, 33, 34],
       [40, 41, 42, 43, 44],
       [50, 51, 52, 53, 54]])
>>> c.shape
(6, 5)

由于这种广播计算很常用，因此numpy提供了一个快速产生如上面a,b数组的方法： ogrid对象：

ogrid是一个很有趣的对象，它像一个多维数组一样，用切片组元作为下标进行存取，返回的是一组可以用来广播计算的数组。其切片下标有两种形式：

• 开始值:结束值:步长，和np.arange(开始值, 结束值, 步长)类似

• 开始值:结束值:长度j，当第三个参数为虚数时，它表示返回的数组的长度，和np.linspace(开始值, 结束值, 长度)类似：

>>> x,y = np.ogrid[0:5,0:5]        #开始值:结束值:步长
>>> x
array([[0],
       [1],
       [2],
       [3],
       [4]])
>>> y
array([[0, 1, 2, 3, 4]])
>>> x, y = np.ogrid[0:1:4j, 0:1:3j]#开始值:结束值:长度j
>>> x
array([[ 0.        ],
       [ 0.33333333],
       [ 0.66666667],
       [ 1.        ]])
>>> y
array([[ 0. ,  0.5,  1. ]])

注意：ogrid不是函数！！！

根据Python的语法，只有在中括号中才能使用用冒号隔开的切片语法，如果ogrid是函数的话，那么这些切片必须使用slice函数创建，这显然会增加代码的长度。

ufunc的方法

ufunc函数本身还有些方法，这些方法只对两个输入一个输出(如加法)的ufunc函数有效，其它的ufunc对象调用这些方法时会抛出ValueError异常。

reduce 方法和Python的reduce函数类似，语法为：

<op>.reduce (array=, axis=0, dtype=None)
#沿着axis轴对array进行操作，将<op>运算符插入到沿axis轴的所有子数组或者元素当中

>>> np.add.reduce([1,2,3]) # 1 + 2 + 3
6
>>> np.add.reduce([[1,2,3],[4,5,6]], axis=1) # 1,4 + 2,5 + 3,6
array([ 6, 15])

accumulate 方法和reduce方法类似，只是它返回的数组和输入的数组的shape相同，保存所有的中间计算结果：

>>> np.add.accumulate([1,2,3])
array([1, 3, 6])
>>> np.add.accumulate([[1,2,3],[4,5,6]], axis=1)
array([[ 1,  3,  6],
       [ 4,  9, 15]])

reduceat 方法计算多组reduce的结果，通过indices参数指定一系列reduce的起始和终了位置。

>>> a = np.array([1,2,3,4])
>>> result = np.add.reduceat(a,indices=[0,1,0,2,0,3,0])
>>> result
array([ 1,  2,  3,  3,  6,  4, 10])

outer 方法，<op>.outer(a,b)方法的计算等同于如下程序：

>>> np.multiply.outer([1,2,3,4,5],[2,3,4])#乘法表最终是通过广播的方式计算出来的
array([[ 2,  3,  4],
       [ 4,  6,  8],
       [ 6,  9, 12],
       [ 8, 12, 16],
       [10, 15, 20]])
#>>> np.multiply([1,2,3,4,5],[2,3,4])#ValueError