Skip to content

Commit 9013ccf

Browse files
rhettingerrhettinger
authored
bpo-36546: Add statistics.quantiles() (#12710)
1 parent d437012 commit 9013ccf

5 files changed

Lines changed: 251 additions & 7 deletions

File tree

‎Doc/library/statistics.rst‎

Lines changed: 51 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -48,6 +48,7 @@ or sample.
4848
:func:`median_grouped` Median, or 50th percentile, of grouped data.
4949
:func:`mode` Single mode (most common value) of discrete or nominal data.
5050
:func:`multimode` List of modes (most common values) of discrete or nomimal data.
51+
:func:`quantiles` Divide data into intervals with equal probability.
5152
======================= ===============================================================
5253

5354
Measures of spread
@@ -499,6 +500,53 @@ However, for reading convenience, most of the examples show sorted sequences.
499500
:func:`pvariance` function as the *mu* parameter to get the variance of a
500501
sample.
501502

503+
.. function:: quantiles(dist, *, n=4, method='exclusive')
504+
505+
Divide *dist* into *n* continuous intervals with equal probability.
506+
Returns a list of ``n - 1`` cut points separating the intervals.
507+
508+
Set *n* to 4 for quartiles (the default). Set *n* to 10 for deciles. Set
509+
*n* to 100 for percentiles which gives the 99 cuts points that separate
510+
*dist* in to 100 equal sized groups. Raises :exc:`StatisticsError` if *n*
511+
is not least 1.
512+
513+
The *dist* can be any iterable containing sample data or it can be an
514+
instance of a class that defines an :meth:`~inv_cdf` method.
515+
Raises :exc:`StatisticsError` if there are not at least two data points.
516+
517+
For sample data, the cut points are linearly interpolated from the
518+
two nearest data points. For example, if a cut point falls one-third
519+
of the distance between two sample values, ``100`` and ``112``, the
520+
cut-point will evaluate to ``104``. Other selection methods may be
521+
offered in the future (for example choose ``100`` as the nearest
522+
value or compute ``106`` as the midpoint). This might matter if
523+
there are too few samples for a given number of cut points.
524+
525+
If *method* is set to *inclusive*, *dist* is treated as population data.
526+
The minimum value is treated as the 0th percentile and the maximum
527+
value is treated as the 100th percentile. If *dist* is an instance of
528+
a class that defines an :meth:`~inv_cdf` method, setting *method*
529+
has no effect.
530+
531+
.. doctest::
532+
533+
# Decile cut points for empirically sampled data
534+
>>> data = [105, 129, 87, 86, 111, 111, 89, 81, 108, 92, 110,
535+
... 100, 75, 105, 103, 109, 76, 119, 99, 91, 103, 129,
536+
... 106, 101, 84, 111, 74, 87, 86, 103, 103, 106, 86,
537+
... 111, 75, 87, 102, 121, 111, 88, 89, 101, 106, 95,
538+
... 103, 107, 101, 81, 109, 104]
539+
>>> [round(q, 1) for q in quantiles(data, n=10)]
540+
[81.0, 86.2, 89.0, 99.4, 102.5, 103.6, 106.0, 109.8, 111.0]
541+
542+
>>> # Quartile cut points for the standard normal distibution
543+
>>> Z = NormalDist()
544+
>>> [round(q, 4) for q in quantiles(Z, n=4)]
545+
[-0.6745, 0.0, 0.6745]
546+
547+
.. versionadded:: 3.8
548+
549+
502550
Exceptions
503551
----------
504552

@@ -606,7 +654,7 @@ of applications in statistics.
606654
<http://www.iceaaonline.com/ready/wp-content/uploads/2014/06/MM-9-Presentation-Meet-the-Overlapping-Coefficient-A-Measure-for-Elevator-Speeches.pdf>`_
607655
between two normal distributions, giving a measure of agreement.
608656
Returns a value between 0.0 and 1.0 giving `the overlapping area for
609-
two probability density functions
657+
the two probability density functions
610658
<https://www.rasch.org/rmt/rmt101r.htm>`_.
611659

612660
Instances of :class:`NormalDist` support addition, subtraction,
@@ -649,8 +697,8 @@ of applications in statistics.
649697
For example, given `historical data for SAT exams
650698
<https://blog.prepscholar.com/sat-standard-deviation>`_ showing that scores
651699
are normally distributed with a mean of 1060 and a standard deviation of 192,
652-
determine the percentage of students with scores between 1100 and 1200, after
653-
rounding to the nearest whole number:
700+
determine the percentage of students with test scores between 1100 and
701+
1200, after rounding to the nearest whole number:
654702

655703
.. doctest::
656704

‎Doc/whatsnew/3.8.rst‎

Lines changed: 4 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -337,6 +337,10 @@ Added :func:`statistics.geometric_mean()`
337337
Added :func:`statistics.multimode` that returns a list of the most
338338
common values. (Contributed by Raymond Hettinger in :issue:`35892`.)
339339

340+
Added :func:`statistics.quantiles` that divides data or a distribution
341+
in to equiprobable intervals (e.g. quartiles, deciles, or percentiles).
342+
(Contributed by Raymond Hettinger in :issue:`36546`.)
343+
340344
Added :class:`statistics.NormalDist`, a tool for creating
341345
and manipulating normal distributions of a random variable.
342346
(Contributed by Raymond Hettinger in :issue:`36018`.)

‎Lib/statistics.py‎

Lines changed: 53 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -7,9 +7,9 @@
77
Calculating averages
88
--------------------
99
10-
================== =============================================
10+
================== ==================================================
1111
Function Description
12-
================== =============================================
12+
================== ==================================================
1313
mean Arithmetic mean (average) of data.
1414
geometric_mean Geometric mean of data.
1515
harmonic_mean Harmonic mean of data.
@@ -19,7 +19,8 @@
1919
median_grouped Median, or 50th percentile, of grouped data.
2020
mode Mode (most common value) of data.
2121
multimode List of modes (most common values of data).
22-
================== =============================================
22+
quantiles Divide data into intervals with equal probability.
23+
================== ==================================================
2324
2425
Calculate the arithmetic mean ("the average") of data:
2526
@@ -78,7 +79,7 @@
7879
7980
"""
8081

81-
__all__ = [ 'StatisticsError', 'NormalDist',
82+
__all__ = [ 'StatisticsError', 'NormalDist', 'quantiles',
8283
'pstdev', 'pvariance', 'stdev', 'variance',
8384
'median', 'median_low', 'median_high', 'median_grouped',
8485
'mean', 'mode', 'multimode', 'harmonic_mean', 'fmean',
@@ -562,6 +563,54 @@ def multimode(data):
562563
maxcount, mode_items = next(groupby(counts, key=itemgetter(1)), (0, []))
563564
return list(map(itemgetter(0), mode_items))
564565

566+
def quantiles(dist, *, n=4, method='exclusive'):
567+
'''Divide *dist* into *n* continuous intervals with equal probability.
568+
569+
Returns a list of (n - 1) cut points separating the intervals.
570+
571+
Set *n* to 4 for quartiles (the default). Set *n* to 10 for deciles.
572+
Set *n* to 100 for percentiles which gives the 99 cuts points that
573+
separate *dist* in to 100 equal sized groups.
574+
575+
The *dist* can be any iterable containing sample data or it can be
576+
an instance of a class that defines an inv_cdf() method. For sample
577+
data, the cut points are linearly interpolated between data points.
578+
579+
If *method* is set to *inclusive*, *dist* is treated as population
580+
data. The minimum value is treated as the 0th percentile and the
581+
maximum value is treated as the 100th percentile.
582+
'''
583+
# Possible future API extensions:
584+
# quantiles(data, already_sorted=True)
585+
# quantiles(data, cut_points=[0.02, 0.25, 0.50, 0.75, 0.98])
586+
if n < 1:
587+
raise StatisticsError('n must be at least 1')
588+
if hasattr(dist, 'inv_cdf'):
589+
return [dist.inv_cdf(i / n) for i in range(1, n)]
590+
data = sorted(dist)
591+
ld = len(data)
592+
if ld < 2:
593+
raise StatisticsError('must have at least two data points')
594+
if method == 'inclusive':
595+
m = ld - 1
596+
result = []
597+
for i in range(1, n):
598+
j = i * m // n
599+
delta = i*m - j*n
600+
interpolated = (data[j] * (n - delta) + data[j+1] * delta) / n
601+
result.append(interpolated)
602+
return result
603+
if method == 'exclusive':
604+
m = ld + 1
605+
result = []
606+
for i in range(1, n):
607+
j = i * m // n # rescale i to m/n
608+
j = 1 if j < 1 else ld-1 if j > ld-1 else j # clamp to 1 .. ld-1
609+
delta = i*m - j*n # exact integer math
610+
interpolated = (data[j-1] * (n - delta) + data[j] * delta) / n
611+
result.append(interpolated)
612+
return result
613+
raise ValueError(f'Unknown method: {method!r}')
565614

566615
# === Measures of spread ===
567616

‎Lib/test/test_statistics.py‎

Lines changed: 142 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -3,6 +3,7 @@
33
44
"""
55

6+
import bisect
67
import collections
78
import collections.abc
89
import copy
@@ -2038,6 +2039,7 @@ def test_compare_to_variance(self):
20382039
expected = math.sqrt(statistics.variance(data))
20392040
self.assertEqual(self.func(data), expected)
20402041

2042+
20412043
class TestGeometricMean(unittest.TestCase):
20422044

20432045
def test_basics(self):
@@ -2126,6 +2128,146 @@ def test_special_values(self):
21262128
with self.assertRaises(ValueError):
21272129
geometric_mean([Inf, -Inf])
21282130

2131+
2132+
class TestQuantiles(unittest.TestCase):
2133+
2134+
def test_specific_cases(self):
2135+
# Match results computed by hand and cross-checked
2136+
# against the PERCENTILE.EXC function in MS Excel.
2137+
quantiles = statistics.quantiles
2138+
data = [120, 200, 250, 320, 350]
2139+
random.shuffle(data)
2140+
for n, expected in [
2141+
(1, []),
2142+
(2, [250.0]),
2143+
(3, [200.0, 320.0]),
2144+
(4, [160.0, 250.0, 335.0]),
2145+
(5, [136.0, 220.0, 292.0, 344.0]),
2146+
(6, [120.0, 200.0, 250.0, 320.0, 350.0]),
2147+
(8, [100.0, 160.0, 212.5, 250.0, 302.5, 335.0, 357.5]),
2148+
(10, [88.0, 136.0, 184.0, 220.0, 250.0, 292.0, 326.0, 344.0, 362.0]),
2149+
(12, [80.0, 120.0, 160.0, 200.0, 225.0, 250.0, 285.0, 320.0, 335.0,
2150+
350.0, 365.0]),
2151+
(15, [72.0, 104.0, 136.0, 168.0, 200.0, 220.0, 240.0, 264.0, 292.0,
2152+
320.0, 332.0, 344.0, 356.0, 368.0]),
2153+
]:
2154+
self.assertEqual(expected, quantiles(data, n=n))
2155+
self.assertEqual(len(quantiles(data, n=n)), n - 1)
2156+
self.assertEqual(list(map(float, expected)),
2157+
quantiles(map(Decimal, data), n=n))
2158+
self.assertEqual(list(map(Decimal, expected)),
2159+
quantiles(map(Decimal, data), n=n))
2160+
self.assertEqual(list(map(Fraction, expected)),
2161+
quantiles(map(Fraction, data), n=n))
2162+
# Invariant under tranlation and scaling
2163+
def f(x):
2164+
return 3.5 * x - 1234.675
2165+
exp = list(map(f, expected))
2166+
act = quantiles(map(f, data), n=n)
2167+
self.assertTrue(all(math.isclose(e, a) for e, a in zip(exp, act)))
2168+
# Quartiles of a standard normal distribution
2169+
for n, expected in [
2170+
(1, []),
2171+
(2, [0.0]),
2172+
(3, [-0.4307, 0.4307]),
2173+
(4 ,[-0.6745, 0.0, 0.6745]),
2174+
]:
2175+
actual = quantiles(statistics.NormalDist(), n=n)
2176+
self.assertTrue(all(math.isclose(e, a, abs_tol=0.0001)
2177+
for e, a in zip(expected, actual)))
2178+
2179+
def test_specific_cases_inclusive(self):
2180+
# Match results computed by hand and cross-checked
2181+
# against the PERCENTILE.INC function in MS Excel
2182+
# and against the quaatile() function in SciPy.
2183+
quantiles = statistics.quantiles
2184+
data = [100, 200, 400, 800]
2185+
random.shuffle(data)
2186+
for n, expected in [
2187+
(1, []),
2188+
(2, [300.0]),
2189+
(3, [200.0, 400.0]),
2190+
(4, [175.0, 300.0, 500.0]),
2191+
(5, [160.0, 240.0, 360.0, 560.0]),
2192+
(6, [150.0, 200.0, 300.0, 400.0, 600.0]),
2193+
(8, [137.5, 175, 225.0, 300.0, 375.0, 500.0,650.0]),
2194+
(10, [130.0, 160.0, 190.0, 240.0, 300.0, 360.0, 440.0, 560.0, 680.0]),
2195+
(12, [125.0, 150.0, 175.0, 200.0, 250.0, 300.0, 350.0, 400.0,
2196+
500.0, 600.0, 700.0]),
2197+
(15, [120.0, 140.0, 160.0, 180.0, 200.0, 240.0, 280.0, 320.0, 360.0,
2198+
400.0, 480.0, 560.0, 640.0, 720.0]),
2199+
]:
2200+
self.assertEqual(expected, quantiles(data, n=n, method="inclusive"))
2201+
self.assertEqual(len(quantiles(data, n=n, method="inclusive")), n - 1)
2202+
self.assertEqual(list(map(float, expected)),
2203+
quantiles(map(Decimal, data), n=n, method="inclusive"))
2204+
self.assertEqual(list(map(Decimal, expected)),
2205+
quantiles(map(Decimal, data), n=n, method="inclusive"))
2206+
self.assertEqual(list(map(Fraction, expected)),
2207+
quantiles(map(Fraction, data), n=n, method="inclusive"))
2208+
# Invariant under tranlation and scaling
2209+
def f(x):
2210+
return 3.5 * x - 1234.675
2211+
exp = list(map(f, expected))
2212+
act = quantiles(map(f, data), n=n, method="inclusive")
2213+
self.assertTrue(all(math.isclose(e, a) for e, a in zip(exp, act)))
2214+
# Quartiles of a standard normal distribution
2215+
for n, expected in [
2216+
(1, []),
2217+
(2, [0.0]),
2218+
(3, [-0.4307, 0.4307]),
2219+
(4 ,[-0.6745, 0.0, 0.6745]),
2220+
]:
2221+
actual = quantiles(statistics.NormalDist(), n=n, method="inclusive")
2222+
self.assertTrue(all(math.isclose(e, a, abs_tol=0.0001)
2223+
for e, a in zip(expected, actual)))
2224+
2225+
def test_equal_sized_groups(self):
2226+
quantiles = statistics.quantiles
2227+
total = 10_000
2228+
data = [random.expovariate(0.2) for i in range(total)]
2229+
while len(set(data)) != total:
2230+
data.append(random.expovariate(0.2))
2231+
data.sort()
2232+
2233+
# Cases where the group size exactly divides the total
2234+
for n in (1, 2, 5, 10, 20, 50, 100, 200, 500, 1000, 2000, 5000, 10000):
2235+
group_size = total // n
2236+
self.assertEqual(
2237+
[bisect.bisect(data, q) for q in quantiles(data, n=n)],
2238+
list(range(group_size, total, group_size)))
2239+
2240+
# When the group sizes can't be exactly equal, they should
2241+
# differ by no more than one
2242+
for n in (13, 19, 59, 109, 211, 571, 1019, 1907, 5261, 9769):
2243+
group_sizes = {total // n, total // n + 1}
2244+
pos = [bisect.bisect(data, q) for q in quantiles(data, n=n)]
2245+
sizes = {q - p for p, q in zip(pos, pos[1:])}
2246+
self.assertTrue(sizes <= group_sizes)
2247+
2248+
def test_error_cases(self):
2249+
quantiles = statistics.quantiles
2250+
StatisticsError = statistics.StatisticsError
2251+
with self.assertRaises(TypeError):
2252+
quantiles() # Missing arguments
2253+
with self.assertRaises(TypeError):
2254+
quantiles([10, 20, 30], 13, n=4) # Too many arguments
2255+
with self.assertRaises(TypeError):
2256+
quantiles([10, 20, 30], 4) # n is a positional argument
2257+
with self.assertRaises(StatisticsError):
2258+
quantiles([10, 20, 30], n=0) # n is zero
2259+
with self.assertRaises(StatisticsError):
2260+
quantiles([10, 20, 30], n=-1) # n is negative
2261+
with self.assertRaises(TypeError):
2262+
quantiles([10, 20, 30], n=1.5) # n is not an integer
2263+
with self.assertRaises(ValueError):
2264+
quantiles([10, 20, 30], method='X') # method is unknown
2265+
with self.assertRaises(StatisticsError):
2266+
quantiles([10], n=4) # not enough data points
2267+
with self.assertRaises(TypeError):
2268+
quantiles([10, None, 30], n=4) # data is non-numeric
2269+
2270+
21292271
class TestNormalDist(unittest.TestCase):
21302272

21312273
# General note on precision: The pdf(), cdf(), and overlap() methods

‎Misc/NEWS.d/next/Library/2019-04-06-14-23-00.bpo-36546.YXjbyY.rst‎

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1 @@
1+
Add statistics.quantiles()

0 commit comments

Comments
 (0)