使用python拟合“多模态”对数正态分布到数据
使用python拟合“多模态”对数正态分布到数据
提问于 2019-11-25 23:02:50
我在实验室里用仪器测量了以下数据。由于仪器根据直径将不同大小的颗粒收集在垃圾箱中，因此测量结果基本上是“入库”的：
import numpy as np
import matplotlib.pylab as plt
from lmfit import models
y = np.array([196, 486, 968, 2262, 3321, 4203, 15072, 46789, 95201, 303494, 421484, 327507, 138931, 27973])
bins = np.array([0.0150, 0.0306, 0.0548, 0.0944, 0.1540, 0.2560, 0.3830, 0.6050, 0.9510, 1.6400, 2.4800, 3.6700, 5.3800, 9.9100, 15])
bin_width=np.diff(bins)
x_plot = np.add(bins[:-1],np.divide(bin_width,2))
x=x_plot
y=y
在这里绘制的是数据的外观。在x尺度的单位中，有一个模式在0.1左右，另一个模式在2左右。
在这个研究领域内，通常将“多峰”对数正态分布拟合到这样的数据:鉴于此，我使用LMFIT拟合了大约2的模型：
model = models.LognormalModel()
params = model.make_params(center=1.5, sigma=0.6, amplitude=2214337)
result = model.fit(y, params, x=x)
print(result.fit_report())
plt.plot(x, y, label='data')
plt.plot(x, result.best_fit, label='fit')
plt.xscale("log")
plt.yscale("log")
plt.legend()
plt.show()
正如预期的那样，这将导致第二个模式在2左右很好地拟合。我的问题是，我如何也拟合0.1左右的第二个模式(本质上两个模式的总和应该符合数据)？我意识到也可以争论三种模式会更好，但我假设一旦我了解了如何使用两种模式，增加第三种模式应该是微不足道的。
回答 2
Stack Overflow用户
回答已采纳
发布于 2019-11-26 03:05:59
可以将lmfit.Models添加到一起，如：
model = (models.LognormalModel(prefix='p1_') +
models.LognormalModel(prefix='p2_') +
models.LognormalModel(prefix='p3_') )
params = model.make_params(p1_center=0.3, p1_sigma=0.2, p1_amplitude=1e4,
p2_center=1.0, p2_sigma=0.4, p2_amplitude=1e6,
p3_center=1.5, p3_sigma=0.6, p3_amplitude=2e7)在复合模型中，model的每个组件都有自己的“前缀”(任何字符串)，这些前缀位于参数名称的前面。你可以在拟合后得到一个模型组件的字典：
components = result.eval_components()
# {'p1_': Array, 'p2_': Array, 'p3_': Array}
for compname, comp in components.keys():
plt.plot(x, comp, label=compname)要拟合在半对数或对数对数图上表示的数据，可以考虑将模型拟合到log(y)。否则，在非常低的y值下，拟合将不会对不拟合非常敏感。
请注意，lmfit模型和参数支持边界。您可能想要或发现您需要放置边界，例如
params['p1_amplitude'].min = 0
params['p1_sigma'].min = 0
params['p1_center'].max = 0.5
params['p3_center'].min = 1.0Stack Overflow用户
发布于 2019-11-25 23:30:04
这是你想要拟合的对数正态混合分布。你可以简单地取你的数据的对数，然后拟合一个高斯混合：
import numpy as np
from sklearn.mixture import GaussianMixture
# Make data from two log-normal distributions
# NOTE: 2d array of shape (n_samples, n_features)
n = 10000
x = np.zeros((n,1))
x[:n//2] = np.random.lognormal(0,1, size=(n//2,1))
x[n//2:] = np.random.lognormal(2,0.5, size=(n//2,1))
# Log transform the data
x_transformed = np.log(x)
# Make gaussian mixture model.
# n_init makes multiple initial guesses and
# depending on data, 1 might be good enough
# Decrease tolerance for speedup or increase for better precision
m = GaussianMixture(n_components=2, n_init=10, tol=1e-6)
# Fit the model
m.fit(x_transformed)
# Get the fitted parameters
# NOTE: covariances are stdev**2
print(m.weights_) # [0.50183897 0.49816103]
print(m.means_) # [1.99866785, -0.00528186]
print(m.covariances_) # [0.25227372,0.99692494]页面原文内容由Stack Overflow提供。腾讯云小微IT领域专用引擎提供翻译支持
原文链接：
https://stackoverflow.com/questions/59034668
复制相关文章
相似问题
腾讯云开发者
