如何基于耦合的粒子位置数组X和速度数组V计算各分箱的速度方差

实现方法 #

你可以通过np.digitize先建立每个位置样本到分箱的映射关系，再基于映射统计每个分箱的速度方差，以下是两种可直接运行的实现方案：

方案1：纯Numpy实现（无额外依赖） #

不需要安装其他第三方库，完全兼容你现有代码逻辑：

import numpy as np
import matplotlib.pyplot as plt

X = np.random.randn(1000)
V = 3*np.random.randn(1000) + 40

bins = np.arange(int(X.min()) - 0.5, int(X.max())+1.5, 1)
# 1. 获取每个X元素对应的分箱索引，即你需要的元素-分箱映射关系
bin_indices = np.digitize(X, bins)
# 过滤出在分箱范围内的有效样本
valid_mask = (bin_indices > 0) & (bin_indices < len(bins))
valid_bin_idx = bin_indices[valid_mask]
valid_V = V[valid_mask]

# 2. 计算每个分箱的速度方差，空分箱默认填充nan
v_var = np.full(len(bins)-1, np.nan)
for i in range(1, len(bins)):
    bin_v = valid_V[valid_bin_idx == i]
    # 至少2个样本才计算方差，可根据需求调整规则
    if len(bin_v) >= 2:
        v_var[i-1] = np.var(bin_v)

# 3. 原有直方图绘制逻辑不变
plt.hist(X, bins=bins, facecolor = '#2ab0ff', edgecolor='#169acf', linewidth=0.7)
plt.show()

# 可直接输出每个分箱中心对应的方差结果
bin_centers = (bins[:-1] + bins[1:]) / 2
for center, var in zip(bin_centers, v_var):
    print(f"分箱中心{center:.0f}，速度方差：{var if not np.isnan(var) else '无有效样本'}")

方案2：Scipy简洁实现 #

如果允许引入scipy依赖，可直接调用封装好的分箱统计函数，一行完成计算：

import numpy as np
import matplotlib.pyplot as plt
from scipy.stats import binned_statistic

X = np.random.randn(1000)
V = 3*np.random.randn(1000) + 40
bins = np.arange(int(X.min()) - 0.5, int(X.max())+1.5, 1)

# 直接按X分箱，统计对应V的方差，返回结果顺序和直方图分箱完全对齐
v_var, bins, bin_mapping = binned_statistic(X, V, statistic='var', bins=bins)

plt.hist(X, bins=bins, facecolor = '#2ab0ff', edgecolor='#169acf', linewidth=0.7)
plt.show()

内容的提问来源于stack exchange，提问作者Oo7