1
$\begingroup$

I am trying to compare various PDFs for goodness of fit using 100,000 points of simulated data. I use some of Mathematica's built-in PDFs but I also need to generate my own. I have run into a couple of problems with the statistical package that I can't make sense of. For instance, when I use DistributionFitTest, even with one of the built in PDF's, the results don't always make sense. Here is an example where I have data that obviously fits the test PDF pretty well but the P-value is low and furthermore, the p-value does not change when I vary the level of significance. Here is the code and results along with a smooth histogram plot of the data compared to the test PDF:

testdata = Flatten[Sort[
Import["C:\\Users\\Missy\\Documents\\Wolfram \Mathematica\\PDFs\\WOSData\\Cn2_ 5.00E-014_path _lngth _ 100m_no _of \_scrns _ 3_sptl _smpl _\\3PixelAperturefull.csv", "Data"]]]
H1 = DistributionFitTest[testdata, LogNormalDistribution[2.18651, 0.138736], "HypothesisTestData", SignificanceLevel -> .01];
H2 = DistributionFitTest[testdata, LogNormalDistribution[2.18651, 0.138736], "HypothesisTestData", SignificanceLevel -> .1];
H3 = DistributionFitTest[testdata,LogNormalDistribution[2.18651, 0.138736], "HypothesisTestData", SignificanceLevel -> .99];
{H1["TestDataTable", All], H2["TestDataTable", All], H3["TestDataTable", All]}

Now the results and plots: Goodness of Fit

Another thing that doesn't make sense to me is when I define a PDF and then try to see the CDF of it. Below is an example. I define the PDF, then plot the PDF, Mathematica calculated CDF, and correct CDF which is integral of original function. Why doesn't it work?

GG[a_, b_, x_] := (2. (a b)^((a + b)/2))/(Gamma[a] Gamma[b])x^((a + b)/2 - 1) BesselK[a - b, 2. Sqrt[a b x]]
ClearAll[EstGG];
EstGG := ProbabilityDistribution[GG[7.174, 1.233, x], {x, 0., Infinity}]
Plot[{PDF[EstGG, x], CDF[EstGG, x], NIntegrate[GG[7.174, 1.233, xp], {xp, 0, x}]}, {x, 0, 10}, PlotRange -> Full, PlotStyle -> {Blue, Green, Red},PlotLegends -> {"PDF", "Mathematica CDF",  "Numerical Integration for CDF"}]

And here is the plot: CDF

$\endgroup$
  • $\begingroup$ For your second issue: use EstGG = instead of EstGG :=. $\endgroup$ – JimB May 15 at 17:01
  • 2
    $\begingroup$ For your first issue: No, it is not obvious that the fit is good and will result in a non-significant P-value. You have 100,000 data points. If there are even small deviations from the hypothesized distribution, a large sample size will tend to result in small P-values as the test concerns deviations from an exact match rather than what you or I might consider as a "practically the same" fit. You need to know what kind and size of differences are important. Such practical differences (especially the size of the difference) are not really considered in a statistical test. $\endgroup$ – JimB May 15 at 17:11
2
$\begingroup$

Consider the following:

  1. Changing the level of significance does not change the observed level of significance as those are different critters.
  2. A P-value is the probability of obtaining at least as extreme a value given that the sample comes from the hypothesized distribution.
  3. A P-value is NOT the probability that the sample comes from the hypothesized distribution.
  4. A P-value can become small (1) by chance, (2) if the data comes from
    some other very different distribution, and/or (3) if the data comes
    from a nearly identical distribution but with a huge sample size so that deviations very small by subject matter considerations result in statistical significance.
  5. Statistical significance is not subject matter significance.

OK. Enough sermonizing.

Suppose the data comes from a mixture distribution where 99% of the time it comes from the log normal distribution that you used and 1% of the time it comes from a slightly different log normal distribution.

d = MixtureDistribution[{0.99, 0.01}, {LogNormalDistribution[2.18651, 0.138736], 
    LogNormalDistribution[2, .3]}];

By design we now have two different log normal distributions but they look very similar (no data yet, these are the actual pdf's):

Plot[{PDF[LogNormalDistribution[2.18651, 0.138736], z], PDF[d, z]}, {z, 0, 25},
  PlotLegends -> {"Log normal distribution", "Mixture distribution"},
  PlotStyle -> {Thickness[0.02], {Orange, Dashed}}, 
  PlotRange -> {{0, 25}, {0, 0.35}}]

True pdfs for log normal and mixture distributions

Now take a sample from the mixture distribution and test for differences from the log normal distribution:

SeedRandom[12345];
testdata = RandomVariate[d, 100000];
Show[SmoothHistogram[testdata, PlotRange -> {{0, 25}, {0, 0.35}}, PlotStyle -> Thickness[0.02]], 
 Plot[PDF[LogNormalDistribution[2.18651, 0.138736], z], {z, 0, 20}, 
  PlotStyle -> {Orange, Dashed}], PlotRange -> All]

Estimated density vs hypothesized density

H1 = DistributionFitTest[testdata, 
   LogNormalDistribution[2.18651, 0.138736], "HypothesisTestData", 
   SignificanceLevel -> .01];
H1["TestDataTable", All]

Table of P-values

So we know the distributions are definitely different but they look very similar. And with a sample size of 100,000 we get very small P-values.

Relying solely on P-values is not a recommended practice to be able to make good decisions. Some subject matter criteria needs to be involved at some point. (Otherwise why would we need subject matter experts? You want me making medical decisions just because I can calculate P-values from datasets I don't understand?)

$\endgroup$
  • $\begingroup$ Thank you. These comments have been very helpful. I was naively expecting the p-value to remain fairly consistent as I changed the number of samples - but it did not and with this explanation and some research I now understand why. Also, my expectations of changing the level of significance were incorrect. The comment regarding changing from := to = in the PDF definition was spot on. Thanks again for response. Long time follower of stack exchange but first time posting. $\endgroup$ – Melissa May 17 at 9:34

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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