利用python进行数据分析之numpy基础(三)

利用python进行数据分析之numpy基础(三)

目录

4.1 NumPy的narray:一种多维数组对象

4.1.1 创建ndarray

4.1.2 ndarray 的数据类型

4.1.3 数组和标量之间的运算

4.1.4 基本的索引和切片

4.1.5 切片索引

4.1.6 布尔型索引

4.1.7 花式索引

4.1.8 数组转置和轴对换

4.2 通用函数:快速的元素级数组函数

4.3 利用数组进行数据处理

4.3.1 将条件逻辑表述为数组运算

4.3.2 数学和统计方法

4.3.3 用于布尔型数组的方法

4.3.4 排序

4.3.5 唯一化以及其他的集合逻辑

4.4 用于数组的文件输入和输出

4.4.1 数组以二进制格式保存到磁盘

4.4.2 存取文本文件

4.5 线性代数

4.6 随机数生成

4.7 范例:随机漫步

4.7.1 一次模拟多个随机漫步



本章主要内容:

  • numpy的adarray:一种多维数组对象
  • 通用函数:快速的元素级数组函数
  • 利用数组进行数据处理
  • 用于数组的文件输入和输出
  • 线性代数
  • 随机数生成
  • 示例:随机漫步

NumPy(Numerical Python)是高性能科学计算和数据分析的基础包。是本书几乎所有的高级工具的构建基础。部分功能如下:

  • 快速高效的多维数组对象ndarray
  • 读写硬盘上基于数组的数据集的工具
  • 线性代数运算、傅里叶变换、以及随机数的生产

大于大部分数据分析而言,最关注的几点:

  • 用于数据整理和清理、子集构造和过滤、转换等快速的矢量化数组运算。
  • 常用的数组算法,如排序、唯一化、集合运算等
  • 将条件逻辑表述为数组表达式(而不是if-else分支的循环)
  • 高效的描述统计和数据聚合/摘要运算

4.1 NumPy的narray:一种多维数组对象

NumPy 最重要的一个特点就是其N维数组对象(即ndarray),该对象是一个快速而灵活的大数据集容器。

In [15]: from numpy import *
In [16]: data=array([[0.645,-0.454,0.234],[0.234,0.765,0.456]])
In [17]: data
Out[17]:
array([[ 0.645, -0.454, 0.234],
[ 0.234, 0.765, 0.456]])
In [18]: data*9
Out[18]:
array([[ 5.805, -4.086, 2.106],
[ 2.106, 6.885, 4.104]])

ndarray是一个通用的同构数据多维容器,也就是说,其他所有元素必须是相同类型的。

每个数组都有一个shape(一个表示各维度大小的元组)和dtype(一个用于说明数组数据类型的对象)

In [19]: data.shape
Out[19]: (2L, 3L)
In [20]: data.dtype
Out[20]: dtype('float64')

4.1.1 创建ndarray

创建数组的最简单的办法就是使用array函数。,接受一切序列型的对象(包括其他数组),然后产生一个新的含有传入数据的NumPy的数组。

以一个列表的转换为例:

In [21]: data1=[6,7.5,8,0,1]
In [23]: arr1=np.array(data1)
In [24]: arr1
Out[24]: array([ 6. , 7.5, 8. , 0. , 1. ])

嵌套序列(比如一组等长列表组成的列表)将会被转换为一个多维数组

In [25]: data2=[[1,2,3,4],[5,6,7,8]]
In [26]: arr2=np.array(data2)
In [27]: arr2
Out[27]:
array([[1, 2, 3, 4],
[5, 6, 7, 8]])
In [28]: arr2.ndim
Out[28]: 2
In [29]: arr2.shape
Out[29]: (2L, 4L)

na.array 会尝试为新建的这个数组推断出一个较为合适的数据类型。

还有,zeros、ones创建指定长度或形状的全0或全1的数组。

empty可以创建一个没有任何具体值的数组。

要用以上的方法创建数组,需传入一个表示形状的元组即可。

In [31]: np.zeros(5)
Out[31]: array([ 0., 0., 0., 0., 0.])
In [32]: np.zeros((2,4))
Out[32]:
array([[ 0., 0., 0., 0.],
[ 0., 0., 0., 0.]])
In [34]: np.empty((2,3,2))
Out[34]:
array([[[ 0.00000000e+000, 6.36598737e-314],
[ 0.00000000e+000, 1.27319747e-313],
[ 1.27319747e-313, 1.27319747e-313]],
[[ 1.27319747e-313, 1.27319747e-313],
[ 0.00000000e+000, 4.44659081e-323],
[ 2.54639495e-313, 6.42285340e-323]]])

arange是python内置函数range的数组版:

In [37]: np.arange(10)
Out[37]: array([0, 1, 2, 3, 4, 5, 6, 7, 8, 9])

由于NumPy关注的是数值计算。所以一般数据类型都指的是float64(浮点数)

4.1.2 ndarray 的数据类型

dtype(数据类型)是一个特殊的对象,它含有ndarray将一块内存为特定数据类型所需要的信息:

In [5]: arr1=np.array([1,2,3],dtype=np.float64)
In [6]: arr2=np.array([1,2,3],dtype=np.int32)
In [7]: arr1.dtype
Out[7]: dtype('float64')
In [8]: arr2.dtype
Out[8]: dtype('int32')

dtype是NumPy如此强大和灵活的原因之一。

数值型dtype的命名方式相同:一个类型名(float或int),后面跟一个用于表示各元素位长的数字。

通常只需知道你所处理的数据的大量类型是浮点数、复数、整数、布尔值、字符串还是普通的python对象即可。

当你需要控制数据在内存中和磁盘中的存储方式的时候(尤其是大数据集),就需要了解如何控制和存储类型。

下表给出numpy的数据类型

可以通过ndarray的astype方法显式地转换其dtype:

下面的例子将整数转换为浮点数。

In [6]: arr=np.array([1,2,3,4,5])
In [7]: arr.dtype
Out[7]: dtype('int32')
In [8]: float_arr=arr.astype(np.float64)
In [9]: float_arr.dtype
Out[9]: dtype('float64')

下面的例子将浮点数转换为整数,小数部分被截断了。

In [10]: arr=np.array([3.7,-1.2,-2.6,0.5,12.9,10.1])
In [11]: arr
Out[11]: array([ 3.7, -1.2, -2.6, 0.5, 12.9, 10.1])
In [12]: arr.astype(np.int32)
Out[12]: array([ 3, -1, -2, 0, 12, 10])

如果某字符串数组表示全是数字,也可以用astype将其转换为数值形式。

In [13]: numeric_strings=np.array(['1.25','-9.6','42'],dtype=np.string_)
In [14]: numeric_strings.astype(float)
Out[14]: array([ 1.25, -9.6 , 42. ])

数组的dtype还有一个用法:

In [15]: int_array=np.arange(10)
In [16]: calibers=np.array([.22,.270,.357,.380,.44,.50],dtype=np.float64)
In [17]: int_array.astype(calibers.dtype)
Out[17]: array([ 0., 1., 2., 3., 4., 5., 6., 7., 8., 9.])

还可以用简洁的类型来表示dtype:

In [19]: empty_uint32=np.empty(8,dtype='u4')
In [20]: empty_uint32
Out[20]:
array([ 0, 1074790400, 0, 1075576832, 0,
3222536192, 0, 1074266112], dtype=uint32)

调用astype无论如何都会创建出一个新的数组(原始数据的一份拷贝)

4.1.2 数组和标量之间的运算

数组很重要,因为数组能使你不写循环就可以对数据执行批量运算。这通常叫做矢量化。

大小相等的数组之间的任何算术运算都会将运算应用到元素级:

In [22]: arr=np.array([[1.,2.,3.],[4.,5.,6.]])
In [23]: arr
Out[23]:
array([[ 1., 2., 3.],
[ 4., 5., 6.]])
In [24]: arr*arr
Out[24]:
array([[ 1., 4., 9.],
[ 16., 25., 36.]])
In [25]: arr-arr
Out[25]:
array([[ 0., 0., 0.],
[ 0., 0., 0.]])
In [26]: arr/arr
Out[26]:
array([[ 1., 1., 1.],
[ 1., 1., 1.]])
In [27]: arr**0.5
Out[27]:
array([[ 1. , 1.41421356, 1.73205081],
[ 2. , 2.23606798, 2.44948974]])
In [28]: 1/arr
Out[28]:
array([[ 1. , 0.5 , 0.33333333],
[ 0.25 , 0.2 , 0.16666667]])

数组与标量的算术运算也会将那个标量传播到各个元素。

不同大小的数组之间的运算叫做广播。

4.1.3 基本的索引和切片

numpy数组的索引是一个内容丰富的主题。

因为选取数据子集或单个元素的方式有很多。一维数组很简单。从表面上看,跟python列表的功能差不多:

In [29]: arr=np.arange(10)
In [30]: arr
Out[30]: array([0, 1, 2, 3, 4, 5, 6, 7, 8, 9])
In [31]: arr[5]
Out[31]: 5
In [32]: arr[5:8]
Out[32]: array([5, 6, 7])
In [33]: arr[5:8]
Out[33]: array([5, 6, 7])
In [34]: arr[5:8]=12
In [35]: arr
Out[35]: array([ 0, 1, 2, 3, 4, 12, 12, 12, 8, 9])

将一个标量赋值给一个切片时(如arr[5:8]=12),该值会自动传播到整个选区。

跟列表最重要的区别在于,数组切片是原始数组的视图。这意味着数据不会被复制,视图上任何修改都会直接反映到源数组上:

In [16]: arr_slice=arr[5:8]
In [17]: arr_slice[1]=12345
In [18]: arr
Out[18]: array([ 0, 1, 2, 3, 4, 12, 12345, 12, 8, 9])
In [19]: arr_slice[:]=64
In [20]: arr
Out[20]: array([ 0, 1, 2, 3, 4, 64, 64, 64, 8, 9])

若你想得到的是ndarray切片的一份副本而非视图,就需要显式的进行复制操作,例如:

arr[5:8].copy()

在一个二维数组中,各索引位置上的元素不再是标量而是一维数组:

In [22]: arr2d=np.array([[1,2,3],[4,5,6],[7,8,9]])
In [23]: arr2d[2]
Out[23]: array([7, 8, 9])
传入一个以逗号隔开的索引列表来选取单个元素。下面两种方式等价:
In [24]: arr2d[0][2]
Out[24]: 3
In [25]: arr2d[0,2]
Out[25]: 3

在多维数组中,若省略了后面的索引,则返回对象是一个维度低一点的ndarray(含有高一维度上的所有数据)。

在2*2*3 数组的arr3中:

In [28]: arr3d=np.array([[[1,2,3],[4,5,6]],[[7,8,9],[10,11,12]]])
In [29]: arr3d
Out[29]:
array([[[ 1, 2, 3],
[ 4, 5, 6]],
[[ 7, 8, 9],
[10, 11, 12]]])

注:括号外的“维度”是一维、二维、三维之类的意思,而括号里面理解为“轴”。

指的是返回的低维数组含有原始高维数组某条轴上的所有数据。

标量值和数组都可以被赋值给arr3d[0]

In [30]: arr3d[0]
Out[30]:
array([[1, 2, 3],
[4, 5, 6]])
In [31]: old_values=arr3d[0].copy()
In [32]: arr3d[0]=42
In [33]: arr3d
Out[33]:
array([[[42, 42, 42],
[42, 42, 42]],
[[ 7, 8, 9],
[10, 11, 12]]])
In [34]: arr3d[0]=old_values
In [35]: arr3d
Out[35]:
array([[[ 1, 2, 3],
[ 4, 5, 6]],
[[ 7, 8, 9],
[10, 11, 12]]])

以此类推,arr3d[1,0]可以访问索引以(1,0)开头的那些值:

In [36]: arr3d[1,0]
Out[36]: array([7, 8, 9])

4.1.4 切片索引

In [6]: arr
Out[6]: array([0, 1, 2, 3, 4, 5, 6, 7, 8, 9])
In [7]: arr[1:6]
Out[7]: array([1, 2, 3, 4, 5])

高维度的对象花样较多,可以在一个或多个轴上进行切片,也可以跟整数索引混合使用。

In [8]: arr2d=np.array([[1,2,3],[4,5,6],[7,8,9]])
In [9]: arr2d[:2]
Out[9]:
array([[1, 2, 3],
[4, 5, 6]])

也可以一次传入多个切片,就像传入多个索引一样:

In [10]: arr2d[:2,1:]
Out[10]:
array([[2, 3],
[5, 6]])

通过将整数的索引和切片混合,可以得到低维度的切片:

In [11]: arr2d[1,:2]
Out[11]: array([4, 5])
In [12]: arr2d[2,:1]
Out[12]: array([7])
In [13]: arr2d[:,:1]
Out[13]:
array([[1],
[4],
[7]])
In [14]: arr2d[:2,:1]
Out[14]:
array([[1],
[4]])

4.1.5 布尔型索引

我们使用numpy.random 中的randn函数生成一些正态分布的随机数据:

In [15]: names=np.array(['Bob','Joe','Will','Bob','Will','Joe','Joe'])
In [25]: data=np.random.randn(7,4)
In [26]: names
Out[26]:
array(['Bob', 'Joe', 'Will', 'Bob', 'Will', 'Joe', 'Joe'],
dtype='|S4')
In [27]: data
Out[27]:
array([[-0.70649745, 0.25280456, 0.33512776, 0.94623182],
[ 0.36423875, 0.71782951, -0.20608715, 0.50924173],
[ 0.23875182, -0.17342783, 0.22927433, -0.6337787 ],
[ 1.42096692, -0.47960948, 1.06856931, 0.13755213],
[ 0.81513134, -1.31711409, -0.05354505, 0.14838148],
[-0.61441209, -0.1091356 , 1.1222585 , -0.34686519],
[ 0.34995035, 0.3710245 , 0.9741113 , 0.42007426]])

我们假设想选出对应于名字“Bob”的所有行。数组的比较运算(如==)也是矢量化的。

因此,对names和字符串“Bob”的比较元素将会产生一个布尔型数组:

In [28]: names=='Bob'
Out[28]: array([ True, False, False, True, False, False, False], dtype=bool)

这个布尔型数组可用于数组索引:

In [29]: data[names=='Bob']
Out[29]:
array([[-0.70649745, 0.25280456, 0.33512776, 0.94623182],
[ 1.42096692, -0.47960948, 1.06856931, 0.13755213]])

布尔型数组长度必须跟被索引的轴长度一致。

可以将布尔数组跟切片、整数混合使用:

In [30]: data[names=='Bob',2:]
Out[30]:
array([[ 0.33512776, 0.94623182],
[ 1.06856931, 0.13755213]])
In [32]: data[names=='Bob',3]
Out[32]: array([ 0.94623182, 0.13755213])

要选择除了“Bob”以外的其他值,既可以使用不等于符号(!=),也可以通过负号(-)对条件进行否定

In [37]: names!='Bob'
Out[37]: array([False, True, True, False, True, True, True], dtype=bool)
In [38]: data[-(names=='Bob')]
H:/Anaconda2-4.3.0.1/Scripts/ipython-script.py:1: DeprecationWarning: numpy boolean negative, the `-` operator, is deprecated, use the `~` operator or the logical_not function instead.
if __name__ == '__main__':
Out[38]:
array([[ 0.36423875, 0.71782951, -0.20608715, 0.50924173],
[ 0.23875182, -0.17342783, 0.22927433, -0.6337787 ],
[ 0.81513134, -1.31711409, -0.05354505, 0.14838148],
[-0.61441209, -0.1091356 , 1.1222585 , -0.34686519],
[ 0.34995035, 0.3710245 , 0.9741113 , 0.42007426]])

选取这三个名字中的两个需要组合应用多个布尔条件,使用&(和)、|(或)之类的布尔算术运算符即可:

In [34]: mask=(names=='Bob')|(names=='Will')
In [35]: mask
Out[35]: array([ True, False, True, True, True, False, False], dtype=bool)
In [36]: data[mask]
Out[36]:
array([[-0.70649745, 0.25280456, 0.33512776, 0.94623182],
[ 0.23875182, -0.17342783, 0.22927433, -0.6337787 ],
[ 1.42096692, -0.47960948, 1.06856931, 0.13755213],
[ 0.81513134, -1.31711409, -0.05354505, 0.14838148]])

注意:python中的关键字and和or在布尔型数组中无效。

为了将data中的所有负值都设置为0,我们只需要:

In [39]: data[data<0]=0
In [40]: data
Out[40]:
array([[ 0. , 0.25280456, 0.33512776, 0.94623182],
[ 0.36423875, 0.71782951, 0. , 0.50924173],
[ 0.23875182, 0. , 0.22927433, 0. ],
[ 1.42096692, 0. , 1.06856931, 0.13755213],
[ 0.81513134, 0. , 0. , 0.14838148],
[ 0. , 0. , 1.1222585 , 0. ],
[ 0.34995035, 0.3710245 , 0.9741113 , 0.42007426]])

通过一维布尔数组设置整行或列的值也很简单:

In [41]: data[names!='Joe']=7
In [42]: data
Out[42]:
array([[ 7. , 7. , 7. , 7. ],
[ 0.36423875, 0.71782951, 0. , 0.50924173],
[ 7. , 7. , 7. , 7. ],
[ 7. , 7. , 7. , 7. ],
[ 7. , 7. , 7. , 7. ],
[ 0. , 0. , 1.1222585 , 0. ],
[ 0.34995035, 0.3710245 , 0.9741113 , 0.42007426]])

4.1.6 花式索引

花式索引(fancy indexing)是一个numpy术语,指的是利用整数数组进行索引。

假设有一个8*4的数组:

In [2]: arr=np.empty((8,4))
In [3]: for i in range(8):
...: arr[i]=i
...:
In [4]: arr
Out[4]:
array([[ 0., 0., 0., 0.],
[ 1., 1., 1., 1.],
[ 2., 2., 2., 2.],
[ 3., 3., 3., 3.],
[ 4., 4., 4., 4.],
[ 5., 5., 5., 5.],
[ 6., 6., 6., 6.],
[ 7., 7., 7., 7.]])

为了特定顺序选取子集,只需要传入一个用于指定顺序的整数列表或ndarray即可:

In [5]: arr[[4,3,0,6]]
Out[5]:
array([[ 4., 4., 4., 4.],
[ 3., 3., 3., 3.],
[ 0., 0., 0., 0.],
[ 6., 6., 6., 6.]])

使用负数索引将会从末尾开始选取行:

In [6]: arr[[-3,-5,-7]]
Out[6]:
array([[ 5., 5., 5., 5.],
[ 3., 3., 3., 3.],
[ 1., 1., 1., 1.]])

4.1.7 数组转置和轴对换

转置(transpose)是重塑的一种特殊方式。返回的是源数据的视图。

数组不仅有transpose方法,还有一个特殊的T属性:

In [7]: arr=np.arange(15).reshape((3,5))
In [8]: arr
Out[8]:
array([[ 0, 1, 2, 3, 4],
[ 5, 6, 7, 8, 9],
[10, 11, 12, 13, 14]])
In [9]: arr.T
Out[9]:
array([[ 0, 5, 10],
[ 1, 6, 11],
[ 2, 7, 12],
[ 3, 8, 13],
[ 4, 9, 14]])

进行矩阵计算的时候,经常需要用到该操作,

比如利用np.dot计算矩阵内积 XTX

In [10]: arr=np.random.randn(6,3)
In [11]: np.dot(arr.T,arr)
Out[11]:
array([[ 4.34145127, -0.33698457, 0.07195685],
[-0.33698457, 2.17753847, 0.34911861],
[ 0.07195685, 0.34911861, 0.86871981]])

对于高纬数组,transpose需要一个由轴编号组成的元组才能对这些轴进行转置

In [16]: arr=np.arange(16).reshape((2,2,4))
In [17]: arr
Out[17]:
array([[[ 0, 1, 2, 3],
[ 4, 5, 6, 7]],
[[ 8, 9, 10, 11],
[12, 13, 14, 15]]])
In [18]: arr.transpose((1,0,2))
Out[18]:
array([[[ 0, 1, 2, 3],
[ 8, 9, 10, 11]],
[[ 4, 5, 6, 7],
[12, 13, 14, 15]]])

简单的转置可以使用.T,其实就是进行轴对换而已。

ndarray还有一个swapaxes方法,需要接受一对轴编号:

In [19]: arr
Out[19]:
array([[[ 0, 1, 2, 3],
[ 4, 5, 6, 7]],
[[ 8, 9, 10, 11],
[12, 13, 14, 15]]])
In [20]: arr.swapaxes(1,2)
Out[20]:
array([[[ 0, 4],
[ 1, 5],
[ 2, 6],
[ 3, 7]],
[[ 8, 12],
[ 9, 13],
[10, 14],
[11, 15]]])

4.2 通用函数:快速的元素级数组函数

通用函数(即ufunc)是一种对ndarray中的数据执行元素级运算的函数,

可以将其看做简单的函数(接受一个或多个标量值,并产生一个或多个标量值)的矢量化包装器。

许多ufunc都是简单的元素级变体,如sqrt和exp:

In [21]: arr=np.arange(10)
In [22]: np.sqrt(arr)
Out[22]:
array([ 0. , 1. , 1.41421356, 1.73205081, 2. ,
2.23606798, 2.44948974, 2.64575131, 2.82842712, 3. ])
In [23]: np.exp(arr)
Out[23]:
array([ 1.00000000e+00, 2.71828183e+00, 7.38905610e+00,
2.00855369e+01, 5.45981500e+01, 1.48413159e+02,
4.03428793e+02, 1.09663316e+03, 2.98095799e+03,
8.10308393e+03])

这些都是一元(unary)ufunc

另外一些(如add或maximum)接受2个数组(因此也叫二元(binary)ufunc),并返回一个结果数组:

In [25]: x=np.random.randn(8)
In [26]: y=np.random.randn(8)
In [27]: x
Out[27]:
array([-1.47208377, 1.77425071, 1.60968406, 2.51201825, -0.42125923,
-0.1886985 , 0.22943125, -1.33311273])
In [28]: y
Out[28]:
array([-0.79881692, -0.79049443, 0.4106048 , -0.65989715, 1.19432883,
2.66407523, -1.08633441, -0.86861555])
In [30]: np.maximum(x,y) #元素级的最大值
Out[30]:
array([-0.79881692, 1.77425071, 1.60968406, 2.51201825, 1.19432883,
2.66407523, 0.22943125, -0.86861555])

通用函数(ufunc)的确可以返回多个数组。modf就是一个例子,是python内置函数divmod的矢量化版本,用于浮点数数组的小数和整数部分。

In [31]: arr=np.random.randn(7)*5
In [32]: np.modf(arr)
Out[32]:
(array([-0.6184109 , -0.69467218, 0.61799246, 0.51558266, -0.26625316,
0.25942319, 0.15090508]),
array([ -5., -0., 1., 3., -6., 11., 6.]))

一元ufunc

二元ufunc

4.3 利用数组进行数据处理

numpy数组使你可以将许多数据处理表达任务表述为简洁的数组表达式。(否则需要编写循环)。

用数组表达式代替循环的做法,通常被称之为 矢量化。

一般来说,矢量化的数组要比等价的python方式上快的多,尤其是各种数值的计算。

假如,我们想在一组值(网格型)上计算sqrt(x^2+y^2) 。

np.meshgrid函数接受两个一维数组,并产生两个二维矩阵(对应于两个数组中所有的(x,y)对)

In [33]: points=np.arange(-5,5,0.01) #1000间隔相等的点
In [34]: xs,ys=np.meshgrid(points,points)
In [35]: ys
Out[35]:
array([[-5. , -5. , -5. , ..., -5. , -5. , -5. ],
[-4.99, -4.99, -4.99, ..., -4.99, -4.99, -4.99],
[-4.98, -4.98, -4.98, ..., -4.98, -4.98, -4.98],
...,
[ 4.97, 4.97, 4.97, ..., 4.97, 4.97, 4.97],
[ 4.98, 4.98, 4.98, ..., 4.98, 4.98, 4.98],
[ 4.99, 4.99, 4.99, ..., 4.99, 4.99, 4.99]])

把这两个数组当做两个浮点数那样编写表达式即可:

In [36]: import matplotlib.pyplot as plt
In [37]: z=np.sqrt(xs**2+ys**2)
In [38]: z
Out[38]:
array([[ 7.07106781, 7.06400028, 7.05693985, ..., 7.04988652,
7.05693985, 7.06400028],
[ 7.06400028, 7.05692568, 7.04985815, ..., 7.04279774,
7.04985815, 7.05692568],
[ 7.05693985, 7.04985815, 7.04278354, ..., 7.03571603,
7.04278354, 7.04985815],
...,
[ 7.04988652, 7.04279774, 7.03571603, ..., 7.0286414 ,
7.03571603, 7.04279774],
[ 7.05693985, 7.04985815, 7.04278354, ..., 7.03571603,
7.04278354, 7.04985815],
[ 7.06400028, 7.05692568, 7.04985815, ..., 7.04279774,
7.04985815, 7.05692568]])
In [39]: plt.imshow(z,cmap=plt.cm.gray);plt.colorbar()
Out[39]: <matplotlib.colorbar.Colorbar at 0x893b208>
In [40]: plt.title("Image plot of $\sqrt{x^2+y^2}$ for a grid of values")
Out[40]: <matplotlib.text.Text at 0x87367b8>
In [41]: plt.show()

4.3.1 将条件逻辑表述为数组运算

numpy.where 函数是 三元表达式 x if condition else y 矢量化版本。

假设我们有一个布尔数组和两个值数组:

In [42]: xarr=np.array([1.1,1.2,1.3,1.4,1.5])
In [43]: yarr=np.array([2.1,2.2,2.3,2.4,2.5])
In [44]: cond=np.array([True,False,True,True,False])

假设我们想要根据cond中的值选取array 和yarr的值:当cond 中的值为True时,选取xarr的值,否则从yarr中选取。

列表推倒式的写法如下表示:

In [42]: xarr=np.array([1.1,1.2,1.3,1.4,1.5])
In [43]: yarr=np.array([2.1,2.2,2.3,2.4,2.5])
In [44]: cond=np.array([True,False,True,True,False])

假设我们想根据cond 中的值选取xarr 和 yarr的值,当cond中的值为True时,选取xarr的值,否则从yarr中选取。列表推导式的写法应该如下:

In [46]: result=[(x if c else y) for x,y,c in zip(xarr,yarr,cond)]
In [47]: result
Out[47]: [1.1000000000000001, 2.2000000000000002, 1.3, 1.3999999999999999, 2.5]

但是:第一,它对大数组的处理速度不是很快;第二,无法用于多维数组。

如果使用np.where,则可以将该功能写得非常简洁:

In [49]: result=np.where(cond,xarr,yarr)
In [50]: result
Out[50]: array([ 1.1, 2.2, 1.3, 1.4, 2.5])

np.where的第二个和第三个参数不必是数组,他们都是标量值。

where通常用于根据另一个数组而产生一个新的数组。假设有一个随机数据组成的矩阵,你希望将所有的正值替换为2,将所有的负值替换为-2.

In [51]: arr=np.random.randn(4,4)
In [52]: arr
Out[52]:
array([[-0.80514356, 1.54340728, 0.69401056, 0.62621073],
[-1.82695429, -0.95471409, 1.15567363, 0.14769847],
[ 1.17490595, 0.78044466, 1.05243103, -0.36894073],
[-0.52798532, 0.21361986, -0.33805238, -0.20382236]])
In [53]: np.where(arr>0,2,-2)
Out[53]:
array([[-2, 2, 2, 2],
[-2, -2, 2, 2],
[ 2, 2, 2, -2],
[-2, 2, -2, -2]])
In [54]: np.where(arr>0,2,arr) #只将正值设置为2
Out[54]:
array([[-0.80514356, 2. , 2. , 2. ],
[-1.82695429, -0.95471409, 2. , 2. ],
[ 2. , 2. , 2. , -0.36894073],
[-0.52798532, 2. , -0.33805238, -0.20382236]])

传递给where的数组大小可以不相等,甚至可以是标量值。

例子:

两个布尔型数组cond1和cond2,希望根据4种不同的布尔值组合实现不同的赋值操作:

for i in range(n):
...: if cond1[i] and cond2[i]:
...: result.append(0)
...: elif cond1[i]:
...: result.append(1)
...: elif cond2[i]:
...: result.append(2)
...: else:
...: result.append(3)
...:

这个for循环可以被改成一个嵌套的where的表达式:

np.where(cond1 & cond2,0,
np.where(cond1,1,
np.where(cond2,2,3)))

还可以利用“布尔值在计算过程中可以被当做0或1处理”,所以还能写成:

result=1*(cond1-cond2)+2*(cond2 &-cond1)+3*(cond1|cond2)

4.3.2 数学和统计方法

可以通过数组上的一组数学函数对整个数组或某个轴向的数据进行统计计算。

sum、mean以及标准差std等聚合计算,既可以当做数组的实例方法调用,也可以当做顶级numpy函数使用:

In [59]: arr=np.random.randn(5,4)
In [63]: arr
Out[63]:
array([[ 1.32400628, 0.66470491, -0.1503368 , -0.02298636],
[ 1.27163137, 0.0088513 , 2.02317321, 1.31306532],
[ 2.24286294, -0.1156468 , 1.51364621, -0.25358941],
[-0.84898582, -0.47816911, -0.50558895, 0.23782362],
[-1.97356203, -0.31476145, -0.9898961 , -0.52539505]])
In [60]: arr.mean()	#正态分布的数据
Out[60]: 0.22104236430646385
In [61]: np.mean(arr)
Out[61]: 0.22104236430646385
In [62]: arr.sum()
Out[62]: 4.4208472861292769

mean和sum这类的函数可以接受一个axis的参数(用于计算该轴向上的统计值),最终结果是一个少一维的数组:

In [64]: arr.mean(axis=1)
Out[64]: array([ 0.45384701, 1.1541803 , 0.84681824, -0.39873006, -0.95090366])
In [65]: arr.sum(0)
Out[65]: array([ 2.01595275, -0.23502115, 1.89099757, 0.74891811])
其他他cumsum和cumprod之类的方法不聚合,而产生一个由中间结果组成的数组:
In [69]: arr=np.array([[0,1,2],[3,4,5],[6,7,8]])
In [70]: arr.cumsum(0)
Out[70]:
array([[ 0, 1, 2],
[ 3, 5, 7],
[ 9, 12, 15]])
In [71]: arr.cumprod(0)
Out[71]:
array([[ 0, 1, 2],
[ 0, 4, 10],
[ 0, 28, 80]])
In [72]: arr.cumprod(1)
Out[72]:
array([[ 0, 0, 0],
[ 3, 12, 60],
[ 6, 42, 336]])

基本的统计方法:

sum 对数组中全部或某轴向的元素求和。

mean 算术平均数

std、var 分别为标准差和方差,自由度可调

min、max 最大值和最小值

argmin、argmax 分别为最大和最小元素的索引

cumsum 所有元素的累计和

cumprod 所有元素的累计和

4.3.3 用于布尔型数组的方法

上述的方法中,布尔值会被强制转换为1(True)和 0(False)。

因此,sum经常会被用来对布尔值数组中的True值计数:

In [73]: arr=np.random.randn(100)
In [74]: (arr>0).sum() #正值的个数
Out[74]: 48

另外还有两个方法any和all,他们对布尔型数组非常有用。any用于测试数组中是否存在一个或多个True,而all则检查数组中的所有值是否都是True:

In [75]: bools=np.array([False,False,True,False])
In [76]: bools.any()
Out[76]: True
In [77]: bools.all()
Out[77]: False

4.3.4 排序

numpy数组可以通过sort方法就地排序:

In [78]: arr=np.random.randn(8)
In [79]: arr
Out[79]:
array([-0.4352837 , -0.43593509, -0.16405127, -2.2530089 , 0.66142323,
-0.68429742, 0.90189568, -0.27002961])
In [80]: arr.sort()
In [82]: arr
Out[82]:
array([-2.2530089 , -0.68429742, -0.43593509, -0.4352837 , -0.27002961,
-0.16405127, 0.66142323, 0.90189568])

多维数组可以在任何一个轴向上进行排序,只需要将轴的编号传给sort即可:

In [84]: arr=np.random.randn(5,3)
In [85]: arr
Out[85]:
array([[ 0.66612646, 0.08637154, 0.43885613],
[-0.61377645, 0.09066706, 0.00999189],
[ 1.04957196, -0.22258295, 0.65275142],
[ 1.30820192, -0.89387398, 0.63250203],
[-0.46729776, -1.50409162, 1.61033952]])
In [86]: arr.sort()
In [87]: arr
Out[87]:
array([[ 0.08637154, 0.43885613, 0.66612646],
[-0.61377645, 0.00999189, 0.09066706],
[-0.22258295, 0.65275142, 1.04957196],
[-0.89387398, 0.63250203, 1.30820192],
[-1.50409162, -0.46729776, 1.61033952]])
In [88]: arr.sort(1)
In [89]: arr
Out[89]:
array([[ 0.08637154, 0.43885613, 0.66612646],
[-0.61377645, 0.00999189, 0.09066706],
[-0.22258295, 0.65275142, 1.04957196],
[-0.89387398, 0.63250203, 1.30820192],
[-1.50409162, -0.46729776, 1.61033952]])

计算数组分位数最简单的办法就是对其进行排序,然后选取特定位置的值:

In [90]: large_arr=np.random.randn(1000)
In [91]: large_arr.sort()
In [92]: large_arr[int(0.05*len(large_arr))]	# 5%的分位数
Out[92]: -1.6226524067201611

4.3.5 唯一化以及其他的集合逻辑

numpy提供了一些针对一维ndarray的基本集合运算,最常用的可能就是np.unique

np.unique 用于找出数组中的唯一值并返回已排序的结果:

In [93]: names=np.array(['Bob','Joe','Will','Bob','Will','Joe','Joe'])
In [94]: np.unique(names)
Out[94]:
array(['Bob', 'Joe', 'Will'],
dtype='|S4')
In [97]: ints=np.array([3,3,3,3,3,3,2,2,2,2,1,1,1,1,4,4])
In [98]: np.unique(ints)
Out[98]: array([1, 2, 3, 4])

另一个函数np.in1d用于测试一个数值中的值和另一个数组中的成员资格,返回一个布尔型数组:

In [99]: values=np.array([6,0,0,3,2,5,6])
In [100]: np.in1d(values,[2,3,6])
Out[100]: array([ True, False, False, True, True, False, True], dtype=bool)

4.4 用于数组的文件输入和输出

numpy能够读写磁盘上的文本数据或二进制数据。

4.4.1 数组以二进制格式保存到磁盘

np.save 和 np.load 是读写磁盘数组数据的两个主要函数。

默认状态下,数组是以未压缩的原始二进制格式保存到扩展名为.npy的文件中。

In [3]: arr=np.arange(10)
In [4]: np.save('some_array',arr)
In [7]: np.load('some_array.npy')
Out[7]: array([0, 1, 2, 3, 4, 5, 6, 7, 8, 9])

通过np.savez可以将多个数组保存到一个压缩文件中,将数组以关键字参数的形式传入即可:

In [8]: np.savez('array_archive.npz',a=arr,b=arr)

加载 .npz文件时,会得到一个类似于字典的对象,该对象会对各个数组进行延迟加载。

In [9]: arch=np.load('array_archive.npz')
In [10]: arch['b']
Out[10]: array([0, 1, 2, 3, 4, 5, 6, 7, 8, 9])

4.4.2 存取文本文件

从文本中加载一个非常标准的任务。

我们主要介绍pandas包中的read_csv和 read_table函数。

我们有时需要用np.loadtxt 或更专门化的 np.genfromtxt将数据加载到普通的数组numpy数组中。

这些函数许多选项都可以提供使用:指定各种分隔符,针对特定的转换器函数,需要跳过的行数等。

以一个简单的逗号分隔文件(CSV)为例子:

In [15]: !type array_ex.txt
0.58,0.17,1.04,1.13
0.19-0.63-0.94,0.12
-0.120.,26-0.6950.047
-1.484,0.004-0.7440.005
2.303,0.2003,1.6702-1.881
In [16]: arr=np.loadtxt('array_ex.txt',delimiter=',')
In [16]: arr
Out[16]:
array{[[0.58,0.17,1.04,1.13],
[0.19-0.63-0.94,0.12],
[-0.120.,26-0.6950.047],
[-1.484,0.004-0.7440.005],
[2.303,0.2003,1.6702-1.881]]}

np.savetxt执行的是相反的操作,将数组写到某种分隔符隔开的文本文件中。genformtxt和loadtxt差不多。

4.4 线性代数

线性代数(如矩阵乘法、矩阵分解、行列式以及其他方阵数学等)是任何数据库的重要组成部分。

numpy提供了一个数组方法也是numpy命名空间中的一个函数:

In [20]: x=np.array([[1.,2.,3.],[4.,5.,6.]])
In [21]: y=np.array([[6.,23.],[-1,-7],[8,9]])
In [22]: x
Out[22]:
array([[ 1., 2., 3.],
[ 4., 5., 6.]])
In [23]: y
Out[23]:
array([[ 6., 23.],
[ -1., -7.],
[ 8., 9.]])
In [25]: x.dot(y)
Out[25]:
array([[ 28., 36.],
[ 67., 111.]])

一个二维数组跟一个大小合适的一维数组的矩阵点积运算之后将会得到一个一维数组:

In [26]: np.dot(x,np.ones(3))	#dot是矩阵乘法
Out[26]: array([ 6., 15.])

numpy.linalg中有一组标准的矩阵分解运算以及诸如求逆和行列式之类的东西。

In [28]: from numpy.linalg import inv,qr
In [31]: X=np.random.randn(5,5)
In [32]: mat=X.T.dot(X)
In [33]: inv(mat)	#inv是矩阵的逆
Out[33]:
array([[ 0.13540336, -0.03570526, -0.03855285, 0.01570608, -0.01468477],
[-0.03570526, 0.65865699, 0.44353126, -0.27194311, -0.53503204],
[-0.03855285, 0.44353126, 0.95741232, -0.54344194, -0.48796139],
[ 0.01570608, -0.27194311, -0.54344194, 0.78882115, 0.81842571],
[-0.01468477, -0.53503204, -0.48796139, 0.81842571, 1.48301536]])
In [34]: mat.dot(inv(mat))
Out[34]:
array([[ 1.00000000e+00, 4.77693942e-18, 2.14638031e-17,
5.52800556e-17, 1.20779876e-17],
[ -9.31195688e-18, 1.00000000e+00, -1.82170005e-16,
1.31835755e-16, -1.40616613e-16],
[ -2.80570276e-18, -7.04587092e-17, 1.00000000e+00,
2.87660800e-17, 3.08691961e-17],
[ -9.24454202e-18, 8.20251278e-17, 4.16132855e-16,
1.00000000e+00, -9.51932954e-17],
[ -2.64100882e-17, 2.06118106e-16, 9.89240316e-17,
-3.91320305e-16, 1.00000000e+00]])
In [35]: q,r=qr(mat)	#计算QR的分解
In [36]: r
Out[36]:
array([[-7.80380924, -1.29182792, 0.46092023, 1.39731869, -1.22123513],
[ 0. , -4.00626236, 3.25373701, 4.54976221, -3.11812236],
[ 0. , 0. , -1.86896416, -2.82967231, 0.86140682],
[ 0. , 0. , 0. , -2.63952078, 1.76534633],
[ 0. , 0. , 0. , 0. , 0.54282638]])

4.5 随机数生成

In [38]: samples
Out[38]:
array([[-1.98098888, 0.15087649, -0.34650381, -0.3328599 ],
[ 1.07834818, -1.76824169, 0.81545647, -0.63545325],
[-0.95969246, -0.40560104, 0.16036069, 1.21139831],
[-1.02963311, 0.1659463 , -0.19833545, -0.57768052]])

python内置的random模块则只能一次生成一个样本值。

从下面的结果可以看出,如果需要产生大量的样本数据,numpy.random快了不止一个数量级。

In [7]: import numpy as np
In [8]: import numpy as py
In [9]: from random import normalvariate
In [10]: N=1000000
In [11]: %timeit samples=[normalvariate(0,1) for _ in xrange(N)]
1 loop, best of 3: 1.07 s per loop
In [14]: %timeit np.random.normal(size=N)
10 loops, best of 3: 50 ms per loop

4.6 范例:随机漫步

从0开始,步长1和-1出现的概率相等。通过内置的random模块以纯python的方式实现1000步的随机漫步:

In [1]: import random
In [2]: position=0
In [3]: walk=[position]
In [4]: steps=1000
In [5]: for i in xrange(steps):
...: step=1 if random.randint(0,1) else -1
...: position += step
...: walk.append(position)
...:

我用np.random模块一次性随机产生1000个“掷硬币的结果,将其分别设置为1或-1,然后计算累计和:

In [6]: nsteps=1000
In [9]: draws=np.random.randint(0,2,size=nsteps)
In [11]: steps=np.where(draws>0,1,-1)
In [12]: walk=steps.cumsum()
In [13]: walk.min()
Out[13]: -15
In [14]: walk.max()
Out[14]: 21

接下来我们想计算:随机漫步过程中的第一次到达某个特定的值的时间。

假设我们知道需要多久才能距离0点至少10步远。

np.ads(walk)>=10可以得到一个布尔型数组,表示的距离是否达到或者超过10。

可以使用argmax来解决这个问题:

In [15]: (np.abs(walk)>=10).argmax()
Out[15]: 139

4.6.1 一次模拟多个随机漫步

希望模拟多个随机漫步过程,只需要对上面的代码做一点点修改即可。只要给numpy.random的函数传入一个二元元组就可以产生一个二维数组。然后一次性计算5000个随机漫步过程的累计和了。

In [25]: nwalks=5000
In [26]: nsteps=1000
In [28]: draws=np.random.randint(0,2,size=(nwalks,nsteps))
In [29]: steps=np.where(draws>0,1,-1)
In [30]: walks=steps.cumsum(1)
In [31]: walks
Out[31]:
array([[ -1, 0, 1, ..., -12, -13, -12],
[ 1, 2, 1, ..., 24, 25, 24],
[ 1, 2, 3, ..., 48, 49, 48],
...,
[ -1, -2, -3, ..., 28, 29, 28],
[ -1, -2, -3, ..., -22, -23, -22],
[ -1, 0, 1, ..., -30, -31, -32]])

再来计算随机漫步的最大值和最小值

In [32]: walks.max()
Out[32]: 125
In [33]: walks.min()
Out[33]: -115

我们来计算30或-30的最小穿越时间。不是所有的数据都到了30,我们用any来处理一下:

In [34]: hits30=(np.abs(walks)>=30).any(1)
In [35]: hits30
Out[35]: array([ True, True, True, ..., True, False, True], dtype=bool)
In [36]: hits30.sum()