Statistik: Kombinationen in Python

122

Ich brauche combinatorials (nCr) in Python zu berechnen , aber die Funktion nicht finden können , das zu tun in math, numpyoder stat Bibliotheken. So etwas wie eine Funktion des Typs:

comb = calculate_combinations(n, r)

Ich brauche die Anzahl der möglichen Kombinationen, nicht die tatsächlichen Kombinationen, also itertools.combinationsinteressiert mich das nicht.

Schließlich möchte ich die Verwendung von Fakultäten vermeiden, da die Zahlen, für die ich die Kombinationen berechnen werde, zu groß werden können und die Fakultäten monströs sein werden.

Dies scheint eine WIRKLICH einfach zu beantwortende Frage zu sein, aber ich bin in Fragen über das Generieren aller tatsächlichen Kombinationen ertrunken, was ich nicht will.

Morlock
quelle

Antworten:

121

Siehe scipy.special.comb (scipy.misc.comb in älteren Versionen von scipy). Wenn exactFalse ist, wird die Gammaln-Funktion verwendet, um eine gute Präzision zu erzielen, ohne viel Zeit in Anspruch zu nehmen. Im genauen Fall wird eine Ganzzahl mit beliebiger Genauigkeit zurückgegeben, deren Berechnung möglicherweise lange dauert.

Jouni K. Seppänen
quelle
5
scipy.misc.combist zugunsten scipy.special.combseit Version veraltet 0.10.0.
Dilawar
119

Warum schreibst du es nicht selbst? Es ist ein Einzeiler oder so:

from operator import mul    # or mul=lambda x,y:x*y
from fractions import Fraction

def nCk(n,k): 
  return int( reduce(mul, (Fraction(n-i, i+1) for i in range(k)), 1) )

Test - Drucken des Pascalschen Dreiecks:

>>> for n in range(17):
...     print ' '.join('%5d'%nCk(n,k) for k in range(n+1)).center(100)
...     
                                                   1                                                
                                                1     1                                             
                                             1     2     1                                          
                                          1     3     3     1                                       
                                       1     4     6     4     1                                    
                                    1     5    10    10     5     1                                 
                                 1     6    15    20    15     6     1                              
                              1     7    21    35    35    21     7     1                           
                           1     8    28    56    70    56    28     8     1                        
                        1     9    36    84   126   126    84    36     9     1                     
                     1    10    45   120   210   252   210   120    45    10     1                  
                  1    11    55   165   330   462   462   330   165    55    11     1               
               1    12    66   220   495   792   924   792   495   220    66    12     1            
            1    13    78   286   715  1287  1716  1716  1287   715   286    78    13     1         
         1    14    91   364  1001  2002  3003  3432  3003  2002  1001   364    91    14     1      
      1    15   105   455  1365  3003  5005  6435  6435  5005  3003  1365   455   105    15     1   
    1    16   120   560  1820  4368  8008 11440 12870 11440  8008  4368  1820   560   120    16     1
>>> 

PS. bearbeitet, um int(round(reduce(mul, (float(n-i)/(i+1) for i in range(k)), 1))) durch zu ersetzen, int(reduce(mul, (Fraction(n-i, i+1) for i in range(k)), 1))damit es nicht für große N / K irrt

Nas Banov
quelle
26
+1 für den Vorschlag, etwas Einfaches zu schreiben, für die Verwendung von
Reduce
6
-1, weil diese Antwort falsch ist: print Fakultät (54) / (Fakultät (54 - 27)) / Fakultät (27) == nCk (54, 27) gibt False.
Robert King
3
@robertking - Ok, du warst sowohl kleinlich als auch technisch korrekt. Was ich tat, war als Beispiel dafür gedacht, wie man seine eigene Funktion schreibt; Ich wusste, dass es aufgrund der Gleitkommapräzision nicht genau genug für N und K ist. Aber wir können das beheben - siehe oben, jetzt sollte es nicht für große Zahlen irren
Nas Banov
9
Dies wäre wahrscheinlich schnell in Haskell, aber leider nicht in Python. Es ist eigentlich ziemlich langsam im Vergleich zu vielen anderen Antworten, zB @Alex Martelli, JF Sebastian und meinen eigenen.
Todd Owen
9
Für Python 3 musste ich auch from functools import reduce.
Velizar Hristov
52

Eine schnelle Suche nach Google-Code ergibt (es wird eine Formel aus der Antwort von @Mark Byers verwendet ):

def choose(n, k):
    """
    A fast way to calculate binomial coefficients by Andrew Dalke (contrib).
    """
    if 0 <= k <= n:
        ntok = 1
        ktok = 1
        for t in xrange(1, min(k, n - k) + 1):
            ntok *= n
            ktok *= t
            n -= 1
        return ntok // ktok
    else:
        return 0

choose()ist 10 mal schneller (getestet an allen 0 <= (n, k) <1e3 Paaren) als scipy.misc.comb()wenn Sie eine genaue Antwort benötigen.

def comb(N,k): # from scipy.comb(), but MODIFIED!
    if (k > N) or (N < 0) or (k < 0):
        return 0L
    N,k = map(long,(N,k))
    top = N
    val = 1L
    while (top > (N-k)):
        val *= top
        top -= 1
    n = 1L
    while (n < k+1L):
        val /= n
        n += 1
    return val
jfs
quelle
Eine schöne Lösung, die kein Paket erfordert
Edward Newell
2
Zu
Ihrer Information
Diese chooseFunktion sollte viel mehr Up-Votes haben! Python 3.8 hat math.comb, aber ich musste Python 3.6 für eine Herausforderung verwenden, und keine Implementierung lieferte genaue Ergebnisse für sehr große Ganzzahlen. Dieser macht und macht es schnell!
Am
42

Wenn Sie genaue Ergebnisse und Geschwindigkeit wünschen , versuchen Sie es mit gmpy - gmpy.combsollte genau das tun, wonach Sie fragen, und es ist ziemlich schnell (natürlich gmpybin ich als Originalautor voreingenommen ;-).

Alex Martelli
quelle
6
In der Tat gmpy2.comb()ist 10-mal schneller als choose()aus meiner Antwort für den Code: for k, n in itertools.combinations(range(1000), 2): f(n,k)Wo f()ist entweder gmpy2.comb()oder choose()auf Python 3.
jfs
Da Sie der Autor des Pakets sind, lasse ich Sie den defekten Link reparieren, damit er auf die richtige Stelle zeigt ...
Selten
@SeldomNeedy, der Link zu code.google.com ist ein richtiger Ort (obwohl sich die Site jetzt im Archivierungsmodus befindet). Von dort aus ist es natürlich einfach, den Github-Speicherort github.com/aleaxit/gmpy und den PyPI- Speicherort pypi.python.org/pypi/gmpy2 zu finden , da er mit beiden verknüpft ist! -)
Alex Martelli
@ AlexMartelli Entschuldigung für die Verwirrung. Auf der Seite wird ein 404 angezeigt, wenn Javascript (selektiv) deaktiviert wurde. Ich denke, das soll Schurken-AIs davon abhalten, archivierte Google Code Project-Quellen so einfach einzubeziehen.
Selten Needy
28

Wenn Sie ein genaues Ergebnis wünschen, verwenden Sie sympy.binomial. Es scheint zweifellos die schnellste Methode zu sein.

x = 1000000
y = 234050

%timeit scipy.misc.comb(x, y, exact=True)
1 loops, best of 3: 1min 27s per loop

%timeit gmpy.comb(x, y)
1 loops, best of 3: 1.97 s per loop

%timeit int(sympy.binomial(x, y))
100000 loops, best of 3: 5.06 µs per loop
Jim Garrison
quelle
22

Eine wörtliche Übersetzung der mathematischen Definition ist in vielen Fällen völlig ausreichend (wobei zu beachten ist, dass Python automatisch eine Arithmetik mit großen Zahlen verwendet):

from math import factorial

def calculate_combinations(n, r):
    return factorial(n) // factorial(r) // factorial(n-r)

Bei einigen von mir getesteten Eingaben (z. B. n = 1000 r = 500) war dies mehr als zehnmal schneller als der eine Liner reduce, der in einer anderen (derzeit am höchsten bewerteten) Antwort vorgeschlagen wurde. Auf der anderen Seite wird es von dem von @JF Sebastian bereitgestellten Snippit übertroffen.

Todd Owen
quelle
11

Ab Python 3.8sofort enthält die Standardbibliothek jetzt die math.combFunktion zur Berechnung des Binomialkoeffizienten:

math.comb (n, k)

Dies ist die Anzahl der Möglichkeiten, k Elemente aus n Elementen ohne Wiederholung auszuwählen
n! / (k! (n - k)!):

import math
math.comb(10, 5) # 252
Xavier Guihot
quelle
10

Hier ist eine andere Alternative. Dieser wurde ursprünglich in C ++ geschrieben, sodass er für eine Ganzzahl mit endlicher Genauigkeit (z. B. __int64) nach C ++ zurückportiert werden kann. Der Vorteil besteht darin, dass (1) nur ganzzahlige Operationen erforderlich sind und (2) das Aufblähen des ganzzahligen Werts durch aufeinanderfolgende Multiplikations- und Divisionspaare vermieden wird. Ich habe das Ergebnis mit Nas Banovs Pascal-Dreieck getestet, es erhält die richtige Antwort:

def choose(n,r):
  """Computes n! / (r! (n-r)!) exactly. Returns a python long int."""
  assert n >= 0
  assert 0 <= r <= n

  c = 1L
  denom = 1
  for (num,denom) in zip(xrange(n,n-r,-1), xrange(1,r+1,1)):
    c = (c * num) // denom
  return c

Begründung: Um die Anzahl der Multiplikationen und Divisionen zu minimieren, schreiben wir den Ausdruck wie folgt um

    n!      n(n-1)...(n-r+1)
--------- = ----------------
 r!(n-r)!          r!

Um einen Multiplikationsüberlauf so weit wie möglich zu vermeiden, werden wir in der folgenden STRICT-Reihenfolge von links nach rechts auswerten:

n / 1 * (n-1) / 2 * (n-2) / 3 * ... * (n-r+1) / r

Wir können zeigen, dass die in dieser Reihenfolge betriebene Ganzzahlarithmetik genau ist (dh kein Rundungsfehler).

Wirawan Purwanto
quelle
5

Bei dynamischer Programmierung beträgt die zeitliche Komplexität Θ (n * m) und die räumliche Komplexität Θ (m):

def binomial(n, k):
""" (int, int) -> int

         | c(n-1, k-1) + c(n-1, k), if 0 < k < n
c(n,k) = | 1                      , if n = k
         | 1                      , if k = 0

Precondition: n > k

>>> binomial(9, 2)
36
"""

c = [0] * (n + 1)
c[0] = 1
for i in range(1, n + 1):
    c[i] = 1
    j = i - 1
    while j > 0:
        c[j] += c[j - 1]
        j -= 1

return c[k]
pantelis300
quelle
4

Wenn Ihr Programm eine Obergrenze für n(sagen wir n <= N) hat und nCr wiederholt berechnen muss (vorzugsweise für >> NZeiten), kann die Verwendung von lru_cache zu einer enormen Leistungssteigerung führen:

from functools import lru_cache

@lru_cache(maxsize=None)
def nCr(n, r):
    return 1 if r == 0 or r == n else nCr(n - 1, r - 1) + nCr(n - 1, r)

Das Erstellen des Caches (was implizit erfolgt) nimmt O(N^2)Zeit in Anspruch. Alle nachfolgenden Anrufe an nCrwerden zurückgegeben O(1).

yzn-pku
quelle
4

Sie können 2 einfache Funktionen schreiben, die sich tatsächlich als 5-8-mal schneller herausstellen als mit scipy.special.comb . Tatsächlich müssen Sie keine zusätzlichen Pakete importieren, und die Funktion ist recht gut lesbar. Der Trick besteht darin, die zuvor berechneten Werte mithilfe der Memoisierung zu speichern und die Definition von nCr zu verwenden

# create a memoization dictionary
memo = {}
def factorial(n):
    """
    Calculate the factorial of an input using memoization
    :param n: int
    :rtype value: int
    """
    if n in [1,0]:
        return 1
    if n in memo:
        return memo[n]
    value = n*factorial(n-1)
    memo[n] = value
    return value

def ncr(n, k):
    """
    Choose k elements from a set of n elements - n must be larger than or equal to k
    :param n: int
    :param k: int
    :rtype: int
    """
    return factorial(n)/(factorial(k)*factorial(n-k))

Wenn wir mal vergleichen

from scipy.special import comb
%timeit comb(100,48)
>>> 100000 loops, best of 3: 6.78 µs per loop

%timeit ncr(100,48)
>>> 1000000 loops, best of 3: 1.39 µs per loop
PyRsquared
quelle
Heutzutage gibt es in functools einen Memoize Decorator namens lru_cache, der Ihren Code vereinfachen könnte.
wahnsinniger Igel
2

Mit Sympy ist das ziemlich einfach.

import sympy

comb = sympy.binomial(n, r)
Bobby
quelle
2

Verwenden Sie nur die mit Python vertriebene Standardbibliothek :

import itertools

def nCk(n, k):
    return len(list(itertools.combinations(range(n), k)))
MarianD
quelle
3
Ich denke nicht, dass seine zeitliche Komplexität (und Speichernutzung) akzeptabel ist.
xmcp
2

Die direkte Formel erzeugt große ganze Zahlen, wenn n größer als 20 ist.

Also noch eine Antwort:

from math import factorial

reduce(long.__mul__, range(n-r+1, n+1), 1L) // factorial(r)

kurz, genau und effizient, da dies große Python-Ganzzahlen vermeidet, indem es bei Longs bleibt.

Es ist genauer und schneller im Vergleich zu scipy.special.comb:

 >>> from scipy.special import comb
 >>> nCr = lambda n,r: reduce(long.__mul__, range(n-r+1, n+1), 1L) // factorial(r)
 >>> comb(128,20)
 1.1965669823265365e+23
 >>> nCr(128,20)
 119656698232656998274400L  # accurate, no loss
 >>> from timeit import timeit
 >>> timeit(lambda: comb(n,r))
 8.231969118118286
 >>> timeit(lambda: nCr(128, 20))
 3.885951042175293
Olivencoder
quelle
Das ist falsch! Wenn n == r, sollte das Ergebnis 1 sein. Dieser Code gibt 0 zurück.
reyammer
Genauer gesagt sollte es range(n-r+1, n+1)statt sein range(n-r,n+1).
Reyammer
1

Dies ist der @ KillerT2333-Code, der den integrierten Memoization Decorator verwendet.

from functools import lru_cache

@lru_cache()
def factorial(n):
    """
    Calculate the factorial of an input using memoization
    :param n: int
    :rtype value: int
    """
    return 1 if n in (1, 0) else n * factorial(n-1)

@lru_cache()
def ncr(n, k):
    """
    Choose k elements from a set of n elements,
    n must be greater than or equal to k.
    :param n: int
    :param k: int
    :rtype: int
    """
    return factorial(n) / (factorial(k) * factorial(n - k))

print(ncr(6, 3))
wahnsinniger Igel
quelle
1

Hier ist ein effizienter Algorithmus für Sie

for i = 1.....r

   p = p * ( n - i ) / i

print(p)

Zum Beispiel nCr (30,7) = Fakt (30) / (Fakt (7) * Fakt (23)) = (30 * 29 * 28 * 27 * 26 * 25 * 24) / (1 * 2 * 3 * 4 * 5 * 6 * 7)

Führen Sie also einfach die Schleife von 1 nach r aus, um das Ergebnis zu erhalten.

kta
quelle
0

Das ist wahrscheinlich so schnell, wie Sie es in reinem Python für relativ große Eingaben tun können:

def choose(n, k):
    if k == n: return 1
    if k > n: return 0
    d, q = max(k, n-k), min(k, n-k)
    num =  1
    for n in xrange(d+1, n+1): num *= n
    denom = 1
    for d in xrange(1, q+1): denom *= d
    return num / denom
Rabih Kodeih
quelle
0

Diese Funktion ist sehr optimiert.

def nCk(n,k):
    m=0
    if k==0:
        m=1
    if k==1:
        m=n
    if k>=2:
        num,dem,op1,op2=1,1,k,n
        while(op1>=1):
            num*=op2
            dem*=op1
            op1-=1
            op2-=1
        m=num//dem
    return m
Santiago Coca Rojas
quelle