请教一下,我有个二维的概率矩阵,理论上想采用Gibbs采样算法,运用马尔可夫链蒙特卡洛模型,但代码报错,是要修改数据类型吗?如何修改代码?
运行报错如下:
sample_freq = np.bincount(samples) / num_samples
^^^^^^^^^^^^^^^^^
ValueError: object too deep for desired array
代码如下:
import numpy as np
# 假设你有一个二维的概率矩阵
prob_matrix = np.array([[0.1, 0.2, 0.7],
[0.3, 0.4, 0.3],
[0.6, 0.3, 0.1]])
# 确保概率矩阵的每一行之和为1
assert np.allclose(prob_matrix.sum(axis=1), 1)
# Gibbs采样函数
def gibbs_sampling(prob_matrix, num_samples):
# 初始化样本
samples = np.empty((num_samples, prob_matrix.shape[0]), dtype=int)
# 随机选择一个初始状态
current_state = np.random.choice(prob_matrix.shape[0], p=prob_matrix[0])
samples[0] = current_state
return samples
# 设置采样次数
num_samples = 10000
# 进行Gibbs采样
samples = gibbs_sampling(prob_matrix, num_samples)
# 打印样本的频率分布,以验证是否接近原始概率矩阵
sample_freq = np.bincount(samples) / num_samples
print("Sample frequencies:", sample_freq)