微信公众号：cha-diary

大家好！

在这一段时间呢，我主要针对自己调研的一个机器学习模型进行了复现。俗话说得好，读论文一时爽，到复现火葬场。我可能现在后悔就后悔在习惯了Py的解释型环境，把C的那一套函数式编程方法忘得一干二净……现在每天面对自己写的shit mountain不知道该说什么好……

当然了，也不完全是一无是处啦！事实上基于统计和优化的机器学习模型，虽然说就我目前的情况来看，要想在效率和规模上胜过深度学习已经比较困难，但这并不代表说，我们不能够从机器学习算法中学到其它有趣的技巧。其实这也是我们这一节的一个主要内容：在数值计算中，我我们有什么办法加速？

当然我们不可能寥寥几句说完数据处理的一些常规操作，所以我们这一节也会补上剩下的这些内容。

那我们开始吧。当然事先请不要忘记导入

import numpy as npimport pandas as pdimport matplotlib.pyplot as plt

这一节我们还需要做数值实验，所以别忘了还有这个

import time

8. 数据排序

(从8开始的原因是第一节是到7结束的……)数据排序也是很常见的操作，比方说在电商中，你需要提取某一天的物品销量前十名，那自然对销量进行排序是重要的。当然还有对索引排序(为了使得数据原始的顺序得以恢复)等。

排序方法在numpy与pandas中各有一种。首先是numpy中的我们有

A = np.array([1,3,4,8,7,6,9,0,5])A.sort()A'''Result: array([0, 1, 3, 4, 5, 6, 7, 8, 9])'''

并不困难，对这个array调用一个sort方法就好。不过如果是在二维及以上的数据，想进行分行分列的排序，应该如何去做呢？

A = np.array([1,3,4,8,7,6,9,0,5]).reshape(3,3)A.sort(axis = 0)A'''Result: array([[1, 0, 4],       [8, 3, 5],       [9, 7, 6]])'''

这里我们多了一个参数axis，它为0即说明是按列排序，如果它为1那么就是按行排序了。默认情况下是按行排序的。

当然，不乏可能有人有一些特别的想法。比方说

Example 1:对数据进行全排序，并按照列重排。比方说如果输入是1-9的一个乱序排列，输出为

array([[1, 4, 7],       [2, 5, 8],       [3, 6, 9]])

实现这个目标并不困难，但是稍显复杂。

A = np.array([1,3,4,8,7,6,9,2,5])A.sort()A = A.reshape(3,3, order = 'F')A'''Result:array([[1, 4, 7],       [2, 5, 8],       [3, 6, 9]])'''

也就是说，我们先正常排序，然后考虑使用reshape方法，这里注意多了一个参数order = 'F'，注意这个'F'并不是与“列名”有关的一个首字母，而是Fortran的意思。具体的解释下面这个官方文档也说的比较清楚。

https://docs.scipy.org/doc/numpy-1.15.1/reference/generated/numpy.reshape.html

虽然这个需求倒不是很多见。不过reshape有一个用法还是挺重要的，我们用一个例子来看

A = np.array([[1, 4, 7],           [2, 5, 8],           [3, 6, 9]])A = A.reshape(-1,)A'''Result: array([1, 4, 7, 2, 5, 8, 3, 6, 9])'''

reshape的含义就是“重铸”，一般都是转换矩阵的行列数。如果你设置成了-1，那么程序会自动计算这一个数。比方说原始有6个数据，你的行列数设置为(2,-1)，这就是(2,3)的意思。如果仅仅只有-1，那么矩阵就会被“拍平”成为一行，相信大家在卷积神经网络中对于这个操作并不陌生。

注意：sort方法是会改变被操作对象的值的。也就是说

A.sort(axis = 0)A = A.sort(axis = 0)

这两句话其实是一个意思。不过一般我还是会写成后者，也是为了阅读方便。

说到这里我们来看一下pandas的情况，和numpy不同的是pandas的DataFrame一般是含表头和索引的数据表，因此我们如果是对它排序，一般都是要根据某一个特征排序。我们还是使用上一节的数据集来举例。希望你没有忘记那个简单的表。

Table = {'date': ['2019/6/1', '2019/7/2', '2019/6/6', '2019/6/17', '2019/7/4', '2019/6/13', '2019/6/14', '2019/6/21', '2019/6/17'], 'order_id': [i+1 for i in range(9)],'commodity_code': ['S1', 'S2', 'S3', 'S5', 'S5', 'S2', 'S9', 'S11', 'S9'], 'commodity_name': ['标准美式','瑞纳冰', '加浓美式', '拿铁', '拿铁', '瑞纳冰', '菠萝卷', '坚果', '菠萝卷'], 'category_name': ['饮品', '饮品', '饮品', '饮品', '饮品', '饮品', '食品', '食品', '食品']}Table = pd.DataFrame(Table)

不如我们就按照commodity_code排序吧，这样的话方法就应该写成这样。

Table = Table.sort_values(by = ['commodity_code'])'''Result:        dateorder_idcommodity_codecommodity_namecategory_name02019/6/11S1标准美式饮品72019/6/218S11坚果食品12019/7/22S2瑞纳冰饮品52019/6/136S2瑞纳冰饮品22019/6/63S3加浓美式饮品32019/6/174S5拿铁饮品42019/7/45S5拿铁饮品62019/6/147S9菠萝卷食品82019/6/179S9菠萝卷食品'''

所以其实你看，数据按照commodity_code的大小顺序(当然这里的大小顺序是字符意义上的)做了一个排序，你可以看到，包括索引在内的所有列都按照commodity_code的顺序重新排列了。这当然是很重要的，你肯定不希望你在对数据的某一列进行排序的时候，其它的列都保持不动，那样的话数据本身的内容就改变了(如果每一行代表一个观测的话，这就相当于人为改变了观测值的某一个特征)，所以常规的排序肯定是不行的。

在我们完成数据分析后，现在自然希望做的是还原数据。在这里我们就是希望能够把索引再重新排回去。这样的话就需要用到另一个排序的方法

Table = Table.sort_index()

这样就可以了。

注意：我们这里的所有排序方法都是升序的，如果你要降序，那么记得方法内部设置一个参数ascending = False

9. 缺失值，异常值，重复值处理

大数据处理中缺失值是很常见的，可能是因为用户不登记信息，也有可能是数据库出现短时间的宕机等。不过它的处理并不是特别困难，因为缺失值在一个正常的数据集里占比通常都非常小，所以一般的方法都够用了。换句话说如果缺失值占了很大一部分，那你应该考虑换一组数据……

下面这个网站通俗易懂的说明了缺失值的几种处理方式，我就不班门弄斧了。

http://www.sohu.com/a/305803505_697896

那么什么叫异常值呢？统计意义上就是大家比较熟悉的outlier，但是实际情况下不仅仅是统计意义上超出常规范围的值才是outlier。举个例子，在一个统计用户性别的数据集中，男生是1，女生是2的话，那么这一列数据集有没有可能出现3？这个答案是有！如果是这种情况，那么大家就需要自己写函数去判断和处理了。

最后，什么叫重复值？在实际的数据分析中，可能因为疏忽或者数据录入的问题会出现这个情况(比方说学生信息表，一个学号可能在数据表中录给了多个人，或者不同的学号录入给了同一个人)。我们拿上面的数据集来展示如何去除重复值。

Table = Table.drop_duplicates(subset=None, keep='first', inplace=False)Table'''Result:        date order_id     commodity_codecommodity_namecategory_name02019/6/11S1标准美式饮品12019/7/22S2瑞纳冰饮品22019/6/63S3加浓美式饮品32019/6/174S5拿铁饮品42019/7/45S5拿铁饮品52019/6/136S2瑞纳冰饮品62019/6/147S9菠萝卷食品72019/6/218S11坚果食品82019/6/179S9菠萝卷食品'''

好像什么都没有发生。这是因为我们设置的subset参数为None，这种情况下，只有一个观测的所有特征出现重复才会被删除。那么比方说我们不希望删除commodity_code这一列中重复的部分，那如何操作呢？

Table = Table.drop_duplicates(subset=['commodity_code'], keep='first', inplace=False)Table'''Result:        date   order_idcommodity_codecommodity_namecategory_name02019/6/11S1标准美式饮品12019/7/22S2瑞纳冰饮品22019/6/63S3加浓美式饮品32019/6/174S5拿铁饮品62019/6/147S9菠萝卷食品72019/6/218S11坚果食品'''

你看，这个时候，commodity_code这一列，出现重复的元素只有顺序的第一个被保留了下来。这是因为我们的keep参数设置为了'first'，当然如果你希望是最后一个，那就改为'last'，一般情况下这够用了。

最后我们提一下这个inplace参数。如果这个参数设置为了True，那么它会改变原对象。就和上面介绍的sort方法一个情况了。

注意：在这个时候，索引会出现“断档”，如果你要拼接数据可能会出现麻烦，因此最好再加一个reset_index方法，具体可以参见上一节。

10. 列表生成式

列表生成式算是一个常见的比较地道的用法，可以很方便的节省代码量。不过数值实验告诉我，它的效率提升不是太显著。比方说下面这个任务

Example 2:找出一组列表中所有为5的数对应的索引。

这个其实列表生成式就可以写得很短

A = np.array([1,1,5,4,9,5,4,3,7,10])B = [i for i in range(A.shape[0]) if A[i]== 5]B'''Result:[2, 5]'''

为了方便理解它的逻辑，我们把它拆成一般的for循环来看

A = np.array([1,1,5,4,9,5,4,3,7,10])B = []for i in range(A.shape[0]):    if(A[i] == 5):        B.append(i)B

所以实际上就是做一个循环，然后每一次循环再做判断，把满足条件的元素放到我们的列表内。

还有一个比较常见的用法是赋二值，具体就是下面这个意思

Example 3:在一列数中修改值为5的数为1，否则为2

代码也不困难

A = np.array([1,1,5,4,9,5,4,3,7,10])A = [1 if A[i] == 5 else 2 for i in range(A.shape[0])]A'''Result: [2, 2, 1, 2, 2, 1, 2, 2, 2, 2]'''

这里我们用到了Py的三元表达式，感兴趣的也可以记一下。

当然，你也可以用它来做二重循环，只需要后面写两个for即可，具体的大家可以百度搜索，这里不再展开。

11. 高阶函数map

其实我们在上一节已经在“数据筛选”那一部分中提到过高阶函数filter，这里提到map是因为在做批量处理的时候，它的效率和轻便性让它得到了很高的评价。实际的业务中，我们经常会碰到文本数据，但是它的前缀往往不是我们想要的内容(比方说我们数据集commodity_code中的开头的S其实就不是很重要)。抽象一下就是下面这个任务

Example 4:删除commodity_code列中的前缀字母'S'

对应代码如下

Table['commodity_code'] = pd.DataFrame(map(lambda x:x[1:],Table['commodity_code']))Table'''Result:         date  order_idcommodity_codecommodity_namecategory_name02019/6/111标准美式饮品12019/7/222瑞纳冰饮品22019/6/633加浓美式饮品32019/6/1745拿铁饮品42019/7/455拿铁饮品52019/6/1362瑞纳冰饮品62019/6/1479菠萝卷食品72019/6/21811坚果食品82019/6/1799菠萝卷食品'''

你可以看出来，其实就是每一个元素只从索引为1(也就是说实际上是第2个)的部分开始提取。

抽象出来写，高阶函数的语法是

map(f, A)

也就是对于A这个对象的每一个元素使用函数f，然后将结果储存。这里要注意的是高阶函数本身是一个迭代器，你需要把它实例化，这里就是把它转换成一个DataFrame再重新拼了进去。

注意：还用到了一个函数就是匿名函数lambda，上一节我们也提到过。这当然是一个函数，你可以理解为它就是下面这个函数

def f(x):    return x[1:]

但是不觉得它实在是太简单了嘛？所以为了防止命名麻烦和这种过于小的函数占用内存，就可以使用匿名函数完成这些很简单的小操作，代码简单又节省空间。

Python加速

在推荐系统中有一种算法叫作协同聚类(co-clustering)，在商品数据中，这样的方法可以用来发现某一部分的人群对于某一种商品的偏好情况。这和一般的聚类还不相同，主要在于这里，同一类人群对于某一种商品的偏好相似，但是他们对于其它的商品的偏好可能完全不同。这种情况下如果是一般的聚类，是不可能被分到一起的。

很多常见的机器学习的算法会考虑到概率分布，而一般情况下我们也会将销量数据抽象成为一个大的离散二位概率分布矩阵。这个时候就会涉及到大量的数值计算，也就是py需要加速的地方了。我们主要用一个问题来引入，分别说明py的一些加速技巧。

12. Python加速(1)：使用numpy的运算方法

要想构造离散概率密度矩阵，很明显，针对每一个元素

你都需要计算出它们对应的值。当三个依赖的矩阵都已经计算好的时候，似乎一个for循环就可以解决问题了，但你不要忘记我们说的：Py中尽量要避免使用循环，因为太慢了……那怎么办？

答案很简单：多考虑一些numpy本身的特性。这里我们写一个数值例子来看

A = np.tile(1, (5000, 5000))B = np.tile(2, (5000, 5000))def Test_Normal(A,B):    C = np.random.randn(5000, 5000)    for i in range(5000):        for j in range(5000):            C[i,j] = A[i,j] * B[i,j]            def Test_Numpy(A,B):    return A*BStart = time.clock()Test_Normal(A,B)Elapsed = time.clock() - Startprint("Step 1:", Elapsed)Start = time.clock()Test_Numpy(A,B)Elapsed = time.clock() - Startprint("Step 2:", Elapsed)'''Result:Step 1: 26.693720300001587Step 2: 0.0829771000026085'''

你可以看到，如果是针对传统的矩阵乘法操作，那么可以直接使用一步A*B即可。这是因为numpy中，如果直接写成A*B，那么默认为按元素相乘，而不是矩阵乘法。你可以发现，这样的操作快了几百倍(当然实际上创建C这个大的array需要一秒多的时间，不过就算去掉这个，常规循环也要25秒多23333)。

如果希望做矩阵乘法，则使用

np.dot(A, B)

即可。

我们加大一点难度。在这个问题中，我们的矩阵都已经事先计算好。那么如果，在计算中不是这样呢？还是上面那个例子，

是条件概率密度，因此按照正常的情况来计算，实际上应该是一个一维的向量。这个时候很明显，如果要使用循环和赋值操作，就会出现重复。因为针对用户

在计算它们的概率的时候都需要涉及到同一个

这就相当于说我们要一遍一遍的把

过一遍，很明显不是我们希望看到的。这个时候怎么办呢？其实很简单，别忘记上一节我们提到过一个东西叫做np.tile()。你可以看出来，我们它在针对“堆叠数据”这个操作上具有非常高效的性能，所以我们并不用害怕太多。

比方说针对一个6*3的矩阵

和一个6*1的向量 ，你希望它们俩做乘法，那你其实只需要用np.tile，把它的列扩充三倍即可。

13. Python加速(2)——减少重复运算

其实，上面这个问题也是减少重复的例子，因为我不希望一遍又一遍的遍历我的向量，所以我把它扩张成为了一个矩阵。那么还有一种情况，我们直接上这个数值实验来说明

def Test1(A, B):    C = np.array([i for i in range(1000) for j in range(1000)]).reshape(1000,1000)    for i in range(1000):        for j in range(1000):            C[i,j] = A.sum(1) * B.sum(1)def Test2(A, B):    Temp1 = A.sum(1)    Temp2 = B.sum(1)    C = np.array([i for i in range(1000) for j in range(1000)]).reshape(1000,1000)    for i in range(1000):        for j in range(1000):            C[i,j] = Temp1 * Temp2

你感觉这两个函数哪一个写的更好一点？很明显是第二个，原因在于它移开了每一步循环中都会使用到的计算。这里就相当于说把求和聚合的部分直接保存到了外面，计算好再拿进来。这个思维不仅仅是py，在任何数值计算的编程中都很重要。

当然在这个问题中，其实我们不太需要关注这个，毕竟py中我们基本上不写for循环对吧……

14. Python加速(3)——pandas，numpy与numba

numba这个东西在某些情况下，那个加速效果杠杠的。还是那个问题，我们再加大一点难度。如何扩张

请注意我说的不是计算，是扩张。因为分的类别数和本身的样本个数很明显不会相同。比方说

是一个2*3的矩阵，但是

却需要是6*4的，这个时候因为每一个商品和用户都会被归到不同的类别，所以np.tile很明显无法再使用了(因为这不是重复计算)。

好的，似乎要使用for循环了，那么是不是一定就火葬场了呢？不一定！这里我们介绍一个库numba。使用的话需要引入包

import numba

那么好，我们可以写我们的数值实验了。

def Test_Normal(A,B):    C = np.random.randn(5000, 5000)    for i in range(5000):        for j in range(5000):            C[i,j] = A[i,j] * B[i,j]            def Test_Numpy(A,B):    C = np.random.randn(5000, 5000)    return A*B@numba.jitdef Test_Jit(A,B):    C = np.random.randn(5000, 5000)    for i in range(5000):        for j in range(5000):            C[i,j] = A[i,j] * B[i,j]Start = time.clock()Test_Jit(A,B)Elapsed = time.clock() - Startprint("Step 1:", Elapsed)Start = time.clock()Test_Normal(A,B)Elapsed = time.clock() - Startprint("Step 2:", Elapsed)Start = time.clock()Test_Numpy(A,B)Elapsed = time.clock() - Startprint("Step 3:", Elapsed)'''Result:Step 1: 1.4499750999998469Step 2: 23.27106929999991Step 3: 1.419889000000012'''

你看，它就是在函数上方添加了一个装饰器numba.jit而已。

注意：

C = np.random.randn(5000, 5000)

这一步会消耗绝大部分的时间，所以实际上你可以看出，如果你事先不考虑这一步运算，其实numba可以把for循环的速度降到和numpy的矩阵运算差不多的等级(也就是快了几百倍)。

并不用感到奇怪，这是因为numba.jit这个装饰器本身就会把常规的python语句转换为numpy的，反过来说需要注意的是，如果你对含有大量numpy语句的函数进行这样的处理，实际上并不会得到什么好的优化。

下面，让我们欢迎pandas登场做小白鼠，如果pandas做索引循环，它的速度怎么样呢？

def Test_Pandas(A,B):    C = np.random.randn(100, 100)    C = pd.DataFrame(C)    for i in range(100):        for j in range(100):            C.iloc[i,j] = A.iloc[i,j] * B.iloc[i,j]def Test_Normal(A,B):    C = np.random.randn(5000, 5000)    for i in range(5000):        for j in range(5000):            C[i,j] = A[i,j] * B[i,j]Start = time.clock()Test_Normal(A,B)Elapsed = time.clock() - Startprint("Step 1:", Elapsed)A = pd.DataFrame(A)B = pd.DataFrame(B)Start = time.clock()Test_Pandas(A,B)Elapsed = time.clock() - Startprint("Step 2:", Elapsed)A = np.array(A)B = np.array(B)'''Result: Step 1: 22.2333093000002Step 2: 3.1992385000003196'''

你可以看出来，论数据量pandas那个比numpy那个少了2500倍，但是时间其实没少太多，算一下速度其实pandas比numpy要慢几百倍。所以这就是为什么我一般都会把创建好的DataFrame转为array，然后到用到pandas的一些方法的时候再转回去……pandas的索引是真的很慢，并且不可以使用numba调优！当然我也不知道为什么，欢迎评论区能够有人给出答案……

15. Python加速(4)——numpy的广播机制

最后，我们来谈谈广播(broadcast)机制。首先自然要解释一下什么是广播，广播你可以理解为是numpy在遇到矩阵大小不一致，但是需要做矩阵运算的时候所采取的处理方法。举个复杂的例子。

简单理解一下，其实就是说每一行的对应元素与每一列的对应元素作相除，然后每一行得到的结果对应求和就可以得到我们的三个结果。你当然可以使用我们上面已经说过的方法，使用np.tile()扩张之后再操作，但是有那个必要吗？没有！事实上非常简单。

A = np.array([i+1 for i in range(4)]).reshape(1,4)B = np.array([i+1 for i in range(12)]).reshape(3,4)C = A/Bnp.sum(C, 1)'''Result: array([4.        , 1.46190476, 0.91717172])'''

你可以看到只有一个A/B，然后按列求和即可。这也就是Python自带的广播机制。具体如下

1. 找到两个目标矩阵的行列最大值，比方说A为(1,3)，B为(2,1)，那么行列最大值就是(2,3)2. 对比目标矩阵的维度，如果对应维度非行列最大值的因数，则无法广播，程序报错。否则按照倍数拉伸为新矩阵，其维度即为行列最大值。3. 按元素对元素的方法运算

比方说这里就是一个1*4的矩阵与一个3*4的矩阵做除法，根据广播原则，我们的1*4的矩阵就会被拉伸为一个3*4的矩阵，然后点对点运算，所以这里只需要一个A/B运算就完成任务了。

那么我们再加大一点难度，考虑这个问题。

虽然只是多了一行，但是这个时候广播就没那么容易了，因为2*4与3*4的矩阵本身是无法广播的。难道要使用我们传统的for循环了吗？也许你会想到使用上面说的numba，但是在实际的代码中，可能一步循环会涉及到多步操作，而且主要都是numpy或者pandas的方法，所以这个方法一般也是不行的。怎么办？

欢迎我们的高维数组登场！没有条件创造条件也要上！既然我们这里的长宽索引分别为(2, 4) 和(3, 4) ，为什么不考虑升维呢？很明显，我们希望做的就是，把矩阵 A对应的维度扩充3倍，而把矩阵 B对应的维度扩充2倍，这样的话结构就完整了，也就可以考虑广播了。

所以我们考虑这样的变换

用reshape方法是可以很轻松的得到这个结果的。

这样根据广播，维数最大分别为(3, 2, 4) ，很明显我们希望得到的是一个 (3, 2)的矩阵，所以自然我们求和就是要对第三个维度进行求和。具体的代码如下

A = np.array([i+1 for i in range(8)]).reshape(2,4)B = np.array([i+1 for i in range(12)]).reshape(3,4)A = A.reshape(1, 2, 4)B = B.reshape(3, 1, 4)C = A/Bnp.sum(C, 2)'''Result:array([[ 4.        , 12.33333333],       [ 1.46190476,  4.        ],       [ 0.91717172,  2.45858586]])'''

你看，是不是我们想要的结果？这样做了一个变换之后，在针对矩阵特别大的情况下，优化效果是非常显著的，差不多也是几千几万倍速度的提升。

当然，实际上还有更好的算法！不过就交给读者思考了

小结

除去比较基本的一些数据分析的操作外，本节也重点关注了一些Python数值计算的加速技巧。但不同于matlab，Python的很多计算不是特别的占性能优势，因此需要一些技巧去优化。numpy和pandas是数值计算的大头，也是我们所重点应用和做数据分析的必备工具，希望各位需要使用Python实现一些数值计算问题的同学们能够从中感受到它们的奇妙。