# Computing the mode of data sampled from a continuous distribution

What are the best methods for fitting the 'mode' of data sampled from a continuous distribution?

Since the mode is technically undefined (right?) for a continuous distribution, I'm really asking 'how do you find the most common value'?

If you assume the parent distribution is gaussian, you could bin the data and find say the mode is the bin location with the greatest counts. However, how do you determine the bin size? Are there robust implementations available? (i.e., robust to outliers). I use python/scipy/numpy, but I can probably translate R without too much difficulty.

-
I'm not sure if the mode is technically defined this way, but the global mode of a continuous distribution is usually taken to mean the point with the highest density. – Macro Dec 18 '11 at 0:50
@Macro - that's helpful. You can then read my question as, 'What are the best methods to determine the (peak) density?' – keflavich Dec 18 '11 at 0:52
Maybe fit a kernel density estimate for your data and estimate the mode as the peak of that? This seems like a reasonable approach but I'm not familiar with the literature on this problem. – Macro Dec 18 '11 at 0:53
If you don't assume the parent distribution is gaussian, is it still possible to bin the data and take the mode to be the bin location with the largest count? Why or why not? On a more serious note, why not find the deciles $x_0=x_{\min},x_1,x_2,\ldots,x_9,x_{10}=x_{\max}$ so that $10\%$ of the samples are in the interval $x_{i+1}-x_i$, and so it is likely that the mode is in the shortest interdecile interval $\min_{1 \leq j \leq 10} x_{j+1}-x_j$? Then take the bin size to be, say, one-fourth of this shortest interdecile interval. – Dilip Sarwate Dec 18 '11 at 1:53
What assumptions can you make about the parent distribution, keflavich? If they are parametric, it's best to estimate the parameters and then estimate the mode from those parameters. (E.g., the sample mean estimates the mode of a normal distribution.) If not, binning can be a poor method. Instead, a sequence of kernel estimators with varying halfwidth can be used to provide a sequence of estimators; typically, if the underlying distribution is unimodal, the modes of the kernel smooths will appear to converge towards a unique mode as the halfwidths get large and that can be your estimate. – whuber Dec 18 '11 at 16:53

In R, applying the method that isn't based on parametric modelling of the underlying distribution and uses the default kernel estimator of density to 10000 gamma distributed variables:

x <- rgamma(10000, 2, 5)
z <- density(x)
plot(z) # always good to check visually
z$x[z$y==max(z$y)]  returns 0.199 which is the value of x estimated to have the highest density (the density estimates are stored as "z$y").

-
The only thing I would do differently from that is use a different bandwidth. The default bandwidth for density() is not particularly good. density(x,bw="SJ") is better. Even better would be to use a bandwidth designed for mode estimation. See sciencedirect.com/science/article/pii/0167715295000240 for some discussion. – Rob Hyndman Feb 17 '12 at 12:02

See the modeest package for R.

-