First up, I'm quite new to Mathematica so any hints on better code would be greatly appreciated.
I have some histogram frequency data from an unknown distribution that I'm trying to fit. Here's the histogram and how I draw it:
Needs["Histograms`"]
(* Data in the form "rank count"*)
raw = ReadList["Desktop/data.txt", {Number, Number}];
(* I have some missing ranks that I want to set as having a count of 0 *)
fd = Table[0, {i, Max[raw[[All, 1]]]}];
Do[fd[[raw[[i, 1]]]] = raw[[i, 2]], {i, Length[raw]}]
Histogram[Log[fd], FrequencyData -> True]
Needs["Histograms"]
seems to be deprecated but I couldn't find a nice way of plotting frequency data without it.
The ranks are actually counts themselves, this is a histogram of the frequencies at which I've observed X number of things. Does that make sense? I'm slightly concerned that I'm confusing myself here :)
Now I have the data and the log plot seems to show a nice continuous curve I thought I could find a line to fit it. I followed these instructions: FindFIt with BinCounts, for using FindFit over frequency data, so far so good. Let's try a power law distribution:
centers = MovingAverage[raw[[All, 2]], 2];
counts = raw[[All, 1]];
centered = Table[{centers[[i]], counts[[i]]}, {i, Length[centers]}];
xmin = 1;
model = ((a - 1)/xmin)*(x/xmin)^(-a)
pars = FindFit[centered, model, {a}, x]
nlm = NonlinearModelFit[centered, model, {a}, x];
nlm[{"BestFit", "ParameterTable"}]
I'm not particularly expecting a power law distribution to work but it demonstrates my method. Here's what I get:
(-1+a) x^-a
{a->1.16913}
| Estimate | Standard Error | t-Statistic | P-Value
----------------------------------------------------------
a | 1.16913 | 122731. | 9.52597*10^-6 | 0.999992
So, my questions. Is what I'm doing correct? If I managed to find a model that describes my data will FindFit estimate the parameters for me? Can anyone help with what that distribution might be? I've tried (grasping at straws):
- Zipf's law
- Power law
- Gumbel
- Laplace
- Frechet
Some have been reported as having a very good fit (P-value ridiculously low) but the lines don't really match up and when I plug in some numbers I get really bad results.
I feel like I'm being hopelessly naive in trying to do this. I've had my head stuck in these numbers for so long I don't really know what's going on :)