蒙特卡洛积分不同样本量误差计算循环报错及优化问询

错误原因分析 #

你遇到的TypeError: only size-1 arrays can be converted to Python scalars和后续的索引错误，根源有三点：

1. 数组转整数失败：int(sample_size)试图将整个长度为10的numpy数组转为整数，而int()仅支持单个数值（size-1数组）的转换，直接触发报错。
2. 列表索引赋值错误：naive和errornaive是空列表，无法通过naive[i]直接赋值（列表索引赋值要求对应位置已存在元素），应该用append()添加新元素。
3. 样本量计算逻辑错误：m = 2**i中，i是样本量本身（如2、4），再次计算2^{i会导致样本量爆炸式增长（如2}2=4、2^4=16），完全偏离你需要的2^1到2^10的需求。

修正后的基础版本代码 #

以下是修复所有错误、满足需求的代码：

import numpy as np

truetheta = 0.2
N = 10
# 直接生成样本量列表：2^1 到 2^10
sample_sizes = [2**(i+1) for i in range(N)]

naive_estimates = []
errors = []

for m in sample_sizes:
    x = np.random.uniform(0, 1, m)
    y = x ** 4
    # 用np.mean替代sum(y)/m，代码更简洁
    estimate = np.mean(y)
    naive_estimates.append(estimate)
    errors.append(abs(estimate - truetheta))

# 格式化打印结果
for size, est, err in zip(sample_sizes, naive_estimates, errors):
    print(f"样本量M={size}: 估计值={est:.6f}, 误差={err:.6f}")

更优实现方案：向量化计算+误差可视化 #

如果需要更高效率或直观观察误差收敛趋势，可采用numpy向量化操作减少Python循环开销，并添加可视化：

import numpy as np
import matplotlib.pyplot as plt

truetheta = 0.2
N = 10
exponents = np.arange(1, N+1)
sample_sizes = 2 ** exponents  # 生成2^1到2^10的样本量数组

# 一次性生成所有随机数，减少多次调用随机函数的开销
max_sample_size = sample_sizes[-1]
x = np.random.uniform(0, 1, (N, max_sample_size))

# 按每个样本量截断数组并计算均值
naive_estimates = np.array([np.mean(x[i, :m]) for i, m in enumerate(sample_sizes)])
errors = np.abs(naive_estimates - truetheta)

# 格式化输出结果
print("样本量 | 估计值    | 误差")
print("-------------------------")
for m, est, err in zip(sample_sizes, naive_estimates, errors):
    print(f"{m:5d} | {est:.6f} | {err:.6f}")

# 绘制误差收敛曲线
plt.figure(figsize=(8, 4))
plt.plot(sample_sizes, errors, 'bo-')
plt.xscale('log', base=2)
plt.xlabel('样本量 M (2^i)')
plt.ylabel('绝对误差')
plt.title('蒙特卡洛积分误差随样本量的收敛情况')
plt.grid(True)
plt.show()

这个方案的优势：

• 一次性生成所有随机数，提升计算效率
• 用列表推导式结合numpy函数，代码更紧凑
• 可视化误差变化，直观观察蒙特卡洛方法的收敛特性

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