Python-Numpy学习

本文是《用Python做科学计算》的学习笔记，并且参考《Python for Data Analysis》，内容稍作整合，代码均在Pycharm中运行通过。

Numpy

NumPy 是 Numerical Python 的简称，是高性能计算和数据分析的基础包。下面是它所能做的一些事情：

• ndarray，快速和节省空间的多维数组，提供数组化的算术运算和高级的 广播 功能。
• ufunc，使用标准数学函数对整个数组的数据进行快速运算，而不需要编写循环。
• 读取/写入磁盘上的阵列数据和操作存储器映像文件的工具。
• 线性代数，随机数生成，和傅里叶变换的能力。
• 集成C，C++，Fortran代码的工具。

  为什么要Numpy

  Python中用列表(list)保存一组值，可以用来当作数组使用，不过由于列表的元素可以是任何对象，因此列表中所保存的是对象的指针。这样为了保存一个简单的[1,2,3]，需要有3个指针和三个整数对象。对于数值运算来说这种结构显然比较浪费内存和CPU计算时间。

此外Python还提供了一个array模块，array对象和列表不同，它直接保存数值，和C语言的一维数组比较类似。但是由于它不支持多维，也没有各种运算函数，因此也不适合做数值运算。

NumPy的诞生弥补了这些不足，提供了两种基本的对象：

• ndarray（N-dimensional array object，下文统一称之为数组）是存储单一数据类型的多维数组
• ufunc（universal function object），是能够对数组进行处理的函数

模块导入

导入NumPy函数库:

import numpy as np

ndarray对象

创建

从python序列创建

通过np.array()，可选参数：

• dtype：指定元素类型，np.float/np.complex/np.
• order：内存顺序，后续会讲

  a = np.array([1, 2, 3, 4])
  b = np.array((1, 2, 3, 4), dtype=np.float)
  c = np.array([[1, 2, 3, 4], [4, 5, 6, 7], [7, 8, 9, 10]])
  print(a, '\n', b, '\n', c)

数组大小可以通过shape属性获得，shape相关：

• n维数组会有n个轴，其shape会有n个元素。一维数组的shape会有一个元素，二位数组的shape会有两个元素。
• 通过改变shape可以改变每个轴的长度，数组元素在内存中的位置不变
• 当某个轴的元素为-1时，将根据数组元素的个数自动计算此轴的长度
• reshape方法，可以创建一个改变了尺寸的新数组，原数组的shape保持不变。并且新旧数组共享内存区域，修改其中一个会改变另一个。

  print(a.dtype, b.dtype, a.shape, c.shape)
  c.shape = 4, 3
  print(c)
  c.shape = 2, -1
  print(c)
  d = a.reshape((2, 2))
  print(d)
  a[1] = 100  # 将数组a的第一个元素改为100
  print(d)    # 注意数组d中的2也被改变了

由numpy函数创建

先创建一个Python序列，然后通过array函数将其转换为数组，这样做显然效率不高。因此NumPy提供了很多专门用来创建数组的函数：

• 数值构造型
  • arange函数：类似于python的range函数，通过指定开始值、终值和步长来创建一维数组，注意数组不包括终值
  • linspace函数：通过指定开始值、终值和元素个数来创建一维数组，可以通过endpoint关键字指定是否包括终值，缺省设置是包括终值
  • logspace函数和linspace类似，不过它创建等比数列，下面的例子产生1(10^0)到100(10^2)、有20个元素的等比数列
  • ones, ones_like 根据提供的shape和dtype产生一个全1的数组。ones_like使用另一歌数组为入参，产生一个shape和dtype都相同的数组。
  • zeros, zeros_like 同ones和ones_like，但是生成全0的数组
  • empty, enpty_like 通过分配新内存来构造新的数组，但不同与ones 和 zeros，不初始任何值。
  • eye, identity 生成一个NxN的单位方阵（对角线上为1，其它地方为0）
• 字节序列型
  • frombuffer
  • fromfile
  • fromstring：Python的字符串实际上是字节序列，每个字符占一个字节。因此如果从字符串s创建一个8bit的整数数组的话，所得到的数组正好就是字符串中每个字符的ASCII编码。如果从字符串s创建16bit的整数数组，那么两个相邻的字节就表示一个整数
• 函数构造型：np.fromfunction，参数：
  • function:构造函数，由元素的index构造元素，其中二维：a(i,j)=fun(i,j)
  • shape:生成的shape

    # 数值型
    print(np.arange(0, 1, 0.1))
    print(np.linspace(0, 1, 12))
    print(np.logspace(0, 2, 20))
    # 字符型
    s = "abcdefgh"
    print(np.fromstring(s, dtype=np.int8))
    print(np.fromstring(s, dtype=np.int16))
    # 函数型
    # 构造一维数组，a(i)=i%4+1
    def func(i):
        return i % 4 + 1
    print(np.fromfunction(func, (10,)))
    # 构造二位数组，a（i,j）=(i+1)*(j+1)
    def func2(i, j):
        return (i+1) * (j+1)
    print(np.fromfunction(func2, (9, 9)))

数据类型

Dtypes是使NumPy如此强大和灵活的一部分。在大多数情况下，它们直接映射到底层的机器表示，这是的很容易地读取和写入二进制流到磁盘上，也能链接低级语言，如C 或Fortran编写的代码

类型	类型码	描述
int8, uint8	i1, u1	有符号和无符号8位（1字节）整数类型
int16, uint16	i2, u2	有符号和无符号16位整数类型
int32, uint32	i4, u4	有符号和无符号32位整数类型
int64, uint64	i8, u8	有符号和无符号64位整数类型
float16	f2	半精度浮点类型
float32	f4 or f	标准精度浮点。与C的 float 兼容
float64, float128	f8 or d	标准双精度浮点。与C的 double 和Python 的 folat 对象兼容
float128	f16 or g	扩展精度浮点
complex64, complex128, complex256	c8, c16, c32	分别使用两个32，64，128位浮点表示的复数
bool	?	布尔值，存储 True 和 False
object	O	Python对象类型
string_	S	定长字符窜类型（每字符一字节）。例如，为了生成长度为10的字符窜，使用 ‘S10’
unicode_	f16 or g	扩展精度浮点（字节书依赖平台）。同 string_ 有相同的语义规范（例如：U10 ）

存取元素

通过下标(切片)

数组元素的存取方法和Python的标准方法相同，可以使用切片。和Python的列表序列不同，通过下标范围获取的新的数组是原始数组的一个视图。它与原始数组共享同一块数据空间。

b = a[3:7]  # 通过下标范围产生一个新的数组b，b和a共享同一块数据空间
print(b)
b[2] = -10  # 将b的第2个元素修改为-10
print(a)  # a的第5个元素也被修改为10

通过整数序列

当使用整数序列对数组元素进行存取时，将使用整数序列中的每个元素作为下标，整数序列可以是列表或者数组。使用整数序列作为下标获得的数组不和原始数组共享数据空间。

x = np.arange(10, 1, -1)
print(x)
print(x[[3, 3, 1, 8]])  # 获取x中的下标为3, 3, 1, 8的4个元素，组成一个新的数组
b = x[np.array([3, 3, -3, 8])]  # 下标可以是负数
b[2] = 100
print(x)  # 由于b和x不共享数据空间，因此x中的值并没有改变
x[[3, 5, 1]] = -1, -2, -3  # 整数序列下标也可以用来修改元素的值
print(x)

使用布尔数组

当使用布尔数组b作为下标存取数组x中的元素时，将收集数组x中所有在数组b中对应下标为True的元素。使用布尔数组作为下标获得的数组不和原始数组共享数据空间，注意这种方式只对应于布尔数组，不能使用布尔列表。

x = np.arange(5, 0, -1)
print(x)
print(x[np.array([True, False, True, False, False])])  # 布尔数组中下标为0，2的元素为True，因此获取x中下标为0,2的元素
print(x[[True, False, True, False, False]])  # 如果是布尔列表，则把True当作1, False当作0，按照整数序列方式获取x中的元素
print(x[np.array([True, False, True, True])])  # 布尔数组的长度不够时，不够的部分都当作False
x[np.array([True, False, True, True])] = -1, -2, -3  # 布尔数组下标也可以用来修改元素
print(x)

多维数组

多维数组的存取和一维数组类似，因为多维数组有多个轴，因此它的下标需要用多个值来表示，NumPy采用组元(tuple)作为数组的下标。

下标存取（切片）

a = np.arange(0, 60, 10).reshape(-1, 1) + np.arange(0, 6)
print(a)
print(a[0, 3:5])
print(a[4:, 4:])
print(a[:, 2])
print(a[2::2, ::2])

整数序列、布尔数组

• a[(0,1,2,3,4),(1,2,3,4,5)] : 用于存取数组的下标和仍然是一个有两个元素的组元，组元中的每个元素都是整数序列，分别对应数组的第0轴和第1轴。从两个序列的对应位置取出两个整数组成下标： a[0,1], a[1,2], …, a[4,5]。
• a[3:, [0, 2, 5]] : 下标中的第0轴是一个范围，它选取第3行之后的所有行；第1轴是整数序列，它选取第0, 2, 5三列。
• a[mask, 2] : 下标的第0轴是一个布尔数组，它选取第0，2，5行；第1轴是一个整数，选取第2列。

  print(a[(0, 1, 2, 3, 4), (1, 2, 3, 4, 5)])   
  print(a[3:, [0, 2, 5]])
  mask = np.array([1, 0, 1, 0, 0, 1], dtype=np.bool)
  print(a[mask, 2])

结构数组

类似于C语言的结构体struct

• name：定义字段名
• formats：定义每个字段的类型
  • S32 : 32个字节的字符串类型，由于结构中的每个元素的大小必须固定，因此需要指定字符串的长度
  • i : 32bit的整数类型，相当于np.int32
  • f : 32bit的单精度浮点数类型，相当于np.float32

    persontype = np.dtype({
        'names': ['name', 'age', 'weight'],
        'formats': ['S32', 'i', 'f']})
    a = np.array([("Zhang", 32, 75.5), ("Wang", 24, 65.2)], dtype=persontype)
    print(a)

通过dtype获取a的类型,类型描述前面为我们添加了 ‘|’, ‘<’ 等字符，这些字符用来描述字段值的字节顺序:

• ‘|’: 忽视字节顺序
• ‘<’: 低位字节在前
• ‘>’: 高位字节在前

  print(a.dtype)
  print(a[0].dtype)

通过a构造出来的子数组，它和数组a共享内存数据，因此可以通过修改它的字段，改变原始数组中的对应字段

c = a[1]
c["name"] = "Li"
print(a[1]["name"])
b = a[:]["age"]
print(b)
b[0] = 40
print(a[0]["age"])

内存对齐

C语言的结构体为了内存寻址方便，会自动的添加一些填充用的字节，这叫做内存对齐。例如如果把下面的name[32]改为name[30]的话，由于内存对齐问题，在name和age中间会填补两个字节，最终的结构体大小不会改变。因此如果numpy中的所配置的内存大小不符合C语言的对齐规范的话，将会出现数据错位。为了解决这个问题，在创建dtype对象时，可以传递参数align=True，这样numpy的结构数组的内存对齐和C语言的结构体就一致了。

内存结构

关于数组的描述信息保存在一个数据结构中，这个结构引用两个对象：

• 一个用于描述元素类型的dtype对象：知道如何将元素的二进制数据转换为可用的值

• 一块用于保存数据的存储区域：保存着数组中所有元素的二进制数据
  

  a = np.array([[0, 1, 2], [3, 4, 5], [6, 7, 8]], dtype=np.float32)
  print(a.dtype)
  print(a.ndim)
  print(a.shape)
  print(a.strides)
  print(a.data)

• strides：保存的是当每个轴的下标增加1时，数据存储区中的指针所增加的字节数。例如图中的strides为12,4，即第0轴的下标增加1时，数据的地址增加12个字节：即a[1,0]的地址比a[0,0]的地址要高12个字节，正好是3个单精度浮点数的总字节数；第1轴下标增加1时，数据的地址增加4个字节，正好是单精度浮点数的字节数。

  b = a[::2, ::2]
  print(b)
  print(b.strides)

由于数组b和数组a共享数据存储区，而b中的第0轴和第1轴都是数组a中隔一个元素取一个，因此数组b的strides变成了24,8，正好都是数组a的两倍。 对照前面的图很容易看出数据0和2的地址相差8个字节，而0和6的地址相差24个字节。

内存顺序

元素在数据存储区中的排列格式有两种，创建时通过参数order指定：

• order=”C”（缺省值），在C语言中，多维数组的第0轴是最上位的，即第0轴的下标增加1时，元素的地址增加的字节数最多
• order=”F”，Fortan语言的多维数组的第0轴是最下位的，即第0轴的下标增加1时，地址只增加一个元素的字节数。

  c = np.array([[0, 1, 2], [3, 4, 5], [6, 7, 8]], dtype=np.float32, order="F")
  print(c.strides)

以一个二维数组a[2][2]为例

• 在C语言中，其在内存中存储为a[0][0] a[0][1] a[1][0] a[1][1]
• 而在Fortran语言中，其顺序为a[0][0] a[1][0] a[0][1] a[1][1]

ufunc运算

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

单输入函数

sin

np.sin对x中的每个元素求正弦值，然后将结果返回。可选参数out：结果输出到out，如果要覆盖x，则将out赋值为x

x = np.linspace(0, 2*np.pi, 10)
y = np.sin(x)
print(y)
t = np.sin(x, x)
print(id(t) == id(x))

与Python标准库速度对比

计算100万次正弦值，所用时间的对比：

• np.sin处理数组(C级别循环)：0.011924890681391243
• python级别循环 + math.sin处理单个值：0.41831748538827573
• python级别循环 + np.sin处理单个值：1.4091675317012233

  import time
  import math
  import numpy as np
  # python级别循环 + math.sin
  x = [i * 0.001 for i in xrange(1000000)]
  start = time.clock()
  for i, t in enumerate(x):
      x[i] = math.sin(t)
  print "math.sin:", time.clock() - start
  # np.sin(C级别循环)
  x = [i * 0.001 for i in xrange(1000000)]
  x = np.array(x)
  start = time.clock()
  np.sin(x,x)
  print "numpy.sin:", time.clock() - start
  # python级别循环 + np.sin处理单个值
  x = [i * 0.001 for i in range(1000000)]
  start = time.clock()
  for i, t in enumerate(x):
      x[i] = np.sin(t)
  print("numpy.sin loop:", time.clock() - start)

结论：numpy.sin处理数组比math.sin+python循环处理数组快几十倍。这得利于numpy.sin在C语言级别的循环计算

但是numpy.sin处理单个值的计算速度比math.sin要慢。这是因为numpy.sin为了同时支持数组和单个值的计算，其C语言的内部实现要比math.sin复杂很多。

多输入函数

add函数

add函数返回一个新的数组，此数组的每个元素都为两个参数数组的对应元素之和。它接受第3个参数指定计算结果所要写入的数组，如果指定的话，add函数就不再产生新的数组。

a = np.arange(0, 4)
b = np.arange(1, 5)
print(a)
print(b)
print(np.add(a, b))   #  a+b
print(np.add(a, b, a))   #  a+=b
print(a)

其他函数

• 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])

优化

数组对象支持这些操作符，极大地简化了算式的编写，不过要注意如果你的算式很复杂，并且要运算的数组很大的话，会因为产生大量的中间结果而降低程序的运算效率。例如：假设a b c三个数组采用算式x=a*b+c计算，那么它相当于:

t = a * b
x = t + c
del t

也就是说需要产生一个数组t保存乘法的计算结果，然后再产生最后的结果数组x。我们可以通过手工将一个算式分解为x = a*b; x += c，以减少一次内存分配。

a = np.arange(0, 100000000)
b = np.arange(0, 100000000)
start = time.clock()
x = a*b + a
print("x = a*b+c time:", time.clock() - start)
start = time.clock()
y = a*b
y += a
print("y = a*b,y += c time:", time.clock() - start)

输出

x = a*b+c time: 0.6211750219831841
y = a*b,y += c time: 0.48815467969828963

frompyfunc

可以用frompyfunc函数将一个计算单个元素的函数转换成ufunc函数。这样就可以方便地用所产生的ufunc函数对数组进行计算了

frompyfunc的调用格式为frompyfunc(func, nin, nout)，其中func是计算单个元素的函数，nin是此函数的输入参数的个数，nout是此函数的返回值的个数。

我们很容易根据上图所示写出如下的计算三角波某点y坐标的函数，有三种方式:

• 可以用下面的方法先使用列表包容(List comprehension)，计算出一个list，然后用array函数将列表转换为数组
• 用一个lambda函数对triangle_wave的参数进行一次包装。这样传入frompyfunc的函数就只有一个参数了
• 通过函数triangle_func包装三角波的三个参数，在其内部定义一个计算三角波的函数trifunc，trifunc函数在调用时会采用triangle_func的参数进行计算。最后triangle_func返回用frompyfunc转换结果。

值得注意的是用frompyfunc得到的函数计算出的数组元素的类型为object，因为frompyfunc函数无法保证Python函数返回的数据类型都完全一致。因此还需要再次 y2.astype(np.float64)将其转换为双精度浮点数组。

def triangle_wave(x, c, c0, hc):
    x = x - int(x) # 三角波的周期为1，因此只取x坐标的小数部分进行计算
    if x >= c: r = 0.0
    elif x < c0: r = x / c0 * hc
    else: r = (c-x) / (c-c0) * hc
    return r

x = np.linspace(0, 2, 1000)
y = np.array([triangle_wave(t, 0.6, 0.4, 1.0) for t in x])

triangle_ufunc = np.frompyfunc( lambda x: triangle_wave(x, 0.6, 0.4, 1.0), 1, 1)
y2 = triangle_ufunc(x)
    
def triangle_func(c, c0, hc):
    def trifunc(x):
        x = x - int(x) # 三角波的周期为1，因此只取x坐标的小数部分进行计算
        if x >= c: r = 0.0
        elif x < c0: r = x / c0 * hc
        else: r = (c-x) / (c-c0) * hc
        return r
    # 用trifunc函数创建一个ufunc函数，可以直接对数组进行计算, 不过通过此函数
    # 计算得到的是一个Object数组，需要进行类型转换
    return np.frompyfunc(trifunc, 1, 1)
y2 = triangle_func(0.6, 0.4, 1.0)(x)

广播

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

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

举例

之前：

• 先创建一个二维数组a，其shape为(6,1)，再创建一维数组b，其shape为(5,)。
• 由于a和b的shape长度(也就是ndim属性)不同，根据规则1，需要让b的shape向a对齐，于是将b的shape前面加1，补齐为(1,5)
• 这样加法运算的两个输入数组的shape分别为(6,1)和(1,5)，根据规则2，可知输出数组的shape为(6,5)
• 由于b的第0轴上的长度为1，而a的第0轴上的长度为6，因此为了让它们在第0轴上能够相加，需要将b在第0轴上的长度扩展为6，相当于b向下复制5行。同理，a向右复制4列。

  a = np.arange(0, 60, 10).reshape(-1, 1)    # [[ 0], [10], [20], [30], [40], [50]]
  b = np.arange(0, 5)    # [0 1 2 3 4]
  print(a.shape, b.shape)    # 输出：(6, 1) (5,)
  print(a + b)
  
  b.shape = 1, 5    # [[0 1 2 3 4]]
  b = b.repeat(6, axis=0)
  a = a.repeat(5, axis=1)
  print(a, b, sep='\n')
  print(a + b)

ogrid对象

当然，numpy在执行a+b运算时，其内部并不会真正将长度为1的轴用repeat函数进行扩展，如果这样做的话就太浪费空间了。numpy提供了一个快速产生如上面a,b数组的方法： ogrid对象

x, y = np.ogrid[0:1:4j, 0:1:3j]
print(x, y)

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

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

ufunc的方法

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

reduce

方法和Python的reduce函数类似，它沿着axis轴对array进行操作，相当于将运算符插入到沿axis轴的所有子数组或者元素当中

a = np.add.reduce([1, 2, 3])  # 1 + 2 + 3
b = np.add.reduce([[1, 2, 3], [4, 5, 6]], axis=0)  # 1+2+3 , 4+5+6
c = np.add.reduce([[1, 2, 3], [4, 5, 6]], axis=1)  # 1+4 , 2+5 ,3+6
print(a, b, c, sep='\n')

accumulate

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

a = np.add.accumulate([1, 2, 3])  # [1 3 6]
b = np.add.accumulate([[1, 2, 3], [4, 5, 6]], axis=0)  # [[1 2 3],[5 7 9]]
c = np.add.accumulate([[1, 2, 3], [4, 5, 6]], axis=1)  # [[ 1  3  6],[ 4  9 15]]
print(a, b, c, sep='\n')

reduceat

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

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

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

a.shape += (1,)*b.ndim
<op>(a,b)
a = a.squeeze()  # 剔除数组a中长度为1的轴

通过outer方法计算的结果是乘法表。如果将这两个数组按照等同程序一步一步的计算的话，就会发现乘法表最终是通过广播的方式计算出来的。

x = np.multiply.outer([1, 2, 3, 4, 5], [2, 3, 4])
print(x)   # [[ 2  3  4],[ 4  6  8],[ 6  9 12],[ 8 12 16],[10 15 20]]

矩阵运算

NumPy和Matlab不一样，对于多维数组的运算，缺省情况下并不使用矩阵运算，如果你希望对数组进行矩阵运算的话，可以调用相应的函数

matrix对象

numpy库提供了matrix类，使用matrix类创建的是矩阵对象，它们的加减乘除运算缺省采用矩阵方式计算，因此用法和matlab十分类似。但是由于NumPy中同时存在ndarray和matrix对象，因此用户很容易将两者弄混。这有违Python的“显式优于隐式”的原则，因此并不推荐在较复杂的程序中使用matrix。

a = np.matrix([[1, 2, 3], [5, 5, 6], [7, 9, 9]])
print(a*a**-1)

计算的是矩阵a和其逆矩阵的乘积，结果是一个单位矩阵

函数

dot

对于两个一维的数组，计算的是这两个数组对应下标元素的乘积和(数学上称之为内积)；对于二维数组，计算的是两个数组的矩阵乘积；对于多维数组，它的通用计算公式如下，即结果数组中的每个元素都是：数组a的最后一维上的所有元素与数组b的倒数第二位上的所有元素的乘积和

inner

和dot乘积一样，对于两个一维数组，计算的是这两个数组对应下标元素的乘积和；对于多维数组，它计算的结果数组中的每个元素都是：数组a和b的最后一维的内积，因此数组a和b的最后一维的长度必须相同

outer

只按照一维数组进行计算，如果传入参数是多维数组，则先将此数组展平为一维数组之后再进行运算。outer乘积计算的列向量和行向量的矩阵乘积

其他

矩阵中更高级的一些运算可以在NumPy的线性代数子库linalg中找到

文件存取

NumPy提供了多种文件操作函数方便我们存取数组内容。文件存取的格式分为：

• 文本
• 二进制
  • NumPy专用的格式化二进制类型
  • 无格式类型。

tofile && fromfile

使用数组的方法函数tofile可以方便地将数组中数据以二进制的格式写进文件。tofile输出的数据没有格式，因此用numpy.fromfile需要在读入的时候设置正确的dtype和shape才能保证数据一致：

a = np.arange(0, 12)
a.shape = 3, 4
print(a)
a.tofile("a.bin")
b = np.fromfile("a.bin", dtype=np.float)  # 按照float类型读入数据
print(b)  # 读入的数据是错误的
b = np.fromfile("a.bin", dtype=np.int32)  # 按照int32类型读入数据
print(b)  # 数据是一维的
b.shape = 3, 4  # 按照a的shape修改b的shape
print(b)  # 这次终于正确了

numpy.load && numpy.save

以NumPy专用的二进制类型保存数据，这两个函数会自动处理元素类型和shape等信息，使用它们读写数组就方便多了，但是numpy.save输出的文件很难和其它语言编写的程序读入：

np.save("a.npy", a)
c = np.load("a.npy")
print(c)

如果你想将多个数组保存到一个文件中的话，可以使用numpy.savez函数。savez函数的第一个参数是文件名，其后的参数都是需要保存的数组，也可以使用关键字参数为数组起一个名字，非关键字参数传递的数组会自动起名为arr_0, arr_1, …。

a = np.array([[1, 2, 3], [4, 5, 6]])
b = np.arange(0, 1.0, 0.1)
c = np.sin(b)
np.savez("result.npz", a, b, sin_array=c)
r = np.load("result.npz")
print(r["arr_0"])  # 数组a
print(r["arr_1"])  # 数组b
print(r["sin_array"])  # 数组c

如果你用解压软件打开result.npz文件的话，会发现其中有三个文件：arr_0.npy， arr_1.npy， sin_array.npy，其中分别保存着数组a, b, c的内容。

numpy.savetxt && numpy.loadtxt

numpy.savetxt和numpy.loadtxt可以读写1维和2维的数组

a = np.arange(0, 12, 0.5).reshape(4, -1)
np.savetxt("a.txt", a)  # 缺省按照'%.18e'格式保存数据，以空格分隔
b = np.loadtxt("a.txt")
np.savetxt("a.txt", a, fmt="%d", delimiter=",")  # 改为保存为整数，以逗号分隔
c = np.loadtxt("a.txt", delimiter=",")  # 读入的时候也需要指定逗号分隔
print(a, b, c, sep='\n')

参考

在线电子书：用Python做科学计算——NumPy-快速处理数据