Numpy100

1

import numpy as np

1
2

print(np.__version__)
np.show_config()

1
2

Z = np.zeros(10)
print(Z)

1
2

Z = np.zeros((10,10))
print("%d bytes" % (Z.size * Z.itemsize))

1

%run `python -c "import numpy; numpy.info(numpy.add)"`

1
2
3

Z = np.zeros(10)
Z[4] = 1
print(Z)

1
2

Z = np.arange(10,50)
print(Z)

1
2
3

Z = np.arange(50)
Z = Z[::-1]
print(Z)

1
2

Z = np.arange(9).reshape(3, 3)
print(Z)

1
2

nz = np.nonzero([1,2,0,0,4,0])
print(nz)

1
2

Z = np.eye(3)
print(Z)

1
2

Z = np.random.random((3,3,3))
print(Z)

1
2
3

Z = np.random.random((10,10))
Zmin, Zmax = Z.min(), Z.max()
print(Zmin, Zmax)

1
2
3

Z = np.random.random(30)
m = Z.mean()
print(m)

1
2
3

Z = np.ones((10,10))
Z[1:-1,1:-1] = 0
print(Z)

1
2
3
4
5
6
7
8

Z = np.ones((5,5))
Z = np.pad(Z, pad_width=1, mode='constant', constant_values=0)
print(Z)

# Using fancy indexing
Z[:, [0, -1]] = 0
Z[[0, -1], :] = 0
print(Z)

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17

0 * np.nan
# nan

np.nan == np.nan
# False

np.inf > np.nan
# False

np.nan - np.nan
# nan

np.nan in set([np.nan])
# True

0.3 == 3 * 0.1
# False

1
2
3
4
5
6
7
8

Z = np.diag(1+np.arange(4), k=-1)
print(Z)

# [[0 0 0 0 0]
#  [1 0 0 0 0]
#  [0 2 0 0 0]
#  [0 0 3 0 0]
#  [0 0 0 4 0]]

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13

Z = np.zeros((8,8),dtype=int)
Z[1::2,::2] = 1
Z[::2,1::2] = 1
print(Z)

# [[0 1 0 1 0 1 0 1]
#  [1 0 1 0 1 0 1 0]
#  [0 1 0 1 0 1 0 1]
#  [1 0 1 0 1 0 1 0]
#  [0 1 0 1 0 1 0 1]
#  [1 0 1 0 1 0 1 0]
#  [0 1 0 1 0 1 0 1]
#  [1 0 1 0 1 0 1 0]]

1
2

print(np.unravel_index(99, (6,7,8)))
# (np.int64(1), np.int64(5), np.int64(3))

1
2
3

# 通过重复 [[0,1],[1,0]] 矩阵来构造
Z = np.tile(np.array([[0,1],[1,0]]), (4,4))
print(Z)

1
2
3
4

Z = np.random.random((5,5))
# 直接减去均值，再除以标准差
Z = (Z-np.mean(Z)) / (np.std(Z))
print(Z)

1
2
3
4

color = np.dtype([("r", np.ubyte),
                  ("g", np.ubyte),
                  ("b", np.ubyte),
                  ("a", np.ubyte)])

1
2
3
4
5
6

Z = np.dot(np.ones((5,3)), np.ones((3,2)))
print(Z)

# Alternative solution, in Python 3.5 and above
Z = np.ones((5,3)) @ np.ones((3,2))
print(Z)

1
2
3
4
5

# Author: Evgeni Burovski

Z = np.arange(11)
Z[(3 < Z) & (Z < 8)] *= -1
print(Z)

1
2
3
4
5
6
7
8

# Author: Jake VanderPlas

print(sum(range(5),-1))
from numpy import *
print(sum(range(5),-1))

# 9，未导入 numpy 的 sum 之前，0+1+2+3+4-1=9
# 10，导入之后，-1 表示的是在某个轴上求和

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17

Z = np.array([1,2,3])
# 数组的元素级幂运算
Z**Z
# 位运算，先左移再右移，作用到每一个元素上
2 << Z >> 2
Z <- Z
# 复数单位 1j 乘以数组的每个元素
1j*Z
Z/1/1
Z<Z>Z

# [ 1  4 27]
# [1 2 4]
# [False False False]
# [0.+1.j 0.+2.j 0.+3.j]
# [1. 2. 3.]
# 报错

1
2
3
4
5
6
7

np.array(0) / np.array(0)
np.array(0) // np.array(0)
np.array([np.nan]).astype(int).astype(float)

# nan
# 0
# [-9.22337204e+18]

 1
 2
 3
 4
 5
 6
 7
 8
 9
10

# Author: Charles R Harris

Z = np.random.uniform(-10,+10,10)
print(np.copysign(np.ceil(np.abs(Z)), Z))
# 先对绝对值向上取整
# 再通过 copysign() 函数将符号设置为 Z 中对应元素的符号

# More readable but less efficient
print(np.where(Z>0, np.ceil(Z), np.floor(Z)))
# 通过 where 判断，如果元素大于零则向上取整，否则向下取整

1
2
3

Z1 = np.random.randint(0,10,10)
Z2 = np.random.randint(0,10,10)
print(np.intersect1d(Z1, Z2))

 1
 2
 3
 4
 5
 6
 7
 8
 9
10

# Suicide mode on
defaults = np.seterr(all="ignore")
Z = np.ones(1) / 0

# Back to sanity
_ = np.seterr(**defaults)

# Equivalently with a context manager
with np.errstate(all="ignore"):
    np.arange(3) / 0

1
2

np.sqrt(-1) == np.emath.sqrt(-1)
# False

1
2
3

yesterday = np.datetime64('today') - np.timedelta64(1)
today     = np.datetime64('today')
tomorrow  = np.datetime64('today') + np.timedelta64(1)

1
2

Z = np.arange('2016-07', '2016-08', dtype='datetime64[D]')
print(Z)

1
2
3
4
5
6

A = np.ones(3)*1
B = np.ones(3)*2
np.add(A, B, out=B)
np.divide(A, 2, out=A)
np.negative(A, out=A)
np.multiply(A, B, out=A)

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12

Z = np.random.uniform(0,10,10)

# Z%1 获取小数部分
print(Z - Z%1)
# Z 整除 1
print(Z // 1)
# Z 向下取整
print(np.floor(Z))
# 直接转为 int 类型
print(Z.astype(int))
# 截断 truncation
print(np.trunc(Z))

1
2
3
4
5
6
7
8

Z = np.zeros((5,5))
Z += np.arange(5)
print(Z)

# without broadcasting
# np.tile() 函数用于重复数组
Z = np.tile(np.arange(0,5), (5,1))
print(Z)

1
2
3
4
5
6
7

def generate():
    for x in range(10):
        yield x


Z = np.fromiter(generate(), dtype=float, count=-1)
print(Z)

1
2

Z = np.linspace(0,1,11, endpoint=False)[1:]
print(Z)

1
2
3

Z = np.random.random(10)
Z.sort()
print(Z)

1
2
3
4
5
6

# Author: Evgeni Burovski

Z = np.arange(10)
np.add.reduce(Z)
# 它从左到右依次将数组中的元素相加（即求和）
# 可以高效地计算大型数组所有元素的总和

 1
 2
 3
 4
 5
 6
 7
 8
 9
10

A = np.random.randint(0,2,5)
B = np.random.randint(0,2,5)

# Assuming identical shape of the arrays and a tolerance for the comparison of values
equal = np.allclose(A, B)
print(equal)

# Checking both the shape and the element values, no tolerance (values have to be exactly equal)
equal = np.array_equal(A, B)
print(equal)

1
2
3

Z = np.zeros(10)
Z.flags.writeable = False
Z[0] = 1

 1
 2
 3
 4
 5
 6
 7
 8
 9
10

Z = np.random.random((10,2))
# 直角坐标 X、Y
X, Y = Z[:,0], Z[:,1]

# 计算距离
R = np.sqrt(X**2+Y**2)
# 计算夹角
T = np.arctan2(Y, X)
print(R)
print(T)

1
2
3

Z = np.random.random(10)
Z[Z.argmax()] = 0
print(Z)

1
2
3
4

Z = np.zeros((5,5), [('x',float),('y',float)])
Z['x'], Z['y'] = np.meshgrid(np.linspace(0,1,5),
                             np.linspace(0,1,5))
print(Z)

1
2
3
4
5
6
7
8

# Author: Evgeni Burovski

X = np.arange(8)
Y = X + 0.5
# np.subtract.outer(X, Y) 计算 X 和 Y 之间的外差（两个数组中的元素两两之间的差值，构成了一个更大的矩阵）
C = 1.0 / np.subtract.outer(X, Y)
# 计算矩阵 C 的行列式
print(np.linalg.det(C))

1
2
3
4
5
6
7
8

for dtype in [np.int8, np.int32, np.int64]:
   print(np.iinfo(dtype).min)
   print(np.iinfo(dtype).max)

for dtype in [np.float32, np.float64]:
   print(np.finfo(dtype).min)
   print(np.finfo(dtype).max)
   print(np.finfo(dtype).eps)

1
2
3

np.set_printoptions(threshold=float("inf"))
Z = np.zeros((40,40))
print(Z)

1
2
3
4
5

Z = np.arange(100)
v = np.random.uniform(0,100)
# 寻找差值最小的下标
index = (np.abs(Z-v)).argmin()
print(Z[index])

1
2
3
4
5
6

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

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13

Z = np.random.random((10,2))
X, Y = np.atleast_2d(Z[:,0], Z[:,1])
D = np.sqrt( (X-X.T)**2 + (Y-Y.T)**2)
print(D)


# Much faster with scipy
import scipy
# Thanks Gavin Heverly-Coulson (#issue 1)
import scipy.spatial
Z = np.random.random((10,2))
D = scipy.spatial.distance.cdist(Z,Z)
print(D)

1
2
3
4
5
6

# Thanks Vikas (https://stackoverflow.com/a/10622758/5989906)
# & unutbu (https://stackoverflow.com/a/4396247/5989906)
Z = (np.random.rand(10)*100).astype(np.float32)
Y = Z.view(np.int32)
Y[:] = Z
print(Y)

1
2
3

1, 2, 3, 4, 5
6,  ,  , 7, 8
 ,  , 9,10,11

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12

from io import StringIO


# Fake file
s = StringIO('''1, 2, 3, 4, 5

                6,  ,  , 7, 8

                 ,  , 9,10,11
''')
Z = np.genfromtxt(s, delimiter=",", dtype=np.int)
print(Z)

1
2
3
4
5
6

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])

1
2
3
4
5
6
7
8

X, Y = np.meshgrid(np.linspace(-1,1,10), np.linspace(-1,1,10))
D = np.sqrt(X*X+Y*Y)

# 均值为 0，标准差为 1
sigma, mu = 1.0, 0.0

G = np.exp(-( (D-mu)**2 / ( 2.0 * sigma**2 ) ) )
print(G)

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11

# Author: Divakar

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

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11

# Author: Warren Weckesser

X = np.random.rand(5, 10)

# Recent versions of numpy
Y = X - X.mean(axis=1, keepdims=True)

# Older versions of numpy
Y = X - X.mean(axis=1).reshape(-1, 1)

print(Y)

1
2
3
4
5

# Author: Steve Tjoa

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

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18

# Author: Warren Weckesser

# 情况1：检查是否存在全为 0 的列
Z = np.random.randint(0,3, (3,10))
# Z.any(axis=0) 会检查每一列，只要该列有一个非0元素，该列结果就为 True
# 取反后，变成：只有当该列所有元素都是0时，结果才为 True
# 最后的 .any() 用于查看是否至少存在一个这样的列
print((~Z.any(axis=0)).any())

# 情况2：检查是否存在全为 np.nan 的列
Z = np.array([
    [0,1,np.nan],
    [1,2,np.nan],
    [4,5,np.nan]
])
# np.isnan(Z) 将数组中的 NaN 转为 True，数字转为 False
# .all(axis=0) 检查每一列是否全为 True（即全为 NaN）
print(np.isnan(Z).all(axis=0))

1
2
3
4
5
6

Z = np.random.uniform(0,1,10)
z = 0.5
# 计算数组每个元素与目标值 z 的绝对差值，找到差值最小的索引
# Z.flat 返回一个一维迭代器，用于安全地提取元素值
m = Z.flat[np.abs(Z - z).argmin()]
print(m)

1
2
3
4
5
6
7
8
9

A = np.arange(3).reshape(3,1)  # 形状 (3,1)
B = np.arange(3).reshape(1,3)  # 形状 (1,3)

# None 表示在迭代器中分配一个可写的空数组，用于存放结果
it = np.nditer([A, B, None])
for x, y, z in it:
    z[...] = x + y  # 利用广播机制逐个元素求和并赋值
# it.operands[2] 即存放结果的数组，形状为 (3,3)
print(it.operands[2])

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14

class NamedArray(np.ndarray):
    def __new__(cls, array, name="no name"):
        # 将输入对象转换为数组，并以 NamedArray 视图的形式返回
        obj = np.asarray(array).view(cls)
        obj.name = name
        return obj

    def __array_finalize__(self, obj):
        # 当数组被切片或由其他数组创建时，此方法会被调用，用于传递属性
        if obj is None: return
        self.name = getattr(obj, 'name', "no name")

Z = NamedArray(np.arange(10), "range_10")
print(Z.name)

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13

# Author: Brett Olsen

Z = np.ones(10)
I = np.random.randint(0, len(Z), 20)  # 可能包含重复值的索引数组

# 方法1：利用 bincount 统计每个索引出现的次数，一次性累加
Z += np.bincount(I, minlength=len(Z))
print(Z)

# 方法2：使用 np.add.at 进行无缓冲的原位加法，能正确处理重复索引
Z = np.ones(10)  # 重置数组
np.add.at(Z, I, 1)
print(Z)

1
2
3
4
5
6
7
8

# Author: Alan G Isaac

X = [1, 2, 3, 4, 5, 6]
I = [1, 3, 9, 3, 4, 1]
# 在索引 1 处累加 X[0] 和 X[5]（1 + 6 = 7）
# 在索引 3 处累加 X[1] 和 X[3]（2 + 4 = 6）
F = np.bincount(I, X)
print(F)

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14

# Author: Fisher Wang

w, h = 256, 256
I = np.random.randint(0, 4, (h, w, 3)).astype(np.ubyte)
# 常规方法：将像素展平成二维数组 (像素数, 3)，沿行查找唯一值
colors = np.unique(I.reshape(-1, 3), axis=0)
n = len(colors)
print(n)

# 更快的方法：利用点积将 (R, G, B) 压缩为一个整数
# 每个通道占 8 位，B 通道乘 1，G 通道乘 256，R 通道乘 65536
I24 = np.dot(I.astype(np.uint32), [1, 256, 65536])
n = len(np.unique(I24))
print(n)

1
2
3
4
5
6
7
8
9

A = np.random.randint(0, 10, (3, 4, 3, 4))

# 方法1：直接传入需要求和的轴元组 (-2, -1) 代表倒数第二和倒数第一轴
sum_val = A.sum(axis=(-2, -1))
print(sum_val)

# 方法2：如果使用的函数不支持轴元组，可先将最后两维展平再求和
sum_val = A.reshape(A.shape[:-2] + (-1,)).sum(axis=-1)
print(sum_val)

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14

# Author: Jaime Fernández del Río

D = np.random.uniform(0, 1, 100)  # 数据数组
S = np.random.randint(0, 10, 100) # 分组标签（0 到 9）

# 分组求和与分组计数
D_sums = np.bincount(S, weights=D)
D_counts = np.bincount(S)
D_means = D_sums / D_counts
print(D_means)

# 用 Pandas 做同样的事情，代码更直观
import pandas as pd
print(pd.Series(D).groupby(S).mean())

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13

# Author: Mathieu Blondel

A = np.random.uniform(0, 1, (5, 5))
B = np.random.uniform(0, 1, (5, 5))

# 慢速版本：计算所有 25 个乘积再取 5 个对角线元素
np.diag(np.dot(A, B))

# 快速版本：利用广播，仅计算 A 的行与 B 的列对应位置的乘积和
np.sum(A * B.T, axis=1)

# 极速版本：使用爱因斯坦求和约定，精确指定运算 `ij,ji->i`
np.einsum("ij,ji->i", A, B)

 1
 2
 3
 4
 5
 6
 7
 8
 9
10

# Author: Warren Weckesser

Z = np.array([1, 2, 3, 4, 5])
nz = 3  # 插入零的个数

# 新数组长度 = 原长度 + (原长度-1) * 插入数
Z0 = np.zeros(len(Z) + (len(Z) - 1) * nz)
# 每隔 nz+1 个位置放入一个原数组元素
Z0[::nz+1] = Z
print(Z0)  # [1. 0. 0. 0. 2. 0. 0. 0. ...]

1
2
3
4

A = np.ones((5, 5, 3))
B = 2 * np.ones((5, 5))
# B[:,:,None] 将形状从 (5,5) 扩展为 (5,5,1)，然后广播到 (5,5,3)
print(A * B[:, :, None])

1
2
3
4
5
6

# Author: Eelco Hoogendoorn

A = np.arange(25).reshape(5, 5)
# A[[0,1]] 选取第 0 行和第 1 行，赋值时位置对调
A[[0, 1]] = A[[1, 0]]
print(A)

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14

# Author: Nicolas P. Rougier

faces = np.random.randint(0, 100, (10, 3))  # 10 个三角形，每行 3 个顶点索引

# 1. 每个三角形有 3 条边：将 [v0,v1,v2] 扩展为边 [v0,v1], [v1,v2], [v2,v0]
F = np.roll(faces.repeat(2, axis=1), -1, axis=1)
# 2. 展平成 (30, 2) 的边数组
F = F.reshape(len(F) * 3, 2)
# 3. 对每条边的两个顶点进行排序（使无向边标准化）
F = np.sort(F, axis=1)
# 4. 将每条边视为一个结构化标量，利用 np.unique 去重
G = F.view(dtype=[('p0', F.dtype), ('p1', F.dtype)])
G = np.unique(G)
print(G)

1
2
3
4
5
6

# Author: Jaime Fernández del Río

C = np.bincount([1, 1, 2, 3, 4, 4, 6])  # 例如 C[1]=2, C[2]=1, C[3]=1...
# 利用 repeat 函数：索引 i 重复 C[i] 次
A = np.repeat(np.arange(len(C)), C)
print(A)  # 输出应与原数组元素排序一致 [1 1 2 3 4 4 6]

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14

# Author: Jaime Fernández del Río

def moving_average(a, n=3):
    ret = np.cumsum(a, dtype=float)  # 计算前缀和
    ret[n:] = ret[n:] - ret[:-n]     # 利用差分获得窗口内的和
    return ret[n - 1:] / n           # 求平均

Z = np.arange(20)
print(moving_average(Z, n=3))

# NumPy >= 1.20.0 的现代写法
from numpy.lib.stride_tricks import sliding_window_view
Z = np.arange(20)
print(sliding_window_view(Z, window_shape=3).mean(axis=-1))

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13

# Author: Joe Kington / Erik Rigtorp
from numpy.lib import stride_tricks

def rolling(a, window):
    shape = (a.size - window + 1, window)
    strides = (a.strides[0], a.strides[0])  # 行步长与元素步长相同，实现错位
    return stride_tricks.as_strided(a, shape=shape, strides=strides)

Z = rolling(np.arange(10), 3)
print(Z)

# 使用内置函数更安全
print(sliding_window_view(np.arange(10), window_shape=3))

1
2
3
4
5
6
7

# Author: Nathaniel J. Smith

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)     # 原地改变符号

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13

def distance(P0, P1, p):
    T = P1 - P0                         # 线段方向向量
    L = (T**2).sum(axis=1)              # 线段长度的平方
    # 计算投影系数 U
    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                    # 垂足指向 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))
print(distance(P0, P1, p))

1
2
3
4
5
6
7

# Author: Italmassov Kuanysh

P0 = np.random.uniform(-10, 10, (10, 2))
P1 = np.random.uniform(-10, 10, (10, 2))
p = np.random.uniform(-10, 10, (10, 2))
# 列表推导式：对 p 中的每个点计算到所有线段的距离
print(np.array([distance(P0, P1, p_i) for p_i in p]))

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28

# Author: Nicolas Rougier

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[tuple(r)] = Z[tuple(z)]  # 将有效区域复制到填充好的结果数组中
print(Z)
print(R)

1
2
3
4
5

# Author: Stefan van der Walt

Z = np.arange(1, 15, dtype=np.uint32)
R = stride_tricks.as_strided(Z, (11, 4), (4, 4))
print(R)

 1
 2
 3
 4
 5
 6
 7
 8
 9
10

# Author: Stefan van der Walt

Z = np.random.uniform(0, 1, (10, 10))
U, S, V = np.linalg.svd(Z)  # 奇异值分解
rank = np.sum(S > 1e-10)    # 大于容差的奇异值数量
print(rank)

# 更直接的函数
rank = np.linalg.matrix_rank(Z)
print(rank)

1
2
3

Z = np.random.randint(0, 10, 50)
# np.bincount 统计非负整数出现次数，argmax 返回次数最多的值
print(np.bincount(Z).argmax())

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11

# Author: Chris Barker

Z = np.random.randint(0, 5, (10, 10))
n = 3
i = 1 + (Z.shape[0] - 3)  # 8
j = 1 + (Z.shape[1] - 3)  # 8
C = stride_tricks.as_strided(Z, shape=(i, j, n, n), strides=Z.strides + Z.strides)
print(C)

# 现代写法
print(sliding_window_view(Z, window_shape=(3, 3)))

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14

# Author: Eric O. Lebigot

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  # S[3, 2] 也会自动变为 42
print(S)

1
2
3
4
5
6
7
8
9

# Author: Stefan van der Walt

p, n = 10, 20
M = np.ones((p, n, n))  # 批量矩阵
V = np.ones((p, n, 1))  # 批量向量

# 对 M 的第 0 轴（批量）和第 2 轴（列），与 V 的第 0 轴（批量）和第 1 轴（行）求和收缩
S = np.tensordot(M, V, axes=[[0, 2], [0, 1]])
print(S.shape)  # (n, 1)

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13

# Author: Robert Kern

Z = np.ones((16, 16))
k = 4
# 先对行分段求和，再对列分段求和
S = np.add.reduceat(np.add.reduceat(Z, np.arange(0, Z.shape[0], k), axis=0),
                                       np.arange(0, Z.shape[1], k), axis=1)
print(S)

# 使用 sliding_window_view 的备选方案
windows = sliding_window_view(Z, (k, k))
S = windows[::k, ::k].sum(axis=(-2, -1))
print(S)

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18

# Author: Nicolas Rougier

def iterate(Z):
    # 计算每个细胞周围的邻居数（利用 3x3 滑动窗口卷积）
    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:])

    # 应用规则
    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)

1
2
3
4
5
6
7
8
9

Z = np.arange(10000)
np.random.shuffle(Z)
n = 5

# 慢速：对全数组排序
print(Z[np.argsort(Z)[-n:]])

# 快速：仅分区，保证前 n 个为最大的 n 个（但不保证顺序）
print(Z[np.argpartition(-Z, n)[:n]])

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15

# Author: Stefan Van der Walt

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])))