前言numpy对python意义重大，在数据分析和机器学习领域为python做出了巨大的贡献。pandas、matplotlib、sklearn等在python中常用于数据分析或机器学习的库，都需要基于numpy构建。毫不夸张地说，如果没有numpy，python在今天的数据分析和机器学习领域只能捉襟见肘。什么是好的数据分析语言数据分析的大部分数据都是二维表。一个好的数据分析语言首先需要能够将二维表直接存储在数据结构中，然后需要配备成熟的类SQL数据操作接口，最后需要一套易于上手的-使用可视化工具。R语言就是一个很好的例子：使用内置的data.frame结构进行数据存储；data.frame本身提供了足够的数据操作能力，而dplyr、tidyr、data.table、plyr、reshape2等其他库提供了更易使用更高效的数据操作能力；在绘图方面，除了基本的绘图功能外，还提供了ggplot2等优雅的绘图语言，并且还通过htmlwidget库与javascript的各种绘图库建立了紧密的联系，让可视化的动态显示效果更上一层楼更远。Excel也是一个很好的例子。它有像单元格一样灵活的结构来支持数据存储，有大量的函数来实现灵活的操作，还有强大的绘图系统。Python目前在数据分析领域具有相当的能力，包括pandas库实现的DataFrame结构，pandas自身提供的数据操作能力，matplotlib提供的数据可视化能力，这些都离不开numpy库。什么是好的机器学习语言一般来说，好的机器学习语言在数据分析方面一定也很受欢迎，因为数据分析往往是机器学习的基础。但是机器学习的要求更高，因为在模型训练阶段往往需要进行复杂的参数估计操作，所以需要语言具有很强的科学计算能力。科学计算能力的核心是矩阵计算能力。关于矩阵计算能力，这篇文章对各种语言都有很好的比较。如果没有numpy，python内部只能用list或者array来表示矩阵。如果用列表表示[1,2,3]，由于列表的元素可以是任意对象，列表中存放的是对象的指针，所以需要有3个指针和3个整型对象，是浪费内存和CPU计算时间。Python的数组不同于列表。直接存值，类似于C语言的一维数组，但不支持多维，表达形式很简单。编写科学的计算算法是很困难的。numpy弥补了这些缺点。它提供的ndarray是一个存储单一数据类型的多维数组，使用预编译的C语言代码，性能也很不错。Python最新的机器学习库sklearn建立在numpy之上，为各种标准机器学习模型提供训练和预测接口。模型训练接口内部实现基于numpy库。例如，对于一个非常常见的线性回归模型，参数估计调用了numpy.linalg.lstsq函数。numpy的核心结构：ndarray。float32)ndarray是numpy的核心数据结构。我们来看看ndarray在内存中是如何存储的：数组的描述信息存储在一个数据结构中，这个数据结构指的是两个对象，一个是存储数据的存储区，一个是描述元素类型的dtype对象。数据存储区保存着数组中所有元素的二进制数据，dtype对象知道如何将元素的二进制数据转换成可用的值。数组的维度和大小等信息存储在ndarray数组对象的数据结构中。strides中存储的是当每个轴的下标都增加1时，数据存储区中的指针增加的字节数。例如图中的strides为12,4，即第0个轴的下标时加1，数据的地址增加12字节：即a[1,0]的地址比a[0,0]的地址高12字节，恰好是3单的总字节数-精度浮点数；当第一个轴的下标增加1时，数据的地址增加4个字节，正好是单精度浮点数的字节数。以下内容摘自Numpy官方文档Numpybasics，关于ndarray的索引方式。有以下几点需要记住：虽然x[0,2]=x0，但是前者比后者效率更高，因为后者是应用了***第一个索引之后，需要先创建一个临时数组，并且然后应用第二个索引最终找到目标值。切片操作不会触发复制操作，而是创建原始ndarray的视图；它们指向的内存是同一个区域，无论是修改原来的ndarray还是修改view，两者的值都会同时改变。索引数组和布尔索引返回副本，而不是视图。至于上面列出的不会触发复制操作的分片操作，我们再讨论。先看numpy的例子：再看R的例子：可以看出numpy和R对矩阵的切片操作有不同的设计理念：在R中切片操作会引起数据复制，而在numpy中不会。其实R的设计理念往往可以用一句话来概括：copyonmodify。一旦数据被修改，就会引起对内存的复制操作。这个操作需要很多时间，所以人们经常抱怨R使用内存而且很慢。因此，我们可以看到，numpy显然要非常重视这件事，根据场景设计不同的策略，而不是简单地采用R的一刀切的做法。当然这也带来了一些学习成本，需要对numpy足够熟悉，以免踩坑。在R社区中，copyonmodify的哲学也受到批评并试图改变。比如data.table和dplyr都是data.frame操作库。data.table的性能远高于dplyr。部分原因是data.table规避了copyonmodify的方式。结构化数组根据numpy的官方文档，定义结构化数组有四种方式。本文采用字典的方式，通过定义一个dtype对象来实现。需要指定的键值包括名称和格式。persontype=np.dtype({'names':['name','age','weight'],'formats':['S32','i','f']})a=np.array([("Zhang",32,75.5),("Wang",24,65.2)],dtype=persontype)下面用IPython的计时函数看看提取数据的效率：%timeita[1]%timeita['name']%timeita[1]['name']%timeita['name'][1]输出结果如下：Theslowestruntook46.83timeslongerthanthefastest.Thiscouldmeanthatanintermediateresultisbeingcached.1000000loops,bestof3:153nsperloopTheslowestruntook34.34timeslongerthanthefastest.Thiscouldmeanthatanintermediateresultisbeingcached.10000000loops,bestof3:174nsperloopTheslowestruntook13.00timeslongerthanthefastest.Thiscouldmeanthatanintermediateresultisbeingcached.1000000loops,bestof3:1.08μsperloop最慢的运行时间比最快的运行时间长9.84倍。这可能意味着一个中间结果被缓存了。1000000loops,bestof3:412nsperloop从上面的结果我们发现获取相同数据的操作有很多，不同操作的性能差别很大。我做了一个猜测，纯属猜测：numpy在构建结构化数组时，会把整个结构连续存储在一起，也就是按行存储，所以a[1]最快；但是为了保证提取列的效率，为a['name']建立了索引，所以a['name']的效率也很高；但是这个索引只对整个a起作用，如果输入的只是a的一部分，还是需要遍历整个a来提取对应的数据输出，所以a[1]['name']的效率要低很多比[‘名字’][1]。例子基于numpy过滤抖动填充时序数据，经常发现两种情况：一种是抖动特别厉害，说明数据不稳定不可靠，支撑这个结果的数据量不够；另一种是不动的直线。这通常是算法填充的默认值，而不是实际值。这些数据是挖矿的噪音，应该被过滤掉。我们使用numpy来完成这项任务。抖动的特点是频繁抖动，即很多一阶差分的绝对值远大于0，那么我们抓取这些抖动点，统计这些点之间的间隔长度。如果间隔长度太小，则认为抖动过大。填充的特点是值长期不变，即一阶差的很多值都是0，那么我们统计连续0的区间长度分布，如果区间长度太长，比如连续填充1小时，或者多次填充30分钟间隔，我们认为是过填充。我们需要定义跳跃点：一阶差分的绝对值超过dev_thresh，一阶差分/max(benchmark1,benchmark2)的绝对值超过ratio_thresh。defjump(speed_array,dev_thresh,ratio_thresh):diff_array=np.diff(speed_array,axis=0)diff_array=diff_array.astype(np.float64)ratio_array=diff_array/np.maxium(speed_array[:-1],speed_array[1:])ret_array=np.zeros(diff_array.size,dtype=np.int8)foriinrange(diff_array.size):ifabs(diff_array[i])>diff_threshandabs(ratio_array[i])>ratio_thresh:ret_array[i]=1returnret_arraydefinterval(jump_array):jump_idx=np.array([0]+[ifori,xinenumerate(jump_array)ifx!=0]+[jump_array.size])interval_size=np.diff(jump_idx)returninterval_sizedefis_jump_too_much(interval_size):flag=0ifnp.mean(interval_size)<=10ornp.max(interval_size)<=30:flag=1returnflagdefis_fill_too_much(interval_size):flag=0bin_array=np.bincount(interval_size)if(len(bin_array)>=30or(len(bin_array)>=11andnp.sum(bin_array[10:])>=4)or(len(bin_array)>=7andnp.sum(bin_array[6:])>=20)):flag=1returnflag基于numpy的局态模拟合线性回归可以得到时间顺序的趋势。defget_ts_trend(ts_array):x=np.arange(0,len(ts_array),1)y=ts_arrayA=np.vstack([x,np.ones(len(x))]).Tm,c=np.linalg.lstsq(A,y)[0]returnm交通数据识别拥堵点比较复杂。它不是纯粹的时间序列问题，而是时空数据。它需要同时考虑时间和空间关系。本节介绍一个经典特征的提取：blockingpointdiscrimination。假设我们的空间有5个链接，2个upstream，1个self，2个downstream；观察5个时间点的拥堵情况。判断当前链路是否为拥塞点——即拓扑中第一个拥塞点；拥塞发生后，拥塞蔓延。defdetect_congest_point（congest_array）：first_congest_flag=Falseddisperse_congest_flag=Trueidx=np.where（congest_array==1）ifidx[1][0]==congest_array.shape[1]/2：first_congest_flag=Truedisperse_dict={}forkinrange（len（idx[0])):ifdisperse_dict.has_key(idx[0][k]):disperse_dict[idx[0][k]].append(idx[1][k])else:disperse_dict[idx[0][k]]=[idx[1][k]]sorted_disperse_list=sorted(disperse_dict.iteritems(),key=lambdad:d[0])foriinrange(1,len(sorted_disperse_list)):ifnotset(sorted_disperse_list[i-1][1])<=set(sorted_disperse_list[i][1]):disperse_congest_flag=Falsereturnfirst_congest_flaganddisperse_congest_flag关于作者：DanChasingBing：数据分析师，编程语言python和R，使用Spark、Hadoop、Storm、ODPS。本文来自单追兵的pytrafficR专栏，转载请注明作者及出处：https://segmentfault.com/blog...
