python中的list和numpy中的矩阵分析

Author : Jasper Yang
School : Bupt

preface

由于之前在做GIbbsLDA++的源码学习，并且将其c++的源码翻译成了python的版本。后来有朋友用我的实现在大数据量的情况下内存跑崩溃了，仔细去网上一查，才发现了python中的list的实现方式是一种很泛化的面对各种类型的数据结构，这个结构用来做二位数组比numpy中的narrays需要占用更多更过的内存，后面我会详细分析。（但是我发现好像并不是如此）

正文

他两在使用上的区别

	l=[1,2L,3.0,'a'] # 这是一个放了不同数据类型的list
	  
	a=np.array([1,2,3]) # list -> array   
	b=np.array([[1,2,3],[4,5,6]])   
	c=list(a)   # array到list的转换  
	print l
	print a,np.shape(a)
	print b,np.shape(b)  
	print c,np.shape(c)

程序输出如下：

[1, 2L, 3.0, 'a']
[1 2 3] (3,)
[[1 2 3]
[4 5 6]] (2, 3)
[1, 2, 3] (3,)

可以看到list都有逗号来隔开，因为list中为每一个数据分配了一个指针，正因为这种实现方式所以可以面对不同的数据类型。但因如此，l用了4个指针和4个数据的内存空间。

	# numpy 定义矩阵
	array=([[1,2,3],
        [4,5,6],
        [7,8,9]])
        
	print array[1]
	print array[1][1]

5
[4, 5, 6]

下面我想看看tuple，list，array占用内存的情况，用到了sys里的函数

	import sys
	a = (1,2,3,4,5)
	b = [1,2,3,4,5]
	c = np.array([1,2,3,4,5])
	d = np.array((1,2,3,4,5))
	print(sys.getsizeof(a))
	print(sys.getsizeof(b))
	print(sys.getsizeof(c))
	print(sys.getsizeof(d))

output:

bytes 为单位

96
112
136
136
这么一看，有点太奇怪了，怎么会是用了array之后比list占用的内存还要大~？

再试试

	s1 = ()
	s2 = []
	s3 = zero(0) # 创建0矩阵，用array的方法
	print(sys.getsizeof(s1))
	print(sys.getsizeof(s2))
	print(sys.getsizeof(s3))

output:

56
72
96

这里发现tuple初始化费用最小，list次之，array的最大，那么说明了array并没有比占用更小的内存空间，但是，没错，凡事都有但是。array在矩阵操作的性能上必然是大大大于list的，见网上总结如下。

Python 的 list 是动态类型，可以包含不同类型的元素，所以没有支持诸如点乘等数学函数，因为要为 list 实现这些操作会牺牲性能。
Numpy 数组是静态类型并且齐次。元素类型在数组创建的时候就已经确定了。
Numpy 数组节约内存。(这点我不能赞同，难道是getsizeof统计不到指针？)
由于是静态类型，对其数学操作函数（如矩阵乘法，矩阵加法）的实现可以使用 C 或者 Fortran 完成。

而且通过array的itemsize属性发现每次增加数据，array的开销是更多的，也就是说数组越大，array占用的内存比list多越多。

我从Python中优化NumPy包使用性能的教程这里看到了一个很有意思的解释。

为什么NumPy数组如此高效？

一个NumPy数组基本上是由元数据（维数、形状、数据类型等）和实际数据构成。数据存储在一个均匀连续的内存块中，该内存在系统内存（随机存取存储器，或RAM）的一个特定地址处，被称为数据缓冲区。这是和list等纯Python结构的主要区别，list的元素在系统内存中是分散存储的。这是使NumPy数组如此高效的决定性因素。

为什么这会如此重要？主要原因是：

低级语言比如C，可以很高效的实现数组计算（NumPy的很大一部分实际上是用C编写）。例如，知道了内存块地址和数据类型，数组计算只是简单遍历其中所有的元素。但在Python中使用list实现，会有很大的开销。(这里也是提到了很大开销，我很想知道怎么才能看到那很大的开销在哪？)
内存访问模式中的空间位置访问会产生显著地性能提高，尤其要感谢CPU缓存。事实上，缓存将字节块从RAM加载到CPU寄存器。然后相邻元素就能高效地被加载了（顺序位置，或引用位置）。
数据元素连续地存储在内存中，所以NumPy可以利用现代CPU的矢量化指令，像英特尔的SSE和AVX，AMD的XOP等。例如，为了作为CPU指令实现的矢量化算术计算，可以加载在128,256或512位寄存器中的多个连续的浮点数。

此外，说一下这样一个事实：NumPy可以通过Intel Math Kernel Library (MKL)与高度优化的线性代数库相连，比如BLAS和LAPACK。NumPy中一些特定的矩阵计算也可能是多线程，充分利用了现代多核处理器的优势。

总之，将数据存储在一个连续的内存块中，根据内存访问模式，CPU缓存和矢量化指令，可以确保以最佳方式使用现代CPU的体系结构。

在逛知乎里，我又发现了很多关于为什么numpy这么快的讨论，很有意思。

numpy的许多函数不仅是用C实现了，还使用了BLAS（一般Windows下link到MKL的，Linux下link到OpenBLAS）。基本上那些BLAS实现在每种操作上都进行了高度优化，例如使用AVX向量指令集，甚至能比你自己用C实现快上许多，更不要说和用Python实现的比。。

作者：Zhipeng
链接：https://www.zhihu.com/question/30823702/answer/49696394
来源：知乎
著作权归作者所有。商业转载请联系作者获得授权，非商业转载请注明出处。

numpy底层使用BLAS做向量，矩阵运算。像求平均值这种vector operation，很容易使用multi-threading或者vectorization来加速。比如MKL就有很多优化。

p.s. 一个小小的experiment可以看看（摘抄--大大的橙子）

	from sys import getsizeof
	   
	class A(object): pass
	class B: pass
	
	for x in (None, 1, 1L, 1.2, 'c', [], (), {}, set(), B, B(), A, A()):
	  print "{0:20s}\t{1:d}".format(type(x).__name__, sys.getsizeof(x))
	
	
	
	NoneType                16
	int                     24
	long                    28
	float                   24
	str                     34
	list                    64
	tuple                   48
	dict                    272
	set                     224
	classobj                96
	instance                64
	type                    896
	A                       56

一份numpy常用函数总结

arange

create a range

x = arange(0, 10, 1) # arguments: start, stop, step
x

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

x = arange(-1, 1, 0.1)
x
output:

array([ -1.00000000e+00, -9.00000000e-01, -8.00000000e-01,
-7.00000000e-01, -6.00000000e-01, -5.00000000e-01,
-4.00000000e-01, -3.00000000e-01, -2.00000000e-01,
-1.00000000e-01, -2.22044605e-16, 1.00000000e-01,
2.00000000e-01, 3.00000000e-01, 4.00000000e-01,
5.00000000e-01, 6.00000000e-01, 7.00000000e-01,
8.00000000e-01, 9.00000000e-01])
mgrid

x, y = mgrid[0:5, 0:5]
print(x)
print(y)
output:

array([[0, 0, 0, 0, 0],
[1, 1, 1, 1, 1],
[2, 2, 2, 2, 2],
[3, 3, 3, 3, 3],
[4, 4, 4, 4, 4]])

array([[0, 1, 2, 3, 4],
[0, 1, 2, 3, 4],
[0, 1, 2, 3, 4],
[0, 1, 2, 3, 4],
[0, 1, 2, 3, 4]])
linspace 和 logspace

linspace(0, 10, 20)
logspace(0, 10, 10, base=e)
output:

array([ 0. , 0.52631579, 1.05263158, 1.57894737,
2.10526316, 2.63157895, 3.15789474, 3.68421053,
4.21052632, 4.73684211, 5.26315789, 5.78947368,
6.31578947, 6.84210526, 7.36842105, 7.89473684,
8.42105263, 8.94736842, 9.47368421, 10. ])

array([ 1.00000000e+00, 3.03773178e+00, 9.22781435e+00,
2.80316249e+01, 8.51525577e+01, 2.58670631e+02,
7.85771994e+02, 2.38696456e+03, 7.25095809e+03,
2.20264658e+04])
random data

random.rand(4,4)
random.randn(4,4)
output:

array([[ 0.2659394 , 0.59798318, 0.16067613, 0.21704906],
[ 0.98081628, 0.10926225, 0.86342814, 0.83130784],
[ 0.83563301, 0.18313372, 0.73531414, 0.67861801],
[ 0.38322499, 0.02943103, 0.43471297, 0.87099786]])

array([[ 0.18834898, 0.28862928, 0.58291415, 0.57712703],
[ 0.75071525, -0.39247518, -0.35748584, -0.42839121],
[-1.12789581, -0.67267265, -0.05525579, -0.89712592],
[-0.14731484, -0.72237449, -0.16594984, 0.62914291]])
diag

diag([1,2,3])
diag([1,2,3], k=1)
output

array([[1, 0, 0],
[0, 2, 0],
[0, 0, 3]])

array([[0, 1, 0, 0],
[0, 0, 2, 0],
[0, 0, 0, 3],
[0, 0, 0, 0]])
zeros 和ones

zeros((2,2))
ones((2,2))

output

array([[ 0., 0.],
[ 0., 0.]])

array([[ 1., 1.],
[ 1., 1.]])

做了上面这份总结之后使用各种数据类型的创建就不用总是查找了，虽然都很简单。

numpy的索引很有意思，基本分成了以下三类

花式索引（arange）
索引列表（take） -----> 最快
布尔掩码（compress）

一些优化，在 stackoverflow 上看的，为了更好的记住这些trick我决定记在这篇博客里。

	import numpy
	x = numpy.array([0] * 1000000)
	for i in range(1,len(x)):
	  x[i] = x[i-1] + i
	  
	a[0] += 1232234234234324353453453
	OverflowError: Python int too large to convert to C long

可以变成

import numpy as np
x = np.arange(1000000).cumsum()
同时上面出错数字超出了范围可以如下解决

a = np.array([0], dtype=object)
a[0] += 1232234234234324353453453

使用 cumsum() 这样做叫做 vectorized operations 。

全文完，希望能帮到你学到姿势~
另外，给大家推荐一篇我的关于python中list的实现的博客 ---> python中list的实现

补充 -----> 浅拷贝和深拷贝以及广播

numpy中的矩阵属于浅拷贝

a = np.zeros((2,2))
b = a
b[1][1] = 1
print(a)

a的值会被改变，这是因为这是浅拷贝，b获得的是a的指针，数据区域还是一样的

list则是深拷贝

a = [1,2,3,4]
b = a
b[2] = 6
print(a)

a的值不会被改变，这是因为深拷贝，将数据区域包括指针都复制了一遍过来

还有一个numpy矩阵的广播，这个是一个很神奇的使用方式

a = np.arange(0, 50, 10).reshape(-1, 1)
b = np.arange(0, 4)
print a
print b

这时的输出

[[ 0]
[10]
[20]
[30]
[40]]
[0 1 2 3]

当我们使用广播

np.add(a,b)

输出变成了

array([[ 0, 1, 2, 3, 4],
[10, 11, 12, 13, 14],
[20, 21, 22, 23, 24],
[30, 31, 32, 33, 34],
[40, 41, 42, 43, 44],
[50, 51, 52, 53, 54]])

是不是很神奇呢，用法一目了然。 see also --> numpy.broadcast

广播用以描述numpy中对两个形状不同的阵列进行数学计算的处理机制。较小的阵列“广播”到较大阵列相同的形状尺度上，使它们对等以可以进行数学计算。广播提供了一种向量化阵列的操作方式，因此Python不需要像C一样循环。广播操作不需要数据复制，通常执行效率非常高。然而，有时广播是个坏主意，可能会导致内存浪费以致计算减慢。

paper done 2017/4/20

posted @ 2017-04-26 21:39 japeryang 阅读(1466) 评论(0) 收藏举报

japeryang

python中的list和numpy中的矩阵分析