Fitting a histogram with Vavilov function

chetna · January 10, 2011, 9:17am

Hi ROOTers,
I want to fit a histogram with Vavilov function.But always It is seen that the peak of a Vavilov function is at zero,no matter whatever initial value I give to the parameters. But my histo has a peak at 25, ranged from 5 to 35,(I have plotted in the range 0-40,40 bins).Attached here with the data file for the histo.
Please help.
data.txt (72 KB)

moneta · January 11, 2011, 10:15am

Hi,

you need to add some extra parameters to the Vavilov ( a location and a scale parameters).
Attached is an example doing the fit using the ROOT::Math::Vavilov class from MathMore (available in ROOT 5.28)

Best Regards

Lorenzo
fit_Vavilov.C (925 Bytes)

chetna · January 12, 2011, 12:19pm

Hi,
Thanks a lot.
Its working fine now.

chetna · March 16, 2011, 6:54am

Hi ROOTers,
The Vavilov function is fitting my distributions well.But I could not directly extract the mean and variance with errors.
Please Help.

  Chetna.

moneta · March 16, 2011, 8:51am

Hi,

See the attached updated macro which retrieve also the parameters from the fit,

Best Regards

Lorenzo
fit_Vavilov.C (1.21 KB)

chetna · March 16, 2011, 10:46am

Hi,
the parameters p2 and p3 are not the actual mean and sigma of the distribution.I have attached a plot from which this thing becomes clear.The mean and sigma (i.e. p2 and p3) shown in the plot are not true.

Now see,
If p0 is kappa and p1 is beta2,then
mean=gamma-1-ln(p0)-p1.

Now since I am using the locators p2 and p3,so the actual mean should be,
MEAN = mean*p3+p2.

Similarly,
actual sigma should be,
SIGMA = sigmap3, where, sigma=sqrt[(2-p1)/(2p0)]

I have calculated this for several datasets and the result is consistent.

Now my question is,how to find the error in MEAN and SIGMA.Can it be extracted from ROOT directly?
plot.pdf (14.6 KB)

moneta · March 18, 2011, 7:58am

Hi,

Yes, you are right. p2 and p3 are just some location and scale parameter but do not correspond to the mean and standard deviation of the distribution.
Since you have the relation of the real mean and sigma as function of p0,p1,p2,p3 you just to need to apply the error propagation using the covariance matrix returned from the fit.
You get the matrix by doing :

TFitResultPtr r =  h1->Fit(f1,"S");
TMatrixDSym cov = r->GetCovarianceMatrix();

Otherwise you can also try to re-parametrize your fit function in terms of the real mean and sigma and then fit for those. You should get approximately the same results.

Best Regards

Lorenzo

chetna · March 18, 2011, 8:34am

Hi Lorenzo,
Thanks.
Using covariance matrix I am externally calculating the errors.But the thing is, the error propagation is very complicated. So I guess reparametrization is a better idea.

However,I have not been sucessful in reparametrizing the fit function.I was trying it for a few days.Can you give me some hint on it?

cplager · March 18, 2011, 4:21pm

[quote=“moneta”]Hi,

Yes, you are right. p2 and p3 are just some location and scale parameter but do not correspond to the mean and standard deviation of the distribution.
Since you have the relation of the real mean and sigma as function of p0,p1,p2,p3 you just to need to apply the error propagation using the covariance matrix returned from the fit.
You get the matrix by doing :

TFitResultPtr r = h1->Fit(f1,"S"); TMatrixDSym cov = r->GetCovarianceMatrix();[/quote]

Assuming you have a function (i.e., relationship) that converts the parameters above to the mean and sigma you want, what you can do is generate many sets of correlated random numbers from the covariance matrix representing the fit parameters, calculate mean and sigma for each set, and look at the distributions of mean and sigma behave (e.g., “mean” and “rms” of the calculated mean and sigma, or even the correlation between the mean and sigma).

Root doesn’t have this built in anymore as far as I know (although CERN lib did), but looking at this post http://root.cern.ch/phpBB3//viewtopic.php?f=3&t=7572&hilit=+correlated, I have example code for this.

I find that this is an exteremly powerful technique since it may not always be possible to get at the information you want by re-parameterizing your fit function.

Good luck,
Charles

chetna · March 19, 2011, 6:32am

Hi Charles,
Thanks.I have used this technique too.
So far I have used a number of datasets generated with different random seeds,and then have examined the distribution of means and sigmas.The distributions are gaussian and I have used the standard deviation as the uncertainty.
But this is useful for a large data.When there are limited data then I wouldnot get sufficient number of mean and sigmas to check their behaviour.

daljeet · March 19, 2011, 10:54am

hi ,

Iam trying to run the file fit_vavilov.C using root 5.27 but all the time it is showing ERROR.
ERROR:Invalid type ‘VavilovAccurate’ in declaration of ‘pdf’ fit_vavilov.C
ERROR:Symbol Vavilovaccurate pdf is not defined in current scope.
please help me to solve this problem.

Thanks
daljeet
fit_Vavilov.C (925 Bytes)
data.txt (72 KB)

chetna · March 19, 2011, 11:20am

Hi, use ROOT5.28.

Regards,
Chetna.

daljeet · March 20, 2011, 5:56am

hi
thanks for reply.
Iam now using root 5.28 but still it is showing same problem.

daljeet

moneta · March 21, 2011, 8:05am

Hi,

The Vavilov function is defined in libMathMore, you need to have a ROOT build with it.
You can check if you have MathMore, by doing:

root-config --has-mathmore

If you are building from sources, to get MathMore, you need to have the GSL library installed in your system.

Best Regards

Lorenzo

daljeet · June 18, 2011, 6:46am

hi
In the previous post of Vavilov function it is mention that
The mean and sigma (i.e. p2 and p3) shown in the plot are not true.

Now see,
If p0 is kappa and p1 is beta2,then
mean=gamma-1-ln(p0)-p1.

Now since I am using the locators p2 and p3,so the actual mean should be,
MEAN = mean*p3+p2.

Similarly,
actual sigma should be,
SIGMA = sigmap3, where, sigma=sqrt[(2-p1)/(2p0)]

I could not able to understand that what is gamma here?
please help

Thanks

moneta · June 29, 2011, 3:18pm

See the detailed description of the Vavilov function in
project-mathlibs.web.cern.ch/pro … vilov.html

Lorenzo

kolahalb · December 31, 2012, 6:38am

Hi Daljeet, I know it’s long since you asked this question, and by now, you must have known the answer…But for anyone else following this discussion I am writing: gamma is the Euler’s constant = 0.577…Go through the discussions on Vavilov distribution…

roarsmith · May 29, 2015, 8:40am

[quote=“moneta”]Hi,

you need to add some extra parameters to the Vavilov ( a location and a scale parameters).
Attached is an example doing the fit using the ROOT::Math::Vavilov class from MathMore (available in ROOT 5.28)

Best Regards

Lorenzo[/quote]

I run your code exactly the way it is. I am getting the error

error: no
type named ‘VavilovAccurate’ in namespace 'ROOT::Math’
ROOT::Math::VavilovAccurate pdf;

I am using root 6.02/02

moneta · June 9, 2015, 9:20am

See
project-mathlibs.web.cern.ch/pro … ml#details

for the definitions and the given references for more information

Lorenzo