给定一个现有分布，如何绘制标准差为X的大小为N的样本？

Question

我有一个值的现有分布，我想绘制5号样本，但那5个样本的标准差必须在一定公差内。例如，我需要5个样本的std为10（即使总体分布为std =〜32）。

下面的示例代码有些起作用，但是对于大型数据集而言它的运行速度很慢。它对分布进行随机采样，直到找到接近目标std的对象，然后删除这些元素，以使它们不再被绘制。

是否有更聪明的方法来正确，更快地执行此操作？对于某些target_std（6以上），它可以正常工作，但是在6以下是不准确的。

import numpy as np
import matplotlib.pyplot as plt
np.random.seed(23)

# Create a distribution
d1 = np.random.normal(95, 5, 200)
d2 = np.random.normal(125, 5, 200)
d3 = np.random.normal(115, 10, 200)
d4 = np.random.normal(70, 10, 100)
d5 = np.random.normal(160, 5, 200)
d6 = np.random.normal(170, 20, 100)
dist = np.concatenate((d1, d2, d3, d4, d5, d6))
print(f"Full distribution: len={len(dist)}, mean={np.mean(dist)}, std={np.std(dist)}")
plt.hist(dist, bins=100)
plt.title("Full Distribution")
plt.show();


batch_size = 5
num_batches = math.ceil(len(dist)/batch_size)
target_std = 10
tolerance = 1
# how many samples to search
num_samples = 100
result = []

# Find samples of batch_size that are closest to target_std
for i in range(num_batches):
    samples = []
    idxs = np.arange(len(dist))
    for j in range(num_samples):
        indices = np.random.choice(idxs, size=batch_size, replace=False)
        sample = dist[indices]
        std = sample.std()
        err = abs(std - target_std)
        samples.append((sample, indices, std, err, np.mean(sample), max(sample), min(sample)))
        if err <= tolerance:
            # close enough, stop sampling
            break
    # sort by smallest err first, then take the first/best result
    samples = sorted(samples, key=lambda x: x[3])
    best = samples[0] 
    if i % 100 == 0:
        pass
        print(f"{i}, std={best[2]}, err={best[3]}, nsamples={num_samples}")
    result.append(best)
    # remove the data from our source
    dist = np.delete(dist, best[1])

df_samples = pd.DataFrame(result, columns=["sample", "indices", "std", "err", "mean", "max", "min"])

df_samples["err"].plot(title="Errors (target_std - batch_std)")
batch_std = df_samples["std"].mean()
batch_err = df_samples["err"].mean()
print(f"RESULT: Target std: {target_std}, Mean batch std: {batch_std}, Mean batch err: {batch_err}")

Answer 1

由于您的问题不仅限于特定分布，因此我使用通常是随机的分布，但这适用于任何分布。但是，运行时间将取决于总体数量。

population = np.random.randn(1000)*32
std = 10.
tol = 1.
n_samples = 5
samples = list(np.random.choice(population, n_samples))
while True:
    center = np.mean(samples)
    dis = [abs(i-center) for i in samples]
    if np.std(samples)>(std+tol):
        samples.pop(dis.index(max(dis)))
    elif np.std(samples)<(std-tol):
        samples.pop(dis.index(min(dis)))
    else:
        break
    samples.append(np.random.choice(population, 1)[0])

这里是代码的工作方式。首先，绘制n_samples，可能是std不在您想要的范围内，因此我们计算平均值和每个样本与平均值的绝对距离。然后，如果std大于所需值加公差，我们将踢最远的样本并绘制一个新样本，反之亦然。

[请注意，如果花费太多时间来计算数据，则在排除异常值之后，您可以计算总体中应绘制的下一个元素的范围，而不是随机取一个。希望这对您有用。

免责声明：这不再是随机抽签，您应该意识到抽签是有偏见的，并不代表总体。