# What is a good algorithm for estimating the median of a huge read-once data set? [duplicate]

I'm looking for a good algorithm (meaning minimal computation, minimal storage requirements) to estimate the median of a data set that is too large to store, such that each value can only be read once (unless you explicitly store that value). There are no bounds on the data that can be assumed.

Approximations are fine, as long as the accuracy is known.

Any pointers?

• You can try to find a median based on grouped frequency distribution, here is some details – ryfm Jul 20 '10 at 20:12
• Perhaps, asking on Stackoverflow may get better answers. – user28 Jul 21 '10 at 15:33
• @Srikant:> it's a pretty active area of research in statistics :) The solution closest to the lower theoretical bounds in terms of storage involve some pretty clever probability constructs as well. All in all i was surprised when i first looked unto it a couple of month ago; there is more stats here than meets the eye. – user603 Oct 9 '10 at 14:21

Could you group the data set into much smaller data sets (say 100 or 1000 or 10,000 data points) If you then calculated the median of each of the groups. If you did this with enough data sets you could plot something like the average of the results of each of the smaller sets and this woul, by running enough smaller data sets converge to an 'average' solution.

• This is interesting, and where some statistical advice could come in! Assume in total I've got (say) 500,000 i.i.d. points and I look at groups of (say) 1,000 of them, and calculate the median of each group. Now I've got 500 medians. Is there theory that could allow me to calculate a confidence interval for the overall median based on these 500 medians? – PeterR Jul 21 '10 at 17:10
• So, according to a long lost colleague, the best apropoach seems to be Chiranjeeb Buragohain and Subhash Suri. Quantiles on Streams. cs.ucsb.edu/~suri/psdir/ency.pdf I also like Ian's approach, as these medians of smaller data sets will converge to a normal distribution, and so I can form conf intervals for the medians. – PeterR Jul 27 '10 at 15:10

How about something like a binning procedure? Assume (for illustration purposes) that you know that the values are between 1 and 1 million. Set up N bins, of size S. So if S=10000, you'd have 100 bins, corresponding to values [1:10000, 10001:20000, ... , 990001:1000000]

Then, step through the values. Instead of storing each value, just increment the counter in the appropriate bin. Using the midpoint of each bin as an estimate, you can make a reasonable approximation of the median. You can scale this to as fine or coarse of a resolution as you want by changing the size of the bins. You're limited only by how much memory you have.

Since you don't know how big your values may get, just pick a bin size large enough that you aren't likely to run out of memory, using some quick back-of-the-envelope calculations. You might also store the bins sparsely, such that you only add a bin if it contains a value.

Edit:

The link ryfm provides gives an example of doing this, with the additional step of using the cumulative percentages to more accurately estimate the point within the median bin, instead of just using midpoints. This is a nice improvement.

• The problem with the binning approach is that we do not have a good upper bound for the data, and so the midpoint for the largest bin would have to be huge. So, we'd need a huge number of bins (not enough memory for that), or have pretty wide bins (which would then lead to a fairly inaccurate answer.) And the data is not very sparse. – PeterR Jul 21 '10 at 15:13
• Since you are only interested in the median why couldn't you make the bins wider at higher values of your variable? – russellpierce Jul 27 '10 at 0:23
• drknexus - because we don't know what the largest bin should be. – PeterR Jul 27 '10 at 14:58
• Do you have any intuition as to what the range will be? If you're fairly sure that over half of the answers will be below number N, then you can make your last bin as large as you want. Maybe your last bin is all numbers greater than 1 trillion - would that be high enough? With the amount of memory in modern systems you can store a LOT of bins and achieve fairly high resolution. In terms of data structures, we're not talking anything fancy and memory intensive here. – chrisamiller Jul 27 '10 at 22:17
• Any intuition? yes. And your approach could work in general. However, in this case we can not have a lot of memory/computation. It is in a networking application where the device could see tens of thousands of items per second, and have VERY little processing left over for this purpose. Not the ideal/typical scenario, I know, but that is what makes it interesting! – PeterR Jul 28 '10 at 16:43

I re-direct you to my answer to a similar question. In a nutshell, it's a read once, 'on the fly' algorithm with $O(n)$ worst case complexity to compute the (exact) median.

I implemented the P-Square Algorithm for Dynamic Calculation of Quantiles and Histograms without Storing Observations in a neat Python module I wrote called LiveStats. It should solve your problem quite effectively.

• (+1) Thanks for stopping by and providing that link, Sean! – whuber Mar 9 '13 at 22:23

The Rivest-Tarjan-Selection algorithm (sometimes also called the median-of-medians algorithm) will let you compute the median element in linear-time without any sorting. For large data sets this is can be quite a bit faster than log-linear sorting. However, it won't solve your memory storage problem.

I've never had to do this, so this is just a suggestion.

I see two (other) possibilities.

Half data

1. Load in half the data and sort
2. Next read in the remaining values and compare against the your sorted list.
1. If the new value is larger, discard it.
2. else put the value in the sorted list and removing the largest value from that list.

Sampling distribution

The other option, is to use an approximation involving the sampling distribution. If your data is Normal, then the standard error for moderate n is:

1.253 * sd / sqrt(n)

To determine the size of n that you would be happy with, I ran a quick Monte-Carlo simulation in R

n = 10000
outside.ci.uni = 0
outside.ci.nor = 0
N=1000
for(i in 1:N){
#Theoretical median is 0
uni = runif(n, -10, 10)
nor  = rnorm(n, 0, 10)

if(abs(median(uni)) > 1.96*1.253*sd(uni)/sqrt(n))
outside.ci.uni = outside.ci.uni + 1

if(abs(median(nor)) > 1.96*1.253*sd(nor)/sqrt(n))
outside.ci.nor = outside.ci.nor + 1
}

outside.ci.uni/N
outside.ci.nor/N


For n=10000, 15% of the uniform median estimates were outside the CI.

• The data set is potentially too big to read in half of it...it is in a networking context where the device doing the processing can see tens of thousands of items per second, and probably has enough memory to store only a few hundred. Also the data is definitely not Gaussian. In fact it does not fit well to any of the common distributions. – PeterR Jul 21 '10 at 15:17

The Remedian Algorithm (PDF) gives a one-pass median estimate with low storage requirements and well defined accuracy.

The remedian with base b proceeds by computing medians of groups of b observations, and then medians of these medians, until only a single estimate remains. This method merely needs k arrays of size b (where n = b^k)...

If the values you are using are within a certain range, say 1 to 100000, you can efficiently compute the median on an extremely large number of values (say, trillions of entries), with an integer bucket (this code taken from BSD licensed ea-utils/sam-stats.cpp)

class ibucket {
public:
int tot;
vector<int> dat;
ibucket(int max) {dat.resize(max+1);tot=0;}

int operator[] (int n) const {
assert(n < size());
int i;
for (i=0;i<dat.size();++i) {
if (n < dat[i]) {
return i;
}
n-=dat[i];
}
}

void push(int v) {
assert(v<dat.size());
++dat[v];
++tot;
}
};

template <class vtype>
double quantile(const vtype &vec, double p) {
int l = vec.size();
if (!l) return 0;
double t = ((double)l-1)*p;
int it = (int) t;
int v=vec[it];
if (t > (double)it) {
return (v + (t-it) * (vec[it+1] - v));
} else {
return v;
}
}

• Also, this can be extended to using a finite number of bins for real-time medians, etc. – Erik Aronesty Aug 20 '14 at 15:27