-1
$\begingroup$

I want to fit this data. enter image description here

I have the following model functions. Classic gaussian:

    def gauss_model(x, mu, sigma):
    return np.exp(-0.5*((x-mu)/sigma)**2)

And a gauss-hermite parametrization, taking account for skewness (parameter h3) and kurtosis (parameter h4) of the distribution:

    def losvd_param(v, v_rot, v_disp, h3, h4):
    y = np.asarray((np.asarray(v)-v_rot)/(v_disp))
    return (np.exp(-0.5 * y**2) * (1 + h3*((2*np.sqrt(2)*y**3-3*np.sqrt(2)*y)/np.sqrt(6)) + h4*((4*y**4-12*y**2+3)/np.sqrt(24))))

Fitting and plotting:

    p0_1=[1300,250]
    boundslow_1=[1000, 100]
    boundsup_1 = [1400, 300]
    p0_2 = [1300, 250, 0, 0]
    boundslow_2 = [1000, 100, -0.1, -0.1]
    boundsup_2 = [1400, 300, 0.1, 0.1]
    gh_moments = curve_fit(gauss_model, vel_corr_peak, broadening_func/max(broadening_func),p0=p0_1, bounds=((boundslow_1),(boundsup_1)))[0]
    gh_moments_2 = curve_fit(losvd_param, vel_corr_peak, broadening_func/max(broadening_func), p0=p0_2,bounds=((boundslow_2),(boundsup_2)))[0]
    plt.plot(vel_corr_peak, broadening_func/max(broadening_func),linewidth='2', label='data', color='black')
    plt.plot(vel_corr_peak, gauss_model(vel_corr_peak, *gh_moments)/max(gauss_model(vel_corr_peak, *gh_moments)),'--',color='red', label='gauss fit')
    plt.plot(vel_corr_peak, losvd_param(vel_corr_peak, *gh_moments_2)/max(losvd_param(vel_corr_peak, *gh_moments_2)),':', color='green', label='gauss hermite fit')

enter image description here

One can see that the gauss-hermite-fit is a bit closer to the actual data but dont match perfectly.

If somebody wants to try the data: x-values:

[449.99918644287044, 479.9991322057285, 509.9990779685865, 539.9990237314445, 569.9989694943025, 599.9989152571605, 629.9988610200186, 659.9988067828766, 689.9987525457346, 719.9986983085927, 749.9986440714507, 779.9985898343087, 809.9985355971668, 839.9984813600248, 869.9984271228828, 899.9983728857409, 929.9983186485989, 959.998264411457, 989.998210174315, 1019.998155937173, 1049.998101700031, 1079.998047462889, 1109.997993225747, 1139.997938988605, 1169.997884751463, 1199.997830514321, 1229.9977762771791, 1259.9977220400372, 1289.9976678028952, 1319.9976135657532, 1349.9975593286113, 1379.9975050914693, 1409.9974508543273, 1439.9973966171854, 1469.9973423800434, 1499.9972881429014, 1529.9972339057595, 1559.9971796686175, 1589.9971254314755, 1619.9970711943336, 1649.9970169571916, 1679.9969627200496, 1709.9969084829077, 1739.9968542457657, 1769.9968000086237, 1799.9967457714818, 1829.9966915343398, 1859.9966372971978, 1889.9965830600559, 1919.996528822914, 1949.996474585772, 1979.99642034863, 2009.996366111488, 2039.996311874346, 2069.996257637204, 2099.996203400062, 2129.99614916292, 2159.996094925778, 2189.996040688636, 2219.995986451494, 2249.9959322143523]

y-values:

[ 0.48158957 0.51413024 0.57190115 0.65708852 0.77221168 0.91982261 1.10264097 1.323603 1.58568143 1.89166066 2.24396947 2.644586 3.09482201 3.59482773 4.14336694 4.73762031 5.3729276 6.04262429 6.73776387 7.44733358 8.15834164 8.85608227 9.52454256 10.14681488 10.7061444 11.18626665 11.5724327 11.85204538 12.0153712 12.05641108 11.97274977 11.76629065 11.44282113 11.01225276 10.48774599 9.88509761 9.22202859 8.51728797 7.7897047 7.05722644 6.33620579 5.64081325 4.98267552 4.37045974 3.8098805 3.30392459 2.85318343 2.45615909 2.10979453 1.80999646 1.55202573 1.33111531 1.14267783 0.98268788 0.84772486 0.73512314 0.64325468 0.57105697 0.51828786 0.48509208 0.47235733]

Improve this question
$\endgroup$
3
  • 2
    $\begingroup$ Could you include your full code and enough of your data so that someone can try to reproduce your results? $\endgroup$
    – Tyberius
    May 6, 2022 at 0:36
  • $\begingroup$ Please clarify your specific problem or provide additional details to highlight exactly what you need. As it's currently written, it's hard to tell exactly what you're asking. $\endgroup$
    – Community Bot
    May 9, 2022 at 13:47
  • $\begingroup$ Where is this data coming from ? Why a Gaussian or Hermite model ? Just because it is bell-shaped or is there a true rationale ? $\endgroup$
    – Yves Daoust
    May 12, 2022 at 8:20

1 Answer 1

Reset to default
4
$\begingroup$

You're not providing an initial guess for the parameters, and so optimize.curve_fit is defaulting to [1.,1.,1.,1.]. The solver assumes that the true answer is roughly the same order as the initial guess, when in reality they differ by a factor of 1000 (for the first parameter at least). Any step it takes in any direction in parameter space doesn't meaningfully change the error between fit and data; it is effectively lost in the flat plane of the objective landscape, and the actual minimum is beyond the horizon. It tries to move in each direction, makes no progress, and throws its hands up and quits.

The "correct" initial guesses are the mean, standard deviation, skewness, and kurtosis of the data, which you can easily calculate numerically, either by numerical integration (if you have data that looks like your function) or by the expectations (if you have realizations from the distributions). Providing these as the initial guess for the solver gets it to converge just fine for me, or even things that are roughly close to the true parameters (e.g. [1000, 100, 1, 1]).

Improve this answer
$\endgroup$
1
  • $\begingroup$ The problem I keep getting is that even with the initial parameter guesses relatively close the fit data still differs significantly from my real data. Also small changes in my initial guesses give me huge differences in my final fit. f.e. initial guess [1250, 200, 0, 0] gives me the fit [ 929.68956083 1277.90277059 -10.92366431 9.70117469] which is wrong. contrary the initial guess [1000, 100, 1, 1] results in the fit [-387.78490945 720.37366804 7.72216571 5.54750336] which is also totally off. What could be a reason for this ? Is there another workarround $\endgroup$
    – trynerror
    May 10, 2022 at 9:38

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.