How to plot a Histogram and compare with a particular probability density function PDF

Question

I would like to know how to plot a make a comparison between a histogram and a set of random data, I have the following probability density function (PDF) which is given for this problem. I can generate the histogram but I do not know how to add the PDF function and make the plot both together.

I have seen some question and problem using RandomVariate function but it does not work. This is the code.

Clear["Global`*"];  Remove["Global`*"];
 L = 10; Ns = 10^4;
 x = Table[RandomReal[L], {x, 0,Ns}];  
 y = Table[RandomReal[L], {y, 0, Ns}];
 z = Table[RandomReal[L], {z, 0, Ns}]; 
 sf = Abs[Sin[x^2  + y^2] + Cos[x y z]];
 Histogram[sf, {0.25}, "Probability"]

The (PDF) is f = k^2/(k^2 + 1)^(5/2) The goal is to obtain a similar plot like this, , where the line represent in this case f.

Answer 1

Use Show:

Show[
Histogram[sf, {0.25}, "Probability"],
Plot[f, {k, 0, 2}]
]

Answer 2

The OP writes: "It is important to know what is the analytical form of the function (plotted in blue color) given, which is approximated by the Histogram"

It might be difficult to find the analytical form of sf. But we can approximate the analytical form quite easily by noting that your data appears exponentially downward sloping, but bounded on (0,2). So a natural choice would be a truncated $\text{Gamma}(a,b)$ model, truncated above at 2. This model will have pdf

$$g(s) = c s^{a-1} \exp{(-\frac{s}{b})}$$

where constant c = 1/(b^a*Gamma[a]) - Gamma[a, 2/b].

Parameters $a$ and $b$ can then be found using maximum likelihood.

Here is a plot of your histogram and the fitted truncated Gamma model when $a = 0.805$ and $b = 4.47$:

...which seems to work nicely. The fitted density is then approximately something simple like:

0.56 Exp[-0.224 s] s^(-.2)

Answer 3

This is a specific function which can be modified for any PDF you like. (It looks horrible because i wrote it.) Set PrintDistributionAndFitOpt -> True to see the plot. Otherwise it just returns the fit.

Clear[GetGaussianFitToDataDistribution];
Options[GetGaussianFitToDataDistribution] = {
UseNumberBinsOpt -> 50
    ,PrintDistributionAndFitOpt -> False
};
GetGaussianFitToDataDistribution[InData_List, OptionsPattern[]] := Module[{
    UseNumberBins, PrintDistributionAndFit
    , DataScore, MaxBin, MinBin, BinWidth, UseBins
    , wsigma, wmu, x, wFit
},
{UseNumberBins, PrintDistributionAndFit} = OptionValue[{UseNumberBinsOpt, PrintDistributionAndFitOpt}];
If[EvenQ[UseNumberBins], UseNumberBins++;];
If[Length[InData] < UseNumberBins
, Print["Warning in GetGaussianFitToDataDistribution: Less data than bins."];
];
(* Determine the bins *)
MaxBin = Ceiling[Max[InData]];
MinBin = Floor[Min[InData]];
BinWidth = N[(MaxBin - MinBin)/UseNumberBins];
UseBins = Range[MinBin + 0.5*BinWidth, MaxBin, BinWidth];
(* Generate data to be fit using BinCounts
- normalize to 1.0 or the optimization gets confused
*)
DataScore = Transpose[{UseBins
, (N[#] &) /@ BinCounts[InData, {MinBin, MaxBin, BinWidth}]
}];
DataScore = ({#[[1]], #[[2]]/Max[DataScore[[All, 2]]]} &) /@ 
DataScore;
(*Print[DataScore];*)
(* Find the Best Fit 
- Uses Mean and StDev for estimates of the partameters
*)
wFit = NonlinearModelFit[DataScore
, PDF[NormalDistribution[wmu, wsigma], x]/
PDF[NormalDistribution[wmu, wsigma], wmu]
, {{wmu, Mean[InData]}
, {wsigma, StandardDeviation[InData]}}, x];
If[PrintDistributionAndFit, Print[Show[{
ListPlot[DataScore, PlotStyle -> Darker[Red], PlotRange -> All
, PlotLegends -> Placed[{"Original Data"}, {0.85, 0.85}]
]
, Plot[wFit[x], {x, MinBin, MaxBin}
, PlotLegends -> Placed[{"Fitted Gaussian"}, {0.85, 0.75}]
]
}
, PlotLabel -> ToString[{wmu, wsigma} /. wFit["BestFitParameters"]]
    <> "\t Samples = " <> ToString[Length[InData]]
, FrameLabel -> {"Bin Mid-Line", "Relative Frequency"}
]];
];
{wmu, wsigma} /. wFit["BestFitParameters"]
];
(*
TestData=RandomVariate[NormalDistribution[1000.,20.],10^5];
GetGaussianFitToDataDistribution[TestData
,PrintDistributionAndFitOpt\[Rule]True
,UseNumberBinsOpt\[Rule]50
]
*)