Numpy 高斯拟合——在嘈杂和有趣的数据集上

在本文中，我们将介绍如何使用Python中的Numpy包进行高斯拟合。我们将探讨如何处理一些嘈杂的和有趣的数据集。高斯分布在现实世界中有着广泛的应用，一些物理量或者人群的某些特质往往符合高斯分布，因此了解如何进行高斯拟合这一技术是十分重要的。

阅读更多：Numpy 教程

前置知识

在介绍如何进行高斯拟合之前，我们需要了解一些基本的概念。首先，我们需要了解高斯分布的形式。高斯分布的数学表达式是：

$f(x) = \frac{1}{\sigma \sqrt{2\pi }} e^{-\frac{(x-\mu)^2}{2\sigma^2}}$

其中，μ是高斯分布的均值，σ是它的标准差。我们可以使用Numpy包中的random.randn方法生成符合高斯分布的数据，具体使用方法如下：

import numpy as np

mean = 0
std = 1
sample_size = 1000

data = np.random.randn(sample_size)*std + mean

这样，我们就可以生成一组符合高斯分布的随机数据。

嘈杂的数据集

现实世界中的数据往往是嘈杂的。这意味着数据中可能存在一些异常值和误差。这些噪声数据可能会对我们的高斯拟合产生负面影响。

例如，我们想要拟合以下数据集：

import numpy as np
import matplotlib.pyplot as plt

# 生成带有噪声的数据
x = np.linspace(0, 10, 1000)
y = np.sin(x) + np.random.normal(0, 0.1, len(x))

# 绘制数据图像
plt.scatter(x, y, s=2)
plt.show()

可以看到，数据中有很多噪声，这可能会对高斯拟合造成一定的影响。为了使我们的拟合结果更加准确，我们需要先进行数据平滑处理。可以使用numpy的convolve方法进行平滑处理。例如以下代码：

import numpy as np
import matplotlib.pyplot as plt

# 生成带有噪声的数据
x = np.linspace(0, 10, 1000)
y = np.sin(x) + np.random.normal(0, 0.1, len(x))

# 使用convolve方法做数据平滑处理
box = np.ones(10) / 10.0
y_smooth = np.convolve(y, box, mode='same')

# 绘制平滑后的数据图像
plt.scatter(x, y_smooth, s=2)
plt.show()

这时，我们就可以对数据进行高斯拟合了。具体的代码如下：

from scipy.optimize import curve_fit

# 高斯分布函数
def gaussian_func(x, a, b, c):
    return a*np.exp(-(x-b)**2/(2*c**2))

# 对数据进行高斯拟合
popt, pcov = curve_fit(gaussian_func, x, y_smooth)

# 绘制拟合曲线
plt.plot(x, gaussian_func(x, *popt), 'r-', linewidth=2)
plt.scatter(x, y_smooth, s=2)
plt.show()

我们可以看到，拟合效果还是比较好的。虽然我们对数据进行了平滑处理，但由于数据本身的噪声较大，使得拟合结果并非完全准确。

“有趣”的数据集

除了嘈杂的数据集外，我们还将探讨一些“有趣”的数据集。这些数据集通常是指在某个区域或者某个范围内的数据分布较为密集，并且在某些位置上拥有明显的峰值或谷值。

例如，我们想要拟合以下数据集：

import numpy as np
import matplotlib.pyplot as plt

# 生成含有明显峰值的数据
x = np.linspace(-5, 5, 1000)
y = np.exp(-x**2) + 0.5 * np.exp(-(x-4)**2) + 0.3 * np.exp(-(x+2)**2)

# 绘制数据图像
plt.scatter(x, y, s=2)
plt.show()

我们可以看到，数据集在三个位置上拥有明显的峰值。为了更好地进行高斯拟合，我们需要对每个峰值进行单独拟合。

具体的方法是，我们先使用Numpy的peak_detection方法，找出数据集中的峰值位置。例如以下代码：

from scipy.signal import find_peaks

# 找到数据中的峰值位置
peaks, _ = find_peaks(y)

# 绘制峰值图像
plt.scatter(x, y, s=2)
plt.plot(x[peaks], y[peaks], 'x')
plt.show()

我们可以看到，数据集中的三个峰值位置都被正确地找出来了。然后，我们对每个峰值位置上的数据进行高斯拟合即可。具体的代码如下：

from scipy.optimize import curve_fit

# 高斯分布函数
def gaussian_func(x, a, b, c):
    return a*np.exp(-(x-b)**2/(2*c**2))

# 对每个峰值位置上的数据进行高斯拟合
for i in range(len(peaks)):
    start = peaks[i] - 10
    end = peaks[i] + 10
    popt, pcov = curve_fit(gaussian_func, x[start:end], y[start:end])
    plt.plot(x[start:end], gaussian_func(x[start:end], *popt), 'r-', linewidth=2)

# 绘制拟合曲线
plt.scatter(x, y, s=2)
plt.show()

我们可以看到，每个峰值位置上的高斯拟合结果都很准确。