Understanding the Cullen and Frey plot

Question

I would like to figure out which distribution fits my data best.

Here is the histogram of my data :

I used the fitdistrplus package in R to try to find the best fit for my data. To get an idea of which family distribution to fit, I did this:

library(fitdistrplus)
descdist(my_data, discrete=FALSE, boot=500)

I get this skewness-kurtosis plot :

With these summary statistics :

min: 23 max: 1989
median: 184
mean: 228.8346
estimated sd: 165.6273
estimated skewness: 1.706379
estimated kurtosis: 11.31023

So apparently no distribution is a good candidate for the data. How to interpret this plot ? Does this mean that my data is a mixture of several distributions ?

EDIT :

This distribution represents DNA fragment lengths obtained from an experiment. My aim is to be able to simulate the result of this experiment by simuating the resulting fragments. (i.e a simulated fragment is defined by two positions in the genome seperated by a distance D). I assume that the fragment length distribution observed from the real experiment can be described by a density function or a mxiture of densities. I am looking fo the best function from which I could sample values of D for my simulations.

To note, I am using only a sub sample to fit a distribution. We generate millions of fragments. I work on subsample of 500.000 fragments.

Answer 1

This plot used to be commonly called a Pearson plot (it also had several other names), though sometimes with skewness rather than its square being plotted. It was used long before Cullen and Frey wrote about it (a fact they clearly acknowledge in their text, though their own mention of having seen it in a book written in the late 60s still considerably underestimates its age).

The aim of such a plot was to help identify a suitable Pearson distribution.

The Cullen and Frey version of the plot doesn't show all the Pearson family on the plot; you can't see from that plot whether the skewness and kurtosis would correspond to that of a Pearson IV or VI distribution because they leave the dividing line off the plot (which corresponds to a shifted and scaled inverse Gamma)

By transforming (squeezing and rotating) the plot to fit the one here it turns out that it's in the region of a Pearson IV, but you can see from the histogram that the skewness and kurtosis not a sufficient way to summarize the distribution -- no Pearson IV distribution is shaped like that; nor are a couple of other candidates that would correspond to that approximate region.

Another thing to note is that the sample kurtosis tends to underestimate the population kurtosis, and that selection by matching third and fourth cumulants isn't usually an especially good way to go about choosing a model.

Indeed, it is likely that no simple, commonly used distribution will fit very well. You might get an adequate fit with a mixture (as you suggest); I'd expect at least 4-5 components from some suitable family might be required.

However, there are few applications where it's really necessary to identify a distributional form like this -- you'd be much better to explain what you'd be using such a distribution for, because there's probably something better that you can do than this.

Answer 2

Your data looks like a mixture to me: there seems to be one component with a mean fragment length of 100 nt, one around 200 nt and another with a mean fragment length around 300 nt (you can see 'bumps' in the histogram).

Is there anything about how the library was prepared that would explain why there is more than one component in the mixture?

I would fit a mixture of 3 Gaussians to the data. I use the R package mixtools. Edit: to check goodness of fit, you might try this function: https://rdrr.io/cran/AdaptGauss/man/Chi2testMixtures.html