numpy：python数据领域的功臣

前言

numpy对python的意义非凡，在数据分析与机器学习领域为python立下了汗马功劳。现在用python搞数据分析或机器学习经常使用的pandas、matplotlib、sklearn等库，都需要基于numpy构建。毫不夸张地说，没有numpy，python今天在数据分析与机器学习领域只能是捉襟见肘。

什么是一门好的数据分析语言

数据分析面向的数据大多数是二维表。一门好的数据分析语言，首先需要能够直接有个数据结构存下这个二维表，然后要配上一套成熟的类SQL的数据操作接口，***要有一套好用的可视化工具。R语言就是一个极好的典范：用内置的data.frame结构做数据的存储;data.frame本身提供足够强大的数据操作能力，另有dplyr、tidyr、data.table、plyr、reshape2等库提供更好用更高效的数据操作能力;在绘图上，除了基本的plot功能外，还提供了ggplot2这样一套优雅的绘图语言，还通过htmlwidget库与javascript各种绘图库建立了紧密的联系，让可视化的动态展示效果更进一步。Excel也是一个极好的例子，有单元格这种灵活的结构为数据存储做支撑，有大量的函数实现灵活的操作，也有强大的绘图系统。

python目前在数据分析领域也已经具备了相当可观的能力，包括pandas库实现的DataFrame结构，pandas本身提供的数据操作能力，matplotlib提供的数据可视化能力，而这一切都离不开numpy库。

什么是一门好的机器学习语言

一般来讲，一门好的机器学习语言在数据分析上也一定很吃得开，因为数据分析往往是机器学习的基础。但是机器学习的要求更高，因为在模型训练阶段往往需要较为复杂的参数估计运算，因此语言需要具备较强的科学计算能力。科学计算能力，最核心的就是矩阵运算能力。关于矩阵运算能力，这篇文章对各种语言有很好的比较。

如果没有numpy，python内部只能用list或array来表示矩阵。假如用list来表示[1,2,3]，由于list的元素可以是任何对象，因此list中所保存的是对象的指针，所以需要有3个指针和三个整数对象，比较浪费内存和CPU计算时间。python的array和list不同，它直接保存数值，和C语言的一维数组比较类似，但是不支持多维，表达形式很简陋，写科学计算的算法很难受。numpy弥补了这些不足，其提供的ndarray是存储单一数据类型的多维数组，且采用预编译好的C语言代码，性能上的表现也十分不错。

python***的机器学习库sklearn构建在numpy之上，提供了各种标准机器学习模型的训练与预测接口，其中模型训练接口的内部实现是基于numpy库实现的。比如很常见的线性回归模型，参数估计调用的是numpy.linalg.lstsq函数。

numpy的核心结构：ndarray

以下内容摘录自用Python做科学计算

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

ndarray是numpy的核心数据结构。我们来看一下ndarray如何在内存中储存的：关于数组的描述信息保存在一个数据结构中，这个结构引用两个对象，一块用于保存数据的存储区域和一个用于描述元素类型的dtype对象。

 数据存储区域保存着数组中所有元素的二进制数据，dtype对象则知道如何将元素的二进制数据转换为可用的值。数组的维数、大小等信息都保存在ndarray数组对象的数据结构中。

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

以下内容总结自Numpy官方文档Numpy basics

关于ndarray的索引方式，有以下几个重点需要记住：

  • 虽然x[0,2] = x0，但是前者效率比后者高，因为后者在应用***个索引后需要先创建一个temporary array，然后再应用第二个索引，***找到目标值。
  • 分片操作不会引发copy操作，而是创建原ndarray的view;他们所指向的内存是同一片区域，无论是修改原ndarray还是修改view，都会同时改变二者的值。
  • index array和boolean index返回的是copy，不是view。

关于上面列举的分片操作不会引发copy操作，我们来进一步探讨一下。先看一下numpy的例子：

 再来看一下R的例子：

 可以看到numpy和R在矩阵的分片操作有不同的设计理念：在R里分片操作会引起数据的复制，在numpy里不会。事实上，R的设计理念很多时候可以用一句话来概括：copy on modify，一旦对数据有修改就会引起内存上的复制操作，这个操作要花不少时间，因此经常会听到人们抱怨R费内存且速度慢。所以，我们可以看到numpy在处理这件事情上明显要用心很多，根据场景设计了不同的策略，不是简单地采用R的一刀切方式。当然，这也带来了一些学习成本，需要对numpy足够熟悉才能避免踩坑。R社区里对copy on modify的哲学也有诟病并在努力改变，比如同是data.frame操作库的data.table和dplyr，data.table性能比dplyr高很多，部分原因也是data.table规避了copy on modify的方式。

Structured Array

根据numpy的官方文档，定义结构化数组有四种方式。本文采用字典方法，通过定义一个dtype对象实现，需要指定的键值有names和formats。

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

我们用IPython的计时函数看一下提取数据的效率： 

 
 
    
  
  
  1. %timeit a[1] 
    2. 
  
  
  2. %timeit a['name'] 
    3. 
  
  
  3. %timeit a[1]['name'] 
    4. 
  
  
  4. %timeit a['name'][1]  
    5.

输出结果如下：

 
 
    
  
  
  1. The slowest run took 46.83 times longer than the fastest. This could mean that an intermediate result is being cached. 
    2. 
  
  
  2. 1000000 loops, best of 3: 153 ns per loop 
    3. 
  
  
  3. The slowest run took 34.34 times longer than the fastest. This could mean that an intermediate result is being cached. 
    4. 
  
  
  4. 10000000 loops, best of 3: 174 ns per loop 
    5. 
  
  
  5. The slowest run took 13.00 times longer than the fastest. This could mean that an intermediate result is being cached. 
    6. 
  
  
  6. 1000000 loops, best of 3: 1.08 µs per loop 
    7. 
  
  
  7. The slowest run took 9.84 times longer than the fastest. This could mean that an intermediate result is being cached. 
    8. 
  
  
  8. 1000000 loops, best of 3: 412 ns per loop  
    9.

从上面的结果，我们发现，获取相同的数据有多种操作，不同的操作性能差别很大。我做了一个推测，纯粹是瞎猜：numpy在建立结构化数组时，将整个结构体连续存储在一起，即按行存储，因此a[1]的速度最快;但是为了保证提取列的效率，对a['name']建立了索引，因此a['name']的效率也很高;但是这个索引只对整个a起作用，如果输入只有a的一部分，仍然需要遍历整个a，去提取出对应的数据，因此a[1]['name']比a['name'][1]的效率差很多。

实例

基于numpy过滤抖动与填补

时间序列数据经常会发现两种情况：一种是抖得特别厉害，说明数据不稳定不可信，支撑这个结果的数据量不够;另一种是一动不动的一条直线，这往往是算法填充出来的默认值，不是实际值。这些数据对于挖掘来说是噪音，应该过滤掉。我们使用numpy来完成这个任务。抖动的特点是频繁跳动，即一阶差分有很多值绝对值比0大很多，那么我们将这些跳动的点抓出来，统计下这些点之间的区间长度，如果区间长度过小，认为是抖动过多。填补的特点是数值长期不变，即一阶差分有很多值为0，那么我们统计一下连续为0的区间长度分布，如果区间长度过长，比如连续填补了1小时，或者出现多个填补了30分钟的区间，我们认为是填补过多。

我们需要对跳点进行定义：一阶差分的绝对值超过dev_thresh，一阶差分/max(基准1，基准2)的绝对值超过ratio_thresh。 

 
 
    
  
  
  1. def jump(speed_array, dev_thresh, ratio_thresh): 
    2. 
  
  
  2.   diff_array = np.diff(speed_array, axis=0) 
    3. 
  
  
  3.   diff_array = diff_array.astype(np.float64) 
    4. 
  
  
  4.   ratio_array = diff_array/np.maxium(speed_array[:-1], speed_array[1:]) 
    5. 
  
  
  5.   ret_array = np.zeros(diff_array.size, dtype=np.int8) 
    6. 
  
  
  6.   for i in range(diff_array.size): 
    7. 
  
  
  7.     if abs(diff_array[i]) > diff_thresh and abs(ratio_array[i]) > ratio_thresh: 
    8. 
  
  
  8.       ret_array[i] = 1 
    9. 
  
  
  9.   return ret_array 
    10. 
  
  
  10.  
    11. 
  
  
  11. def interval(jump_array): 
    12. 
  
  
  12.   jump_idx = np.array([0] + [i for i,x in enumerate(jump_array) if x != 0] + [jump_array.size]) 
    13. 
  
  
  13.   interval_size = np.diff(jump_idx) 
    14. 
  
  
  14.   return interval_size 
    15. 
  
  
  15.  
    16. 
  
  
  16. def is_jump_too_much(interval_size): 
    17. 
  
  
  17.   flag = 0 
    18. 
  
  
  18.   if np.mean(interval_size) <= 10 or np.max(interval_size) <= 30: 
    19. 
  
  
  19.       flag = 1 
    20. 
  
  
  20.   return flag 
    21. 
  
  
  21.  
    22. 
  
  
  22. def is_fill_too_much(interval_size): 
    23. 
  
  
  23.   flag = 0 
    24. 
  
  
  24.   bin_array = np.bincount(interval_size) 
    25. 
  
  
  25.   if ( len(bin_array) >= 30 or  
    26. 
  
  
  26.        ( len(bin_array) >= 11 and np.sum(bin_array[10:]) >= 4 ) or  
    27. 
  
  
  27.        ( len(bin_array) >= 7 and np.sum(bin_array[6:]) >= 20 ) 
    28. 
  
  
  28.       ): 
    29. 
  
  
  29.       flag = 1 
    30. 
  
  
  30.   return flag  
    31.

基于numpy的局部趋势拟合

用线性回归可以得到时间序列的趋势。 

 
 
    
  
  
  1. def get_ts_trend(ts_array): 
    2. 
  
  
  2.   x = np.arange(0, len(ts_array), 1) 
    3. 
  
  
  3.   y = ts_array 
    4. 
  
  
  4.   A = np.vstack([x, np.ones(len(x))]).T 
    5. 
  
  
  5.   m, c = np.linalg.lstsq(A, y)[0] 
    6. 
  
  
  6.   return m  
    7.

堵点判别

交通数据比较复杂，不纯粹是时间序列问题，而是时空数据，需要同时考虑时间关系和空间关系。本节介绍一个经典特征的提取：堵点判别。

假设我们空间上有5个link，上游2个，自身1个，下游2个;观察5个时间点的拥堵状态。判断当前link是不是堵点——即自身是拓扑中***个发生拥堵的点;发生拥堵后，拥堵是扩散的。 

 
 
    
  
  
  1. def detect_congest_point(congest_array): 
    2. 
  
  
  2.     first_congest_flag = False 
    3. 
  
  
  3.     disperse_congest_flag = True 
    4. 
  
  
  4.     idx = np.where(congest_array == 1) 
    5. 
  
  
  5.     if idx[1][0] == congest_array.shape[1]/2: 
    6. 
  
  
  6.         first_congest_flag = True 
    7. 
  
  
  7.      
    8. 
  
  
  8.     disperse_dict = {} 
    9. 
  
  
  9.     for k in range(len(idx[0])): 
    10. 
  
  
  10.         if disperse_dict.has_key(idx[0][k]): 
    11. 
  
  
  11.             disperse_dict[idx[0][k]].append(idx[1][k]) 
    12. 
  
  
  12.         else: 
    13. 
  
  
  13.             disperse_dict[idx[0][k]]= [idx[1][k]]
    14. 
  
  
  14.     sorted_disperse_list = sorted(disperse_dict.iteritems(), key=lambda d:d[0]) 
    15. 
  
  
  15.     for i in range(1, len(sorted_disperse_list)): 
    16. 
  
  
  16.         if not set(sorted_disperse_list[i-1][1]) <= set(sorted_disperse_list[i][1]): 
    17. 
  
  
  17.             disperse_congest_flag = False 
    18. 
  
  
  18.      
    19. 
  
  
  19.     return first_congest_flag and disperse_congest_flag  
    20.

关于作者：丹追兵：数据分析师一枚，编程语言python和R，使用Spark、Hadoop、Storm、ODPS。本文出自丹追兵的pytrafficR专栏，转载请注明作者与出处：https://segmentfault.com/blog...

文章标题：numpy：python数据领域的功臣
文章地址：http://www.shufengxianlan.com/qtweb/news10/378360.html

网站建设、网络推广公司-创新互联，是专注品牌与效果的网站制作，网络营销seo公司；服务项目有等

广告

声明：本网站发布的内容（图片、视频和文字）以用户投稿、用户转载内容为主，如果涉及侵权请尽快告知，我们将会在第一时间删除。文章观点不代表本网站立场，如需处理请联系客服。电话：028-86922220；邮箱：631063699@qq.com。内容未经允许不得转载，或转载时需注明来源： 创新互联

猜你还喜欢下面的内容

自适应网站知识

分类信息网站

重庆服务器托管　　 　成都全网营销推广　　 　雅安做网站　　 　成都棕树机房　　 　背发光字　　 　成都移动服务器托管　　 　艾瑞格衣柜　　 　成都广告招牌　　 　成都学生服定制　　 　重庆企业网站建设　　 　小程序开发　　 　成都app开发公司　　 　成都微信开发　　 　专业网站设计　　 　成都展柜制作　　 　仿古建筑设计　　 　成都响应式网站建设公司　　 　彭州网站建设　　 　郫县封阳台　　 　雅安电信机房