numpy操作技巧二三事

1、对于数值操作，少用list，尽量使用numpy array。

array是基于矢量化运算的，而list需要for循环操作。python的动态属性导致for loop速度比C慢上千倍（取决于for循环深度和数据量），而矢量化的array速度与C不相上下，因为numpy底层是用C语言优化过的。

要注意，numpy提供了一个叫vectorize的函数，不要被名字迷惑，它实现的还是一个python下的for loop，仅仅为了方便使用，对提高性能毫无作用，其返回一个array。

该函数使得输入的list可以逐元素操作，看起来像array。

用法：

import numpy as np
def myfunc(a, b):
     "Return a-b if a>b, otherwise return a+b"
     if a > b:
         return a - b
     else:
         return a + b

vfunc = np.vectorize(myfunc)
vfunc([1, 2, 3, 4], 2)#返回array([3, 4, 1, 2])

2、numpy数据类型

python内置的数据类型有整型int、浮点型float、复数complex和bool。但是这些数据类型精度取决于计算机和python版本。

numpy提供了非常丰富的数据类型，可以很方便的算出数据量对应的内存。

整型有9种，分别为

inti—平台决定的整型，要么int32或int64\
int8—有符号整数[-2**7， 2**7-1]\
int16—有符号整数[-2**15， 2**15-1]\
int32—有符号整数[-2**31， 2**31-1]\
int64—有符号整数[-2**63， 2**63-1]\
uint8—无符号整数[0， 2**8-1]\
uint16—无符号整数[0， 2**16-1]\
uint32—无符号整数[0， 2**32-1]\
uint64—无符号整数[0， 2**64-1]\

浮点型有3种：

float16—半精度浮点数：1个符号位，5个指数位，10个小数位\
float32—单精度浮点数：1个符号位，8个指数位，23个小数位\
float64—半精度浮点数：1个符号位，11个指数位，52个小数位\

复数有两种：

complex64—复部和实部都是float32\
complex128—复部和实部都是float64\

bool型与python内置的一样。

其中大多数数据类型都可以互相转换，但是复数不能随意转换为别的类型。
数据类型对象是numpy.dtype类的实例。

可用dtype的属性itemsize来检查每个元素的字节数。

a = np.arange(5)
a.dtype.itemsize#返回4（依赖于平台）

3、array的维度操作注意事项

对于一维array，slice操作与list完全一致。

a=np.arange(3)
a[::-1]#返回array([2,1,0])

flatten vs ravel：

flatten操作后的数据是原始数据的copy，所以随后的操作不会影响原始数据；而ravel后的数据与原始数据仍然占据同一内存，所以修改后会影响原始数据。

因为ravel不需要copy，所以性能更优，但是使用时要仔细考虑清楚。

reshape vs resize：

reshape与ravel类似，更改后的数据与原始数据占用同一内存。

resize只能在当前数据上操作，即没法作为赋值对象。

hstack vs vstack vs dstack vs concatenate：

对于2D array，hstack((a,b))与concatenate((a,b), axis=1)产生同样效果。可见hstack是按列堆积，行为等同于铺地砖。

对于2D array，vstack((a,b))与concatenate((a,b), axis=0)产生同样效果。可见vstack是按行堆积，行为等同于盖楼房。

c=dstack((a,b))按第三个维度堆积，所以c[:, :, 0] == a， c[:, :, 1] == b。一般用于图像处理，比如把几个2D图像array堆叠在一起产生新的图像。

对于1D array，要实现按列堆积，需要用column_stack((a,b))，在2D情况下，等同于hstack。\
要实现按行堆积，要用row_stack((a,b))，在2D情况下，等同于vstack。

hsplit vs vsplit vs dsplit vs split:

这些操作分别对应于上面四个stack函数的逆操作，即怎么stack的就怎么split开来。

4、计算一个array占用的内存字节数

任一个元素占用的字节数=a.itemsize，一个array拥有的总元素=a.size，所以一个array占用的总字节数=a.itemsize*a.size。

但是numpy自带的电池a.nbytes做了前面乘法所做的事。

5、flat有啥用

flatten、ravel、flat都将array变成1D，那么要flat何用？flatten返回一个原始数据1D化之后的copy，需要花费额外的内存，尤其对于大的数据量；而ravel虽然不产生copy，但其内存与原始数据指向相同，修改后会改变原始数据；那么有没有一种既不占用内存有不会改变原始数据的方法，flat正好做了这个工作。

flat返回一个生成器。这就令人恍然大悟了，生成器是python路上的利器，可以参看我以前的文章迭代对象、迭代器、生成器浅析。有了flat操作，那么利用循环进行逐元素操作就很easy了。

6、copy vs view

copy将原始数据非配到别的内存块里，所以与原始数据毫无瓜葛；而view的内存与原始数据一致，所以是关联的；可用函数np.may_share_memory(a, b)查看两个array是不是同一个内存地址，如果返回True则说明b是a的view，False则为copy。

何时用copy何时用view需要因时制宜。

7、花式索引 Fancy indexing

一个老掉牙的例子，将lena的两个对角线变成黑色。可用for循环来完成，但Fancy indexing更高效更pythonic。

尴尬的是lena没了，貌似因为版权问题从scipy remove掉了。那就用ones来代替吧。

import numpy as np
import matplotlib.pyplot as plt

lena = np.ones((20,20))
xmax = lena.shape[0]
ymax = lena.shape[1]

# Fancy indexing
# 斜率为正的对角线
lena[range(xmax), range(ymax)] = 0

# 斜率为负的对角线
lena[range(xmax-1,-1,-1), range(ymax)] = 0

plt.imshow(lena)
plt.show()

8、使用bool变量和ix_对array index操作

ix_函数接收N个1D序列（可以是list或array），返回N个N维index。如果要对一个2D array进行操作，则为ix_(xindices, yindices)。

假设原始array为2X3，则ix_将以如下方式index原始图像：

[xindices[0],yindices[0]]
[xindices[0],yindices[1]]
[xindices[0],yindices[2]]
[xindices[1],yindices[0]]
[xindices[1],yindices[1]]
[xindices[1],yindices[2]]

import scipy.misc
import matplotlib.pyplot as plt
import numpy as np

lena = np.ones((20,20))
xmax = lena.shape[0]
ymax = lena.shape[1]

lena[range(xmax), range(ymax)] = 0
lena[range(xmax-1,-1,-1), range(ymax)] = 0

xarr = np.arange(xmax)
np.random.shuffle(xarr)#打乱顺序
yarr = np.arange(ymax)
np.random.shuffle(yarr)#打乱顺序

plt.imshow(lena[np.ix_(xarr,yarr)])#乱序显示图片

用bool值索引array就更简单了，可以直接在index里显式输入True或Fale，也可以通过条件比较生产bool值。

h=np.arange(20)
h[(h>h.max()*0.2) & (h<(h.min()+1)*10)]=0
#大于3小于10的数被置为0

Of course，python还有更pythonic的方法，那就是mask。要注意，满足条件时mask返回True，否则返回False，即被masked的区间为True。

上面的例子可用mask来实现：

h1=np.arange(20)
h1ma = np.ma.masked_where((h1>h1.max()*0.2) & (h1<(h1.min()+1)*10), h1)#生产mask
print(h1ma)
#masked_array(data = [0 1 2 3 -- -- -- -- -- -- 10 11 12 #13 14 15 16 17 18 19],
#             mask = [False False False False  True  True  #True  True  True  True False False
# False False False False False False False False],
#       fill_value = 999999)
h1[h1ma.mask]=0#大于3小于10的数被置为0
np.array_equal(h,h1)#返回True
(h==h1).all()#返回True

~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
（内容同步更新到微信公众号python数学物理，微信号python_math_physics）