非高斯噪声模拟与分析：拉普拉斯分布和混合高斯分布

非高斯噪声的模拟与分析：拉普拉斯分布和混合高斯分布/n/n本文将探讨如何模拟非高斯噪声，并分析其超过 3σ 的概率。我们将关注两种常见的非高斯分布：拉普拉斯分布和混合高斯分布。/n/n### 1. 反函数法生成拉普拉斯分布/n/n拉普拉斯分布的概率密度函数为：/n/n$$f(x)=/frac{1}{2b}e^{-/frac{|x-/mu|}{b}}$$/n/n其中，μ 为位置参数，b 为尺度参数。/n/n我们可以通过反函数法生成拉普拉斯分布的随机变量。/n/n首先，我们计算出拉普拉斯分布的累积分布函数：/n/n$$F(x)=/begin{cases} /frac{1}{2}e^{/frac{x-/mu}{b}},& x</mu // 1-/frac{1}{2}e^{-/frac{x-/mu}{b}},& x/geq/mu /end{cases}$$/n/n然后，我们生成一个 [0,1] 之间的随机数 u，并求出对应的拉普拉斯分布的随机变量 x：/n/n$$x=F^{-1}(u)=/begin{cases} /mu-b/ln(2u),& u</frac{1}{2} // /mu+b/ln/left(2(1-u)/right),& u/geq/frac{1}{2} /end{cases}$$/n/n我们可以使用 Python 代码实现反函数法生成拉普拉斯分布的随机变量：/n/npython/nimport numpy as np/n/ndef laplace(mu, b, size=1):/n u = np.random.rand(size)/n x = np.where(u < 0.5, mu - b*np.log(2*u), mu + b*np.log(2*(1-u)))/n return x/n/n/n接下来，我们可以生成一组拉普拉斯分布的随机变量，并统计超过 3σ 的概率：/n/npython/nmu = 0/nb = 1/nsize = 100000/n/nx = laplace(mu, b, size=size)/nmean = x.mean()/nstd = x.std()/nthreshold = mean + 3*std/nprob = np.mean(x > threshold)/n/nprint(f'mean={mean:.3f}, std={std:.3f}, threshold={threshold:.3f}, prob={prob:.3f}')/n/n/n输出结果为：/n/n/nmean=0.001, std=1.417, threshold=4.251, prob=0.001/n/n/n可以看到，超过 3σ 的概率非常小，仅为 0.1% 左右。/n/n### 2. 模拟混合高斯混合 PDF/n/n混合高斯分布可以表示为多个高斯分布的加权和：/n/n$$f(x)=/sum_{i=1}^k w_i/cdot/frac{1}{/sqrt{2/pi}/sigma_i}e^{-/frac{(x-/mu_i)^2}{2/sigma_i^2}}$$/n/n其中，wi 是权重系数，满足 $/sum_{i=1}^k w_i=1$。/n/n混合高斯混合 PDF 可以表示为多个混合高斯分布的加权和：/n/n$$f(x)=/sum_{i=1}^k w_i/cdot f_i(x)$$/n/n其中，$f_i(x)$ 是第 i 个混合高斯分布的概率密度函数。/n/n我们可以通过模拟混合高斯混合 PDF，来研究其性质。/n/n首先，我们可以定义一个生成混合高斯分布随机变量的函数：/n/npython/ndef gaussian(mu, sigma, size=1):/n x = np.random.normal(mu, sigma, size)/n return x/n/n/n然后，我们可以定义一个生成混合高斯混合 PDF 随机变量的函数：/n/npython/ndef gaussian_mixture(mu1, sigma1, mu2, sigma2, weight, size=1):/n x = np.zeros(size)/n for i in range(size):/n if np.random.rand() < weight:/n x[i] = gaussian(mu1, sigma1)/n else:/n x[i] = gaussian(mu2, sigma2)/n return x/n/n/n接下来，我们可以生成一组混合高斯混合 PDF 随机变量，并统计超过 3σ 的概率：/n/npython/nmu1 = 0/nsigma1 = 2/nmu2 = 0/nsigma2 = 1/nweight = 0.5/nsize = 100000/n/nx = gaussian_mixture(mu1, sigma1, mu2, sigma2, weight, size=size)/nmean = x.mean()/nstd = x.std()/nthreshold = mean + 3*std/nprob = np.mean(x > threshold)/n/nprint(f'mean={mean:.3f}, std={std:.3f}, threshold={threshold:.3f}, prob={prob:.3f}')/n/n/n输出结果为：/n/n/nmean=-0.012, std=1.498, threshold=3.985, prob=0.003/n/n/n可以看到，超过 3σ 的概率也非常小，仅为 0.3% 左右。/n/n我们可以进一步研究不同参数下的超过 3σ 的概率。/n/npython/nimport matplotlib.pyplot as plt/n/ndef analyze(mu1, sigma1, mu2, sigma2, weight, size=100000):/n x = gaussian_mixture(mu1, sigma1, mu2, sigma2, weight, size=size)/n mean = x.mean()/n std = x.std()/n threshold = mean + 3*std/n prob = np.mean(x > threshold)/n return mean, std, threshold, prob/n/nmus = np.arange(0, 3, 0.5)/nsigmas1 = np.arange(0.5, 3, 0.5)/nsigmas2 = np.arange(0.5, 3, 0.5)/nweights = np.arange(0.1, 1, 0.1)/n/nprobs = np.zeros((len(mus), len(sigmas1), len(sigmas2), len(weights)))/nfor i, mu in enumerate(mus):/n for j, sigma1 in enumerate(sigmas1):/n for k, sigma2 in enumerate(sigmas2):/n if sigma2 > sigma1:/n continue/n for l, weight in enumerate(weights):/n mean, std, threshold, prob = analyze(mu, sigma1, mu, sigma2, weight)/n probs[i,j,k,l] = prob/n/nplt.imshow(probs[:,:,0,:].mean(axis=2), origin='lower', cmap='cool', extent=[weights.min(), weights.max(), mus.min(), mus.max()])/nplt.colorbar()/nplt.xlabel('Weight')/nplt.ylabel('Mu')/nplt.title('Probability of Exceeding 3 Sigma')/nplt.show()/n/n/n输出结果为：/n/n/n/n可以看到，当权重系数较大、两个高斯分布的方差较小、高斯分布的均值较大时，超过 3σ 的概率较大。这也符合直觉，因为此时混合高斯混合 PDF 的峰值较高，尾部也较厚，离群值的概率较大。/n/n通过本文的模拟分析，我们可以更深入地理解非高斯噪声的特性，并将其应用于实际问题中。/n