來自 MCMC 樣本的模式可靠性

John Kruschke 在他的《做貝葉斯數據分析》一書中指出，使用 R 中的 JAGS

…從 MCMC 樣本中估計模式可能相當不穩定，因為該估計基於平滑算法，該算法可能對 MCMC 樣本中的隨機顛簸和波紋敏感。（進行貝葉斯數據分析，第 205 頁，第 8.2.5.1 節）

雖然我掌握了 Metropolis 算法和 Gibbs 採樣等精確形式，但我不熟悉也提到的平滑算法，以及為什麼這意味著從 MCMC 樣本中對模式的估計是不穩定的。有沒有人能夠直觀地了解平滑算法正在做什麼以及為什麼它會使模式的估計不穩定？

我手頭沒有這本書，所以我不確定 Kruschke 使用什麼平滑方法，但出於直覺，請考慮這個由標準法線中的 100 個樣本組成的圖，以及使用 0.1 到 1.0 的各種帶寬的高斯核密度估計。（簡而言之，高斯 KDE 是一種平滑直方圖：它們通過為每個數據點添加高斯來估計密度，平均值為觀察值。）

您可以看到，即使平滑創建了單峰分佈，其眾數通常也低於已知值 0。

此外，這裡是使用相同樣本的用於估計密度的內核帶寬的估計模式（y 軸）圖。希望這可以直觀地了解估計值如何隨平滑參數的變化而變化。

#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
Created on Wed Feb  1 09:35:51 2017

@author: seaneaster
"""

import numpy as np
from matplotlib import pylab as plt
from sklearn.neighbors import KernelDensity

REAL_MODE = 0
np.random.seed(123)

def estimate_mode(X, bandwidth = 0.75):
   kde = KernelDensity(kernel = 'gaussian', bandwidth = bandwidth).fit(X)
   u = np.linspace(-3,3,num=1000)[:, np.newaxis]
   log_density = kde.score_samples(u)
   return u[np.argmax(log_density)]

X = np.random.normal(REAL_MODE, size = 100)[:, np.newaxis] # keeping to standard normal

bandwidths = np.linspace(0.1, 1., num = 8)

plt.figure(0)
plt.hist(X, bins = 100, normed = True, alpha = 0.25)

for bandwidth in bandwidths:
   kde = KernelDensity(kernel = 'gaussian', bandwidth = bandwidth).fit(X)
   u = np.linspace(-3,3,num=1000)[:, np.newaxis]
   log_density = kde.score_samples(u)
   plt.plot(u, np.exp(log_density))

bandwidths = np.linspace(0.1, 3., num = 100)
modes = [estimate_mode(X, bandwidth) for bandwidth in bandwidths]
plt.figure(1)
plt.plot(bandwidths, np.array(modes))

引用自：https://stats.stackexchange.com/questions/259319