100 numpy exercises 后50

https://github.com/rougier/numpy-100/blob/master/100%20Numpy%20exercises.md

51.新建一个向量表示（x，y)和(r,g,b)

z = np.zeros(10, [('position', [('x', float, 1),('y', float, 1)]), ('color', [('r', float, 1), ('g', float, 1), ('b', float, 1)])])

52.在一个10*2的矩阵求出两点间距离

z = np.random.random((10,2)) # atleast_2d查看至少2维的数组 x, y = np.atleast_2d(z[:,0]), np.atleast_2d(z[:,1]) D = np.sqrt((x-x.T)**2 + (y-y.T)**2) #更快的版本 import scipy import scipy.spatial z = np.random.random((10, 2)) D = scipy.spatial.distance.cdist(z, z)

53.怎样原地把float的数组转换成int

z = np.arange(10, dtype=np.float32) z = z.astype(np.int.32, copy=False)

54.怎样从文件读取(有问题)

from io import StringIO #假设文件读入IO s = StringIO("""1,2,3,4,5\n 6, , ,7,8\n , ,9,10,11\n""") z = np.genfromtxt(s, delimiter=",", dtype=np.int)

55.和numpy枚举相同的做法

z = np.arange(9).reshape(3,3) for index, value = in np.ndenumerate(z) : print(index, value) for index in np.ndindex(z.shape) : print(index, z[index])

56.生成一个2维高斯分布

x, y = np.meshgrid(np.linspace(-1,1,10), np.linspace(-1,1,10)) D = np.sqrt(x**x + y**y) sigma, mu = 1.0, 0.0 G = np.exp(-((D-mu)**2 / (2.0 *sigma**2)))

57.怎样在2维数组随机替换p个元素

n = 10 p = 3 z = np.zeros((n, n)) np.put(z, np.random.choice(range(n*n), p, replace=False), 1)

58.每行减平均数

y = x - x.mean(axis=1, keepdimas=True)

59.以第n列排序为重新排序矩阵

z = np.random.randint(0,10,(3,3)) z[z[:, 1].argsort()]

60.怎样发现给定的二维矩阵有空列

z = np.random.randint(0,3,(3,10)) print(~z.any(axis=0).any())

61.给定数组里的值找到最近的值

z = np.random.uniform(0,1,10) x = 0.5 m = z.flat[np.abs(z-x).argmin()]

62.如何用生成器计算（1，3）和（3，1）数组的和

a = np.arange(3).reshape(3,1) b = np.arange(3).reshape(1,3) it = np.nditer([a, b, None]) for x, y, z in it: z[...] = x + y print(it.operands[2])

63.新建一个数组类有数组属性

class NameArray(np.ndarray) : def __new__(cls, array, name="no name") : obj = np.asarray(array).view(cls) obj.name = name return obj def __array__finalize__(self, obj) : if obj is None : return self.info = getattr(obj, 'name', "no name") z = NameArray(np.arange(10), "range_10") print(z.name)

64.新建一个向量，如何用第二个向量做下标向第一个元素加一

z = np.ones(10) I = np.random.randint(0, len(z), 20) z += np.bincount(I, minlength=len(z)) #bincount 返回每个数字出现次数

65.根据I，累计X到F

X = [1,2,3,4,5,6] I = [1,3,9,3,4,1] F = np.bincount(I, X) #相当于把X当权重付给I，然后计数

66.给定一个(w,h,3)的图片，计算不一样的颜色个数

w, h = 16, 16 I = np.random,randint(0,2,(h,w,3)).astype(np.ubyte) F = I[...,0]*256*256 + I[...,1]*256 + I[...,2] n = len(np.unique(F))

67.给定一个四维的数组，如果一次性求最后两维的总和

a = np.random.randint(0,10,(3,4,3,4)) #把最后两维前的保留，-1为根据需要压缩 sum = a.reshape(a.shape[:-2] + (-1,)).sum(axis=-1)

68.给定一个一维向量D，怎样用一个一样大小的向量S描述子集的索引计算D子集的平均值

例如 D=[0.1, 0.2, 0.3] S = [0, 1, 0] 即索引值为0的平均值为 (0.1+0.3)/2=0.2

D = np.random.uniform(0,1,100) S = np.random.randint(0,10,100) D_sums = np.bincount(S, weight=D) D_count = np.bincount(S) D_means = D_sums / D_counts

69.怎样获得矩阵乘法结果的对角线

a = np.random.uniform(0,1,(5,5)) b = np.random.uniform(0,1,(5,5)) #慢版本 np.diag(np.dot(a, b)) #快版本 np.sum(a * b.T, axis=1) #更快版本，einsum貌似可以根据字符串来描述加的方法 np.einsum("ij, ji->i", a, b)

70.给定一个[1,2,3,4]，怎样新建一个新的向量两元素间有连续3个零

z = np.array([1,2,3,4,5]) nz = 3 z0 = np.zeros(len(z) + (len(z)-1)*(nz)) z0[::nz+1] = z

71.给定一个(5,5,3)的数组，怎样和一个(5,5)相乘

a = np.ones((5,5,3)) b = 2*np.ones((5,5)) a*b[:,:,None]

72.怎样交换两行

a = np.arange(25).reshape(5,5) # 1赋值0，0赋值1。[0,1]代表两行 a[[0,1]] = a[[1,0]]

73.给定10个3元组来表述10个三角形，线段两两组合

faces = np.random.randint(0,100,(10,3)) # roll 滚动 repeat [1,2] => [1,1,2,2] F = np.roll(faces.repeat(2, axis=1), -1, axis=-1) F = F.reshape(len(F)*3, 2) F = np.sort(F, axis=1) G = F.view(dtype=[('p0', F.dtype), ('p1', F.dtype)]) G = np.unique(G) # 先重复[1,2,3] => [1,1,2,2,3,3] # 再滚动 [1,1,2,2,3,3] => [1,2,2,3,3,1] # 这样就生成[1,2],[2,3],[3,1]

74.给定一个bincount数组C，怎样生成一个A，np.bincount == C

c = np.bincount([1,1,2,3,4,4,6]) [0,1,2,3,4,5,6] 0重复c[0]次 1重复c[1]次 a = np.repeat(np.arange(len(c)), c)

75.用计算滑窗内的平均值

def moving_average(a, n=3) : # cumsum 累计和 ret = np.cumsum(a, dtype=float) ret[n:] = ret[n:] - ret[:-n] return ret[n-1:] / n z = np.arange(20) moving_average(z, n=3)

76.给定一个一维向量，构建一个新二维向量，第一行z[0] z[1] z[2]，后面每行往后移一位

from numpy.lib import stride_tricks def rolling(a, window) : shape = (a.size - window + 1, window) # srtides 代表跳过一行需要多少字节，跳过一列需要多少字节 strides = (a.itemsize, a.itemsize) return stride_tricks.as_strided(a, shape=shapem strides=strides) z = rolling(np.arange(10), 3)

77.布尔去反，或原地改浮点数符号

z = np.random.randint(0,2,100) np.logical_not(z, out=z) z = np.random.uniform(-1.0,1.0,100) np.negative(z, out=z)

78.给定两个点集合P0, P1描述线和一个点p，怎样比较p点到每条线的距离

def distance(P0, P1, p) : T = P1 - P0 L = (T**2).sum(axis=1) U = -((P0[:,0]-p[...,0]) * T[:,0] + (P0[:,1]-p[...,1]) * T[:,1]) / L U = U.reshape(len(U), 1) D = P0 + U*T - p return np.sqrt((D**2).sum(axis=1)) P0 = np.random.uniform(-10,10,(10,2)) P1 = np.random.uniform(-10,10,(10,2)) p = np.random.uniform(-10,10,(1,2)) ditance(P0,P1,p)

79.每个点到每条线的距离

P0 = np.random.uniform(-10,10,(10,2)) P1 = np.random.uniform(-10,10(10,2)) p = np.random.uniform(-10,10,(10,2)) np.array([distance(P0,P1,p_i) for p_i in p])

80.给定一个随机矩阵，写一个函数提取一个子集合，开始点位给定位置，其余用fill元素填充

Z = np.random.randint(0,10,(10,10)) shape = (5,5) fill = 0 position = (1,1) R = np.ones(shape, dtype=Z.dtype)*fill P = np.array(list(position)).astype(int) Rs = np.array(list(R.shape)).astype(int) Zs = np.array(list(Z.shape)).astype(int) R_start = np.zeros((len(shape),)).astype(int) R_stop = np.array(list(shape)).astype(int) Z_start = (P-Rs//2) Z_stop = (P+Rs//2)+Rs%2 R_start = (R_start - np.minimum(Z_start,0)).tolist() Z_start = (np.maximum(Z_start,0)).tolist() R_stop = np.maximum(R_start, (R_stop - np.maximum(Z_stop-Zs,0))).tolist() Z_stop = (np.minimum(Z_stop,Zs)).tolist() r = [slice(start,stop) for start,stop in zip(R_start,R_stop)] z = [slice(start,stop) for start,stop in zip(Z_start,Z_stop)] R[r] = Z[z] print(Z) print(R)

81.给定一个z=[1:15]，怎样生成z[0] z[1] z[2] z[3],每行左移1的新矩阵

z = np.arange(1,15,dtype=np.unit32) r = stride_tricks.as_strided(z, (11,4), (4,4))

82.计算矩阵的秩

z = np.random.uniform(0,1,(10,10)) # 分解矩阵z = u * np.diag(s) * v u, s, v = np.linalg.svd(z) # 大于0 rank = np.sum(s > 1e-10)

83.找众数

z = np.random.randint(0,10,50) np.bincount(z).argmax()

84.在10*10的矩阵内用3*3的滑窗新建一个向量

z = np.random.randint(0,5,(10,10)) n = 3 i = 1+(z.shape[0]-3) j = 1+(z.shape[1]-3) c = stride_tricks.as_strided(z, shape(i,j,n,n), strides=z.strides+z.strides)

85.新建一个2维数组，z[i,j]=z[j,i]，修改一个值能自动更改

class Symetric(np.ndarray): def __setitem__(self, index, value): i,j = index super(Symetric, self).__setitem__((i,j), value) super(Symetric, self).__setitem__((j,i), value) def symetric(Z): return np.asarray(Z + Z.T - np.diag(Z.diagonal())).view(Symetric) S = symetric(np.random.randint(0,10,(5,5))) S[2,3] = 42

86.假设p个(n,n)的矩阵，p个(n,1)的向量，怎样一次性计算和

p, n = 10, 20 M = np.ones((p,n,n)) V = np.ones((p,n,1)) # 0跟0，2跟1相乘 S = np.tensordot(M, V, axes=[[0,2], [0,1]])

87.给定一个16*16的数组，怎样获得模块和（模块大小4*4）

z = np.ones((16,16)) k = 4 s = np.reduceat(np.add.reduceat(z, np.arange(0,z.shape[0], k),axis=0),np.arange(0,z.shape[1],k), axis=1)

88. Game of Life

def iterate(Z): # Count neighbours N = (Z[0:-2,0:-2] + Z[0:-2,1:-1] + Z[0:-2,2:] + Z[1:-1,0:-2] + Z[1:-1,2:] + Z[2: ,0:-2] + Z[2: ,1:-1] + Z[2: ,2:]) # Apply rules birth = (N==3) & (Z[1:-1,1:-1]==0) survive = ((N==2) | (N==3)) & (Z[1:-1,1:-1]==1) Z[...] = 0 Z[1:-1,1:-1][birth | survive] = 1 return Z Z = np.random.randint(0,2,(50,50)) for i in range(100): Z = iterate(Z) print(Z)

89. N大数

z = np.arange(10000) np.random.shuffle(z) n = 5 #slow z[np.argsort(z)[-n:]] #fast argpartition 小于n的放在前面，大于的放在后面 z[np.argpartion(-z,n)[:n]]

90.给定一个随机向量，新建一个笛卡尔组合

def cartesian(arrays): arrays = [np.asarray(a) for a in arrays] shape = (len(x) for x in arrays) ix = np.indices(shape, dtype=int) ix = ix.reshape(len(arrays), -1).T for n, arr in enumerate(arrays): ix[:, n] = arrays[n][ix[:, n]] return ix print (cartesian(([1, 2, 3], [4, 5], [6, 7])))

91.怎样从一个数组新建一个记录数组

z = np.array([("Hello" ,2.5, 3),("World", 3.6, 2)]) r = np.core.records.fromarray(z.T, names='col1, col2, col3', formats='S8,f8,i8')

92.给定一个很大的z，怎样求立方

x = np.random.rand(5e7) %timeit np.power(x,3) %timeit x*x*x %time np.einsum('i,i,i->i', x,x,x)

93.给定A(8,3) B(2,2), 怎样找到A的行包括B每一行的元素不管顺序

A = np.random.randint(0,5,(8,3)) B = np.random.randint(0,5,(2,2)) C = (A[..., np.newaxis, np.newaxis] == B) rows = (C.sum(axis=(1,2,3)) >= B.shape[1]).nonzero()[0] print(rows)

94.给定一个10*3的矩阵，提取有不同值得行

z = np.random.randint(0,5,(10,3)) e = np.logical_and.reduce(z[:,1:]==z[:,:,-1], axis=-1) u = z[~e]

95.转换一个整数的向量，变成二进制表达的矩阵

I = np.array([0, 1, 2, 3, 15, 16, 32, 64, 128]) B = ((I.reshape(-1,1) & (2**np.arange(8))) != 0).astype(int) print(B[:,::-1]) I = np.array([0, 1, 2, 3, 15, 16, 32, 64, 128], dtype=np.uint8) print(np.unpackbits(I[:, np.newaxis], axis=1))

96.给定一个2维数组，怎样提出唯一的行

Z = np.random.randint(0,2,(6,3)) T = np.ascontiguousarray(Z).view(np.dtype((np.void, Z.dtype.itemsize * Z.shape[1]))) _, idx = np.unique(T, return_index=True) uZ = Z[idx] print(uZ)

97.给定两个A，B, 用einsum写出inner outer sum mul功能

a = np.random.uniform(0,1,10) b = np.random.uniform(0,1,10) np.einsum('i->',A) #np.sum(a) np.enisum('i,i->i',a,b) #a * b np.enisum('i,i',a,b) #np.inner(a,b) np.einsum('i,j',a,b) #np.outer(a,b)

98.给定一个用(x,y)描述的路，怎样等距抽样

phi = np.arange(0, 10*np.pi, 0.1) a = 1 x = a*phi*np.cos(phi) y = a*phi*np.sin(phi) dr = (np.diff(x)**2 + np.diff(y)**2)**.5 # segment lengths r = np.zeros_like(x) r[1:] = np.cumsum(dr) # integrate path r_int = np.linspace(0, r.max(), 200) # regular spaced path x_int = np.interp(r_int, r, x) # integrate path y_int = np.interp(r_int, r, y)

99.给定一个整数n和二维数组X，从x中选出符合n分布的行

X = np.asarray([[1.0, 0.0, 3.0, 8.0], [2.0, 0.0, 1.0, 1.0], [1.5, 2.5, 1.0, 0.0]]) n = 4 M = np.logical_and.reduce(np.mod(X, 1) == 0, axis=-1) M &= (X.sum(axis=-1) == n) print(X[M])

100.重新抽样，计算每个样本的平均值，计算平均值的百分位数的值

X = np.random.randn(100) # random 1D array N = 1000 # number of bootstrap samples idx = np.random.randint(0, X.size, (N, X.size)) means = X[idx].mean(axis=1) confint = np.percentile(means, [2.5, 97.5]) print(confint)