如何通过 Gnuplot 拟合双峰分布函数？答案

【问题标题】：How to fit a bimodal distribution function via Gnuplot?如何通过 Gnuplot 拟合双峰分布函数？
【发布时间】：2021-10-15 19:21:55
【问题描述】：

我有以下数据，我喜欢用双峰分布拟合

     <0       0.000000e+00         0          0.000000e+00    0.000000e+00         0          0.000000e+00   
     <=0.3       3.000000e-01         1          8.333333e-02    8.333333e+00         1          8.333333e+00   
     <=0.6       6.000000e-01         0          0.000000e+00    0.000000e+00         1          8.333333e+00   
     <=0.9       9.000000e-01         2          1.666667e-01    1.666667e+01         3          2.500000e+01   
     <=1.2       1.200000e+00         6          5.000000e-01    5.000000e+01         9          7.500000e+01   
     <=1.5       1.500000e+00         3          2.500000e-01    2.500000e+01        12          1.000000e+02   
     <=1.8       1.800000e+00         0          0.000000e+00    0.000000e+00        12          1.000000e+02   
     <=2.1       2.100000e+00         0          0.000000e+00    0.000000e+00        12          1.000000e+02   
     <=2.4       2.400000e+00         0          0.000000e+00    0.000000e+00        12          1.000000e+02   
     <=2.7       2.700000e+00         0          0.000000e+00    0.000000e+00        12          1.000000e+02   
       <=3       3.000000e+00         0          0.000000e+00    0.000000e+00        12          1.000000e+02   
        >3       3.000000e+00         0          0.000000e+00    0.000000e+00        12          1.000000e+02

pl 'data.txt' u 2:($5/100) w lp
gauss(x)=a/(sqrt(2*pi)*sigma)*exp(-(x-mean)**2/(2*sigma**2))
fit gauss(x) 'data.txt' u 2:($5/100) via a,sigma,mean

这是的结果，经过拟合，我只得到 1 个峰值

谁能建议如何进行？是否可以通过 Gnuplot 拟合这些数据？

【问题讨论】：

到目前为止您尝试过什么？你的分布函数应该是什么，高斯、洛伦兹还是别的什么？
我希望分布函数是高斯的
@theozh 我已经更新了我的答案

标签： gnuplot distribution data-fitting

【解决方案1】：

只需定义由两个高斯峰组成的拟合函数。并以一些合理的起始值帮助拟合算法。

代码：

### fitting two gauss peaks
reset session

$Data <<EOD
     <  0       0.000000e+00         0          0.000000e+00    0.000000e+00         0          0.000000e+00   
     <=0.3       3.000000e-01         1          8.333333e-02    8.333333e+00         1          8.333333e+00   
     <=0.6       6.000000e-01         0          0.000000e+00    0.000000e+00         1          8.333333e+00   
     <=0.9       9.000000e-01         2          1.666667e-01    1.666667e+01         3          2.500000e+01   
     <=1.2       1.200000e+00         6          5.000000e-01    5.000000e+01         9          7.500000e+01   
     <=1.5       1.500000e+00         3          2.500000e-01    2.500000e+01        12          1.000000e+02   
     <=1.8       1.800000e+00         0          0.000000e+00    0.000000e+00        12          1.000000e+02   
     <=2.1       2.100000e+00         0          0.000000e+00    0.000000e+00        12          1.000000e+02   
     <=2.4       2.400000e+00         0          0.000000e+00    0.000000e+00        12          1.000000e+02   
     <=2.7       2.700000e+00         0          0.000000e+00    0.000000e+00        12          1.000000e+02   
       <=3       3.000000e+00         0          0.000000e+00    0.000000e+00        12          1.000000e+02   
        >3       3.000000e+00         0          0.000000e+00    0.000000e+00        12          1.000000e+02   
EOD

Gauss(x,a,mean,sigma) = a/(sqrt(2*pi)*sigma)*exp(-(x-mean)**2/(2*sigma**2))
g1(x) = Gauss(x,a1,mean1,sigma1)
g2(x) = Gauss(x,a2,mean2,sigma2)
f(x) = g1(x) + g2(x)

# set some initial values to help the fitting algorithm
mean1 = 0.3
mean2 = 1.2

set fit results nolog
fit f(x) $Data u 2:($5/100) via a1,mean1,sigma1,a2,mean2,sigma2

set samples 200

myResults = sprintf("a1 = %g\nmean1 = %g\nsigma1 = %g\na2 = %g\nmean2 = %g\nsigma2 = %g\n",a1,mean1,sigma1,a2,mean2,sigma2)
set label 1 at graph 0.05, 0.97 myResults

plot $Data u 2:($5/100)  w lp pt 7 dt 3 ti "Data", \
     f(x) w l lc "red"
### end of code

结果：

【讨论】：

非常感谢，我意识到我使用了错误的功能。