probability.py 88 KB

12345678910111213141516171819202122232425262728293031323334353637383940414243444546474849505152535455565758596061626364656667686970717273747576777879808182838485868788899091929394959697989910010110210310410510610710810911011111211311411511611711811912012112212312412512612712812913013113213313413513613713813914014114214314414514614714814915015115215315415515615715815916016116216316416516616716816917017117217317417517617717817918018118218318418518618718818919019119219319419519619719819920020120220320420520620720820921021121221321421521621721821922022122222322422522622722822923023123223323423523623723823924024124224324424524624724824925025125225325425525625725825926026126226326426526626726826927027127227327427527627727827928028128228328428528628728828929029129229329429529629729829930030130230330430530630730830931031131231331431531631731831932032132232332432532632732832933033133233333433533633733833934034134234334434534634734834935035135235335435535635735835936036136236336436536636736836937037137237337437537637737837938038138238338438538638738838939039139239339439539639739839940040140240340440540640740840941041141241341441541641741841942042142242342442542642742842943043143243343443543643743843944044144244344444544644744844945045145245345445545645745845946046146246346446546646746846947047147247347447547647747847948048148248348448548648748848949049149249349449549649749849950050150250350450550650750850951051151251351451551651751851952052152252352452552652752852953053153253353453553653753853954054154254354454554654754854955055155255355455555655755855956056156256356456556656756856957057157257357457557657757857958058158258358458558658758858959059159259359459559659759859960060160260360460560660760860961061161261361461561661761861962062162262362462562662762862963063163263363463563663763863964064164264364464564664764864965065165265365465565665765865966066166266366466566666766866967067167267367467567667767867968068168268368468568668768868969069169269369469569669769869970070170270370470570670770870971071171271371471571671771871972072172272372472572672772872973073173273373473573673773873974074174274374474574674774874975075175275375475575675775875976076176276376476576676776876977077177277377477577677777877978078178278378478578678778878979079179279379479579679779879980080180280380480580680780880981081181281381481581681781881982082182282382482582682782882983083183283383483583683783883984084184284384484584684784884985085185285385485585685785885986086186286386486586686786886987087187287387487587687787887988088188288388488588688788888989089189289389489589689789889990090190290390490590690790890991091191291391491591691791891992092192292392492592692792892993093193293393493593693793893994094194294394494594694794894995095195295395495595695795895996096196296396496596696796896997097197297397497597697797897998098198298398498598698798898999099199299399499599699799899910001001100210031004100510061007100810091010101110121013101410151016101710181019102010211022102310241025102610271028102910301031103210331034103510361037103810391040104110421043104410451046104710481049105010511052105310541055105610571058105910601061106210631064106510661067106810691070107110721073107410751076107710781079108010811082108310841085108610871088108910901091109210931094109510961097109810991100110111021103110411051106110711081109111011111112111311141115111611171118111911201121112211231124112511261127112811291130113111321133113411351136113711381139114011411142114311441145114611471148114911501151115211531154115511561157115811591160116111621163116411651166116711681169117011711172117311741175117611771178117911801181118211831184118511861187118811891190119111921193119411951196119711981199120012011202120312041205120612071208120912101211121212131214121512161217121812191220122112221223122412251226122712281229123012311232123312341235123612371238123912401241124212431244124512461247124812491250125112521253125412551256125712581259126012611262126312641265126612671268126912701271127212731274127512761277127812791280128112821283128412851286128712881289129012911292129312941295129612971298129913001301130213031304130513061307130813091310131113121313131413151316131713181319132013211322132313241325132613271328132913301331133213331334133513361337133813391340134113421343134413451346134713481349135013511352135313541355135613571358135913601361136213631364136513661367136813691370137113721373137413751376137713781379138013811382138313841385138613871388138913901391139213931394139513961397139813991400140114021403140414051406140714081409141014111412141314141415141614171418141914201421142214231424142514261427142814291430143114321433143414351436143714381439144014411442144314441445144614471448144914501451145214531454145514561457145814591460146114621463146414651466146714681469147014711472147314741475147614771478147914801481148214831484148514861487148814891490149114921493149414951496149714981499150015011502150315041505150615071508150915101511151215131514151515161517151815191520152115221523152415251526152715281529153015311532153315341535153615371538153915401541154215431544154515461547154815491550155115521553155415551556155715581559156015611562156315641565156615671568156915701571157215731574157515761577157815791580158115821583158415851586158715881589159015911592159315941595159615971598159916001601160216031604160516061607160816091610161116121613161416151616161716181619162016211622162316241625162616271628162916301631163216331634163516361637163816391640164116421643164416451646164716481649165016511652165316541655165616571658165916601661166216631664166516661667166816691670167116721673167416751676167716781679168016811682168316841685168616871688168916901691169216931694169516961697169816991700170117021703170417051706170717081709171017111712171317141715171617171718171917201721172217231724172517261727172817291730173117321733173417351736173717381739174017411742174317441745174617471748174917501751175217531754175517561757175817591760176117621763176417651766176717681769177017711772177317741775177617771778177917801781178217831784178517861787178817891790179117921793179417951796179717981799180018011802180318041805180618071808180918101811181218131814181518161817181818191820182118221823182418251826182718281829183018311832183318341835183618371838183918401841184218431844184518461847184818491850185118521853185418551856185718581859186018611862186318641865186618671868186918701871187218731874187518761877187818791880188118821883188418851886188718881889189018911892189318941895189618971898189919001901190219031904190519061907190819091910191119121913191419151916191719181919192019211922192319241925192619271928192919301931193219331934193519361937193819391940194119421943194419451946194719481949195019511952195319541955195619571958195919601961196219631964196519661967196819691970197119721973197419751976197719781979198019811982198319841985198619871988198919901991199219931994199519961997199819992000200120022003200420052006200720082009201020112012201320142015201620172018201920202021202220232024202520262027202820292030203120322033203420352036203720382039204020412042204320442045204620472048204920502051205220532054205520562057205820592060206120622063206420652066206720682069207020712072207320742075207620772078207920802081208220832084208520862087208820892090209120922093209420952096209720982099210021012102210321042105210621072108210921102111211221132114211521162117211821192120212121222123212421252126212721282129213021312132213321342135213621372138213921402141214221432144214521462147214821492150215121522153215421552156215721582159216021612162216321642165216621672168216921702171217221732174217521762177217821792180218121822183218421852186218721882189219021912192219321942195219621972198219922002201220222032204220522062207220822092210221122122213221422152216221722182219222022212222222322242225222622272228222922302231223222332234223522362237223822392240224122422243224422452246224722482249225022512252225322542255225622572258225922602261226222632264226522662267226822692270227122722273227422752276227722782279228022812282228322842285228622872288228922902291229222932294229522962297229822992300230123022303230423052306230723082309231023112312231323142315231623172318231923202321232223232324232523262327232823292330233123322333233423352336233723382339234023412342234323442345234623472348234923502351235223532354235523562357235823592360236123622363236423652366236723682369237023712372237323742375237623772378237923802381238223832384238523862387238823892390239123922393239423952396239723982399240024012402240324042405240624072408240924102411241224132414241524162417241824192420242124222423242424252426242724282429243024312432243324342435243624372438243924402441244224432444244524462447244824492450245124522453245424552456245724582459246024612462246324642465246624672468246924702471247224732474247524762477247824792480248124822483248424852486248724882489249024912492249324942495249624972498249925002501250225032504250525062507250825092510251125122513251425152516251725182519252025212522252325242525252625272528252925302531253225332534253525362537253825392540
  1. # -*- coding: utf-8 -*-
  2. # Natural Language Toolkit: Probability and Statistics
  3. #
  4. # Copyright (C) 2001-2019 NLTK Project
  5. # Author: Edward Loper <edloper@gmail.com>
  6. # Steven Bird <stevenbird1@gmail.com> (additions)
  7. # Trevor Cohn <tacohn@cs.mu.oz.au> (additions)
  8. # Peter Ljunglöf <peter.ljunglof@heatherleaf.se> (additions)
  9. # Liang Dong <ldong@clemson.edu> (additions)
  10. # Geoffrey Sampson <sampson@cantab.net> (additions)
  11. # Ilia Kurenkov <ilia.kurenkov@gmail.com> (additions)
  12. #
  13. # URL: <http://nltk.org/>
  14. # For license information, see LICENSE.TXT
  15. """
  16. Classes for representing and processing probabilistic information.
  17. The ``FreqDist`` class is used to encode "frequency distributions",
  18. which count the number of times that each outcome of an experiment
  19. occurs.
  20. The ``ProbDistI`` class defines a standard interface for "probability
  21. distributions", which encode the probability of each outcome for an
  22. experiment. There are two types of probability distribution:
  23. - "derived probability distributions" are created from frequency
  24. distributions. They attempt to model the probability distribution
  25. that generated the frequency distribution.
  26. - "analytic probability distributions" are created directly from
  27. parameters (such as variance).
  28. The ``ConditionalFreqDist`` class and ``ConditionalProbDistI`` interface
  29. are used to encode conditional distributions. Conditional probability
  30. distributions can be derived or analytic; but currently the only
  31. implementation of the ``ConditionalProbDistI`` interface is
  32. ``ConditionalProbDist``, a derived distribution.
  33. """
  34. from __future__ import print_function, unicode_literals, division
  35. import math
  36. import random
  37. import warnings
  38. import array
  39. from collections import defaultdict, Counter
  40. from functools import reduce
  41. from abc import ABCMeta, abstractmethod
  42. from six import itervalues, text_type, add_metaclass
  43. from nltk import compat
  44. from nltk.internals import raise_unorderable_types
  45. _NINF = float('-1e300')
  46. ##//////////////////////////////////////////////////////
  47. ## Frequency Distributions
  48. ##//////////////////////////////////////////////////////
  49. @compat.python_2_unicode_compatible
  50. class FreqDist(Counter):
  51. """
  52. A frequency distribution for the outcomes of an experiment. A
  53. frequency distribution records the number of times each outcome of
  54. an experiment has occurred. For example, a frequency distribution
  55. could be used to record the frequency of each word type in a
  56. document. Formally, a frequency distribution can be defined as a
  57. function mapping from each sample to the number of times that
  58. sample occurred as an outcome.
  59. Frequency distributions are generally constructed by running a
  60. number of experiments, and incrementing the count for a sample
  61. every time it is an outcome of an experiment. For example, the
  62. following code will produce a frequency distribution that encodes
  63. how often each word occurs in a text:
  64. >>> from nltk.tokenize import word_tokenize
  65. >>> from nltk.probability import FreqDist
  66. >>> sent = 'This is an example sentence'
  67. >>> fdist = FreqDist()
  68. >>> for word in word_tokenize(sent):
  69. ... fdist[word.lower()] += 1
  70. An equivalent way to do this is with the initializer:
  71. >>> fdist = FreqDist(word.lower() for word in word_tokenize(sent))
  72. """
  73. def __init__(self, samples=None):
  74. """
  75. Construct a new frequency distribution. If ``samples`` is
  76. given, then the frequency distribution will be initialized
  77. with the count of each object in ``samples``; otherwise, it
  78. will be initialized to be empty.
  79. In particular, ``FreqDist()`` returns an empty frequency
  80. distribution; and ``FreqDist(samples)`` first creates an empty
  81. frequency distribution, and then calls ``update`` with the
  82. list ``samples``.
  83. :param samples: The samples to initialize the frequency
  84. distribution with.
  85. :type samples: Sequence
  86. """
  87. Counter.__init__(self, samples)
  88. # Cached number of samples in this FreqDist
  89. self._N = None
  90. def N(self):
  91. """
  92. Return the total number of sample outcomes that have been
  93. recorded by this FreqDist. For the number of unique
  94. sample values (or bins) with counts greater than zero, use
  95. ``FreqDist.B()``.
  96. :rtype: int
  97. """
  98. if self._N is None:
  99. # Not already cached, or cache has been invalidated
  100. self._N = sum(self.values())
  101. return self._N
  102. def __setitem__(self, key, val):
  103. """
  104. Override ``Counter.__setitem__()`` to invalidate the cached N
  105. """
  106. self._N = None
  107. super(FreqDist, self).__setitem__(key, val)
  108. def __delitem__(self, key):
  109. """
  110. Override ``Counter.__delitem__()`` to invalidate the cached N
  111. """
  112. self._N = None
  113. super(FreqDist, self).__delitem__(key)
  114. def update(self, *args, **kwargs):
  115. """
  116. Override ``Counter.update()`` to invalidate the cached N
  117. """
  118. self._N = None
  119. super(FreqDist, self).update(*args, **kwargs)
  120. def setdefault(self, key, val):
  121. """
  122. Override ``Counter.setdefault()`` to invalidate the cached N
  123. """
  124. self._N = None
  125. super(FreqDist, self).setdefault(key, val)
  126. def B(self):
  127. """
  128. Return the total number of sample values (or "bins") that
  129. have counts greater than zero. For the total
  130. number of sample outcomes recorded, use ``FreqDist.N()``.
  131. (FreqDist.B() is the same as len(FreqDist).)
  132. :rtype: int
  133. """
  134. return len(self)
  135. def hapaxes(self):
  136. """
  137. Return a list of all samples that occur once (hapax legomena)
  138. :rtype: list
  139. """
  140. return [item for item in self if self[item] == 1]
  141. def Nr(self, r, bins=None):
  142. return self.r_Nr(bins)[r]
  143. def r_Nr(self, bins=None):
  144. """
  145. Return the dictionary mapping r to Nr, the number of samples with frequency r, where Nr > 0.
  146. :type bins: int
  147. :param bins: The number of possible sample outcomes. ``bins``
  148. is used to calculate Nr(0). In particular, Nr(0) is
  149. ``bins-self.B()``. If ``bins`` is not specified, it
  150. defaults to ``self.B()`` (so Nr(0) will be 0).
  151. :rtype: int
  152. """
  153. _r_Nr = defaultdict(int)
  154. for count in self.values():
  155. _r_Nr[count] += 1
  156. # Special case for Nr[0]:
  157. _r_Nr[0] = bins - self.B() if bins is not None else 0
  158. return _r_Nr
  159. def _cumulative_frequencies(self, samples):
  160. """
  161. Return the cumulative frequencies of the specified samples.
  162. If no samples are specified, all counts are returned, starting
  163. with the largest.
  164. :param samples: the samples whose frequencies should be returned.
  165. :type samples: any
  166. :rtype: list(float)
  167. """
  168. cf = 0.0
  169. for sample in samples:
  170. cf += self[sample]
  171. yield cf
  172. # slightly odd nomenclature freq() if FreqDist does counts and ProbDist does probs,
  173. # here, freq() does probs
  174. def freq(self, sample):
  175. """
  176. Return the frequency of a given sample. The frequency of a
  177. sample is defined as the count of that sample divided by the
  178. total number of sample outcomes that have been recorded by
  179. this FreqDist. The count of a sample is defined as the
  180. number of times that sample outcome was recorded by this
  181. FreqDist. Frequencies are always real numbers in the range
  182. [0, 1].
  183. :param sample: the sample whose frequency
  184. should be returned.
  185. :type sample: any
  186. :rtype: float
  187. """
  188. n = self.N()
  189. if n == 0:
  190. return 0
  191. return self[sample] / n
  192. def max(self):
  193. """
  194. Return the sample with the greatest number of outcomes in this
  195. frequency distribution. If two or more samples have the same
  196. number of outcomes, return one of them; which sample is
  197. returned is undefined. If no outcomes have occurred in this
  198. frequency distribution, return None.
  199. :return: The sample with the maximum number of outcomes in this
  200. frequency distribution.
  201. :rtype: any or None
  202. """
  203. if len(self) == 0:
  204. raise ValueError(
  205. 'A FreqDist must have at least one sample before max is defined.'
  206. )
  207. return self.most_common(1)[0][0]
  208. def plot(self, *args, **kwargs):
  209. """
  210. Plot samples from the frequency distribution
  211. displaying the most frequent sample first. If an integer
  212. parameter is supplied, stop after this many samples have been
  213. plotted. For a cumulative plot, specify cumulative=True.
  214. (Requires Matplotlib to be installed.)
  215. :param title: The title for the graph
  216. :type title: str
  217. :param cumulative: A flag to specify whether the plot is cumulative (default = False)
  218. :type title: bool
  219. """
  220. try:
  221. import matplotlib.pyplot as plt
  222. except ImportError:
  223. raise ValueError(
  224. 'The plot function requires matplotlib to be installed.'
  225. 'See http://matplotlib.org/'
  226. )
  227. if len(args) == 0:
  228. args = [len(self)]
  229. samples = [item for item, _ in self.most_common(*args)]
  230. cumulative = _get_kwarg(kwargs, 'cumulative', False)
  231. percents = _get_kwarg(kwargs, 'percents', False)
  232. if cumulative:
  233. freqs = list(self._cumulative_frequencies(samples))
  234. ylabel = "Cumulative Counts"
  235. if percents:
  236. freqs = [f / freqs[len(freqs) - 1] * 100 for f in freqs]
  237. ylabel = "Cumulative Percents"
  238. else:
  239. freqs = [self[sample] for sample in samples]
  240. ylabel = "Counts"
  241. # percents = [f * 100 for f in freqs] only in ProbDist?
  242. ax = plt.gca()
  243. ax.grid(True, color = "silver")
  244. if "linewidth" not in kwargs:
  245. kwargs["linewidth"] = 2
  246. if "title" in kwargs:
  247. ax.set_title(kwargs["title"])
  248. del kwargs["title"]
  249. ax.plot(freqs, **kwargs)
  250. ax.set_xticks(range(len(samples)))
  251. ax.set_xticklabels([text_type(s) for s in samples], rotation=90)
  252. ax.set_xlabel("Samples")
  253. ax.set_ylabel(ylabel)
  254. plt.show()
  255. return ax
  256. def tabulate(self, *args, **kwargs):
  257. """
  258. Tabulate the given samples from the frequency distribution (cumulative),
  259. displaying the most frequent sample first. If an integer
  260. parameter is supplied, stop after this many samples have been
  261. plotted.
  262. :param samples: The samples to plot (default is all samples)
  263. :type samples: list
  264. :param cumulative: A flag to specify whether the freqs are cumulative (default = False)
  265. :type title: bool
  266. """
  267. if len(args) == 0:
  268. args = [len(self)]
  269. samples = [item for item, _ in self.most_common(*args)]
  270. cumulative = _get_kwarg(kwargs, 'cumulative', False)
  271. if cumulative:
  272. freqs = list(self._cumulative_frequencies(samples))
  273. else:
  274. freqs = [self[sample] for sample in samples]
  275. # percents = [f * 100 for f in freqs] only in ProbDist?
  276. width = max(len("%s" % s) for s in samples)
  277. width = max(width, max(len("%d" % f) for f in freqs))
  278. for i in range(len(samples)):
  279. print("%*s" % (width, samples[i]), end=' ')
  280. print()
  281. for i in range(len(samples)):
  282. print("%*d" % (width, freqs[i]), end=' ')
  283. print()
  284. def copy(self):
  285. """
  286. Create a copy of this frequency distribution.
  287. :rtype: FreqDist
  288. """
  289. return self.__class__(self)
  290. # Mathematical operatiors
  291. def __add__(self, other):
  292. """
  293. Add counts from two counters.
  294. >>> FreqDist('abbb') + FreqDist('bcc')
  295. FreqDist({'b': 4, 'c': 2, 'a': 1})
  296. """
  297. return self.__class__(super(FreqDist, self).__add__(other))
  298. def __sub__(self, other):
  299. """
  300. Subtract count, but keep only results with positive counts.
  301. >>> FreqDist('abbbc') - FreqDist('bccd')
  302. FreqDist({'b': 2, 'a': 1})
  303. """
  304. return self.__class__(super(FreqDist, self).__sub__(other))
  305. def __or__(self, other):
  306. """
  307. Union is the maximum of value in either of the input counters.
  308. >>> FreqDist('abbb') | FreqDist('bcc')
  309. FreqDist({'b': 3, 'c': 2, 'a': 1})
  310. """
  311. return self.__class__(super(FreqDist, self).__or__(other))
  312. def __and__(self, other):
  313. """
  314. Intersection is the minimum of corresponding counts.
  315. >>> FreqDist('abbb') & FreqDist('bcc')
  316. FreqDist({'b': 1})
  317. """
  318. return self.__class__(super(FreqDist, self).__and__(other))
  319. def __le__(self, other):
  320. if not isinstance(other, FreqDist):
  321. raise_unorderable_types("<=", self, other)
  322. return set(self).issubset(other) and all(
  323. self[key] <= other[key] for key in self
  324. )
  325. # @total_ordering doesn't work here, since the class inherits from a builtin class
  326. __ge__ = lambda self, other: not self <= other or self == other
  327. __lt__ = lambda self, other: self <= other and not self == other
  328. __gt__ = lambda self, other: not self <= other
  329. def __repr__(self):
  330. """
  331. Return a string representation of this FreqDist.
  332. :rtype: string
  333. """
  334. return self.pformat()
  335. def pprint(self, maxlen=10, stream=None):
  336. """
  337. Print a string representation of this FreqDist to 'stream'
  338. :param maxlen: The maximum number of items to print
  339. :type maxlen: int
  340. :param stream: The stream to print to. stdout by default
  341. """
  342. print(self.pformat(maxlen=maxlen), file=stream)
  343. def pformat(self, maxlen=10):
  344. """
  345. Return a string representation of this FreqDist.
  346. :param maxlen: The maximum number of items to display
  347. :type maxlen: int
  348. :rtype: string
  349. """
  350. items = ['{0!r}: {1!r}'.format(*item) for item in self.most_common(maxlen)]
  351. if len(self) > maxlen:
  352. items.append('...')
  353. return 'FreqDist({{{0}}})'.format(', '.join(items))
  354. def __str__(self):
  355. """
  356. Return a string representation of this FreqDist.
  357. :rtype: string
  358. """
  359. return '<FreqDist with %d samples and %d outcomes>' % (len(self), self.N())
  360. ##//////////////////////////////////////////////////////
  361. ## Probability Distributions
  362. ##//////////////////////////////////////////////////////
  363. @add_metaclass(ABCMeta)
  364. class ProbDistI(object):
  365. """
  366. A probability distribution for the outcomes of an experiment. A
  367. probability distribution specifies how likely it is that an
  368. experiment will have any given outcome. For example, a
  369. probability distribution could be used to predict the probability
  370. that a token in a document will have a given type. Formally, a
  371. probability distribution can be defined as a function mapping from
  372. samples to nonnegative real numbers, such that the sum of every
  373. number in the function's range is 1.0. A ``ProbDist`` is often
  374. used to model the probability distribution of the experiment used
  375. to generate a frequency distribution.
  376. """
  377. SUM_TO_ONE = True
  378. """True if the probabilities of the samples in this probability
  379. distribution will always sum to one."""
  380. @abstractmethod
  381. def __init__(self):
  382. """
  383. Classes inheriting from ProbDistI should implement __init__.
  384. """
  385. @abstractmethod
  386. def prob(self, sample):
  387. """
  388. Return the probability for a given sample. Probabilities
  389. are always real numbers in the range [0, 1].
  390. :param sample: The sample whose probability
  391. should be returned.
  392. :type sample: any
  393. :rtype: float
  394. """
  395. def logprob(self, sample):
  396. """
  397. Return the base 2 logarithm of the probability for a given sample.
  398. :param sample: The sample whose probability
  399. should be returned.
  400. :type sample: any
  401. :rtype: float
  402. """
  403. # Default definition, in terms of prob()
  404. p = self.prob(sample)
  405. return math.log(p, 2) if p != 0 else _NINF
  406. @abstractmethod
  407. def max(self):
  408. """
  409. Return the sample with the greatest probability. If two or
  410. more samples have the same probability, return one of them;
  411. which sample is returned is undefined.
  412. :rtype: any
  413. """
  414. @abstractmethod
  415. def samples(self):
  416. """
  417. Return a list of all samples that have nonzero probabilities.
  418. Use ``prob`` to find the probability of each sample.
  419. :rtype: list
  420. """
  421. # cf self.SUM_TO_ONE
  422. def discount(self):
  423. """
  424. Return the ratio by which counts are discounted on average: c*/c
  425. :rtype: float
  426. """
  427. return 0.0
  428. # Subclasses should define more efficient implementations of this,
  429. # where possible.
  430. def generate(self):
  431. """
  432. Return a randomly selected sample from this probability distribution.
  433. The probability of returning each sample ``samp`` is equal to
  434. ``self.prob(samp)``.
  435. """
  436. p = random.random()
  437. p_init = p
  438. for sample in self.samples():
  439. p -= self.prob(sample)
  440. if p <= 0:
  441. return sample
  442. # allow for some rounding error:
  443. if p < 0.0001:
  444. return sample
  445. # we *should* never get here
  446. if self.SUM_TO_ONE:
  447. warnings.warn(
  448. "Probability distribution %r sums to %r; generate()"
  449. " is returning an arbitrary sample." % (self, p_init - p)
  450. )
  451. return random.choice(list(self.samples()))
  452. @compat.python_2_unicode_compatible
  453. class UniformProbDist(ProbDistI):
  454. """
  455. A probability distribution that assigns equal probability to each
  456. sample in a given set; and a zero probability to all other
  457. samples.
  458. """
  459. def __init__(self, samples):
  460. """
  461. Construct a new uniform probability distribution, that assigns
  462. equal probability to each sample in ``samples``.
  463. :param samples: The samples that should be given uniform
  464. probability.
  465. :type samples: list
  466. :raise ValueError: If ``samples`` is empty.
  467. """
  468. if len(samples) == 0:
  469. raise ValueError(
  470. 'A Uniform probability distribution must ' + 'have at least one sample.'
  471. )
  472. self._sampleset = set(samples)
  473. self._prob = 1.0 / len(self._sampleset)
  474. self._samples = list(self._sampleset)
  475. def prob(self, sample):
  476. return self._prob if sample in self._sampleset else 0
  477. def max(self):
  478. return self._samples[0]
  479. def samples(self):
  480. return self._samples
  481. def __repr__(self):
  482. return '<UniformProbDist with %d samples>' % len(self._sampleset)
  483. @compat.python_2_unicode_compatible
  484. class RandomProbDist(ProbDistI):
  485. """
  486. Generates a random probability distribution whereby each sample
  487. will be between 0 and 1 with equal probability (uniform random distribution.
  488. Also called a continuous uniform distribution).
  489. """
  490. def __init__(self, samples):
  491. if len(samples) == 0:
  492. raise ValueError(
  493. 'A probability distribution must ' + 'have at least one sample.'
  494. )
  495. self._probs = self.unirand(samples)
  496. self._samples = list(self._probs.keys())
  497. @classmethod
  498. def unirand(cls, samples):
  499. """
  500. The key function that creates a randomized initial distribution
  501. that still sums to 1. Set as a dictionary of prob values so that
  502. it can still be passed to MutableProbDist and called with identical
  503. syntax to UniformProbDist
  504. """
  505. samples = set(samples)
  506. randrow = [random.random() for i in range(len(samples))]
  507. total = sum(randrow)
  508. for i, x in enumerate(randrow):
  509. randrow[i] = x / total
  510. total = sum(randrow)
  511. if total != 1:
  512. # this difference, if present, is so small (near NINF) that it
  513. # can be subtracted from any element without risking probs not (0 1)
  514. randrow[-1] -= total - 1
  515. return dict((s, randrow[i]) for i, s in enumerate(samples))
  516. def max(self):
  517. if not hasattr(self, '_max'):
  518. self._max = max((p, v) for (v, p) in self._probs.items())[1]
  519. return self._max
  520. def prob(self, sample):
  521. return self._probs.get(sample, 0)
  522. def samples(self):
  523. return self._samples
  524. def __repr__(self):
  525. return '<RandomUniformProbDist with %d samples>' % len(self._probs)
  526. @compat.python_2_unicode_compatible
  527. class DictionaryProbDist(ProbDistI):
  528. """
  529. A probability distribution whose probabilities are directly
  530. specified by a given dictionary. The given dictionary maps
  531. samples to probabilities.
  532. """
  533. def __init__(self, prob_dict=None, log=False, normalize=False):
  534. """
  535. Construct a new probability distribution from the given
  536. dictionary, which maps values to probabilities (or to log
  537. probabilities, if ``log`` is true). If ``normalize`` is
  538. true, then the probability values are scaled by a constant
  539. factor such that they sum to 1.
  540. If called without arguments, the resulting probability
  541. distribution assigns zero probability to all values.
  542. """
  543. self._prob_dict = prob_dict.copy() if prob_dict is not None else {}
  544. self._log = log
  545. # Normalize the distribution, if requested.
  546. if normalize:
  547. if len(prob_dict) == 0:
  548. raise ValueError(
  549. 'A DictionaryProbDist must have at least one sample '
  550. + 'before it can be normalized.'
  551. )
  552. if log:
  553. value_sum = sum_logs(list(self._prob_dict.values()))
  554. if value_sum <= _NINF:
  555. logp = math.log(1.0 / len(prob_dict), 2)
  556. for x in prob_dict:
  557. self._prob_dict[x] = logp
  558. else:
  559. for (x, p) in self._prob_dict.items():
  560. self._prob_dict[x] -= value_sum
  561. else:
  562. value_sum = sum(self._prob_dict.values())
  563. if value_sum == 0:
  564. p = 1.0 / len(prob_dict)
  565. for x in prob_dict:
  566. self._prob_dict[x] = p
  567. else:
  568. norm_factor = 1.0 / value_sum
  569. for (x, p) in self._prob_dict.items():
  570. self._prob_dict[x] *= norm_factor
  571. def prob(self, sample):
  572. if self._log:
  573. return 2 ** (self._prob_dict[sample]) if sample in self._prob_dict else 0
  574. else:
  575. return self._prob_dict.get(sample, 0)
  576. def logprob(self, sample):
  577. if self._log:
  578. return self._prob_dict.get(sample, _NINF)
  579. else:
  580. if sample not in self._prob_dict:
  581. return _NINF
  582. elif self._prob_dict[sample] == 0:
  583. return _NINF
  584. else:
  585. return math.log(self._prob_dict[sample], 2)
  586. def max(self):
  587. if not hasattr(self, '_max'):
  588. self._max = max((p, v) for (v, p) in self._prob_dict.items())[1]
  589. return self._max
  590. def samples(self):
  591. return self._prob_dict.keys()
  592. def __repr__(self):
  593. return '<ProbDist with %d samples>' % len(self._prob_dict)
  594. @compat.python_2_unicode_compatible
  595. class MLEProbDist(ProbDistI):
  596. """
  597. The maximum likelihood estimate for the probability distribution
  598. of the experiment used to generate a frequency distribution. The
  599. "maximum likelihood estimate" approximates the probability of
  600. each sample as the frequency of that sample in the frequency
  601. distribution.
  602. """
  603. def __init__(self, freqdist, bins=None):
  604. """
  605. Use the maximum likelihood estimate to create a probability
  606. distribution for the experiment used to generate ``freqdist``.
  607. :type freqdist: FreqDist
  608. :param freqdist: The frequency distribution that the
  609. probability estimates should be based on.
  610. """
  611. self._freqdist = freqdist
  612. def freqdist(self):
  613. """
  614. Return the frequency distribution that this probability
  615. distribution is based on.
  616. :rtype: FreqDist
  617. """
  618. return self._freqdist
  619. def prob(self, sample):
  620. return self._freqdist.freq(sample)
  621. def max(self):
  622. return self._freqdist.max()
  623. def samples(self):
  624. return self._freqdist.keys()
  625. def __repr__(self):
  626. """
  627. :rtype: str
  628. :return: A string representation of this ``ProbDist``.
  629. """
  630. return '<MLEProbDist based on %d samples>' % self._freqdist.N()
  631. @compat.python_2_unicode_compatible
  632. class LidstoneProbDist(ProbDistI):
  633. """
  634. The Lidstone estimate for the probability distribution of the
  635. experiment used to generate a frequency distribution. The
  636. "Lidstone estimate" is parameterized by a real number *gamma*,
  637. which typically ranges from 0 to 1. The Lidstone estimate
  638. approximates the probability of a sample with count *c* from an
  639. experiment with *N* outcomes and *B* bins as
  640. ``c+gamma)/(N+B*gamma)``. This is equivalent to adding
  641. *gamma* to the count for each bin, and taking the maximum
  642. likelihood estimate of the resulting frequency distribution.
  643. """
  644. SUM_TO_ONE = False
  645. def __init__(self, freqdist, gamma, bins=None):
  646. """
  647. Use the Lidstone estimate to create a probability distribution
  648. for the experiment used to generate ``freqdist``.
  649. :type freqdist: FreqDist
  650. :param freqdist: The frequency distribution that the
  651. probability estimates should be based on.
  652. :type gamma: float
  653. :param gamma: A real number used to parameterize the
  654. estimate. The Lidstone estimate is equivalent to adding
  655. *gamma* to the count for each bin, and taking the
  656. maximum likelihood estimate of the resulting frequency
  657. distribution.
  658. :type bins: int
  659. :param bins: The number of sample values that can be generated
  660. by the experiment that is described by the probability
  661. distribution. This value must be correctly set for the
  662. probabilities of the sample values to sum to one. If
  663. ``bins`` is not specified, it defaults to ``freqdist.B()``.
  664. """
  665. if (bins == 0) or (bins is None and freqdist.N() == 0):
  666. name = self.__class__.__name__[:-8]
  667. raise ValueError(
  668. 'A %s probability distribution ' % name + 'must have at least one bin.'
  669. )
  670. if (bins is not None) and (bins < freqdist.B()):
  671. name = self.__class__.__name__[:-8]
  672. raise ValueError(
  673. '\nThe number of bins in a %s distribution ' % name
  674. + '(%d) must be greater than or equal to\n' % bins
  675. + 'the number of bins in the FreqDist used '
  676. + 'to create it (%d).' % freqdist.B()
  677. )
  678. self._freqdist = freqdist
  679. self._gamma = float(gamma)
  680. self._N = self._freqdist.N()
  681. if bins is None:
  682. bins = freqdist.B()
  683. self._bins = bins
  684. self._divisor = self._N + bins * gamma
  685. if self._divisor == 0.0:
  686. # In extreme cases we force the probability to be 0,
  687. # which it will be, since the count will be 0:
  688. self._gamma = 0
  689. self._divisor = 1
  690. def freqdist(self):
  691. """
  692. Return the frequency distribution that this probability
  693. distribution is based on.
  694. :rtype: FreqDist
  695. """
  696. return self._freqdist
  697. def prob(self, sample):
  698. c = self._freqdist[sample]
  699. return (c + self._gamma) / self._divisor
  700. def max(self):
  701. # For Lidstone distributions, probability is monotonic with
  702. # frequency, so the most probable sample is the one that
  703. # occurs most frequently.
  704. return self._freqdist.max()
  705. def samples(self):
  706. return self._freqdist.keys()
  707. def discount(self):
  708. gb = self._gamma * self._bins
  709. return gb / (self._N + gb)
  710. def __repr__(self):
  711. """
  712. Return a string representation of this ``ProbDist``.
  713. :rtype: str
  714. """
  715. return '<LidstoneProbDist based on %d samples>' % self._freqdist.N()
  716. @compat.python_2_unicode_compatible
  717. class LaplaceProbDist(LidstoneProbDist):
  718. """
  719. The Laplace estimate for the probability distribution of the
  720. experiment used to generate a frequency distribution. The
  721. "Laplace estimate" approximates the probability of a sample with
  722. count *c* from an experiment with *N* outcomes and *B* bins as
  723. *(c+1)/(N+B)*. This is equivalent to adding one to the count for
  724. each bin, and taking the maximum likelihood estimate of the
  725. resulting frequency distribution.
  726. """
  727. def __init__(self, freqdist, bins=None):
  728. """
  729. Use the Laplace estimate to create a probability distribution
  730. for the experiment used to generate ``freqdist``.
  731. :type freqdist: FreqDist
  732. :param freqdist: The frequency distribution that the
  733. probability estimates should be based on.
  734. :type bins: int
  735. :param bins: The number of sample values that can be generated
  736. by the experiment that is described by the probability
  737. distribution. This value must be correctly set for the
  738. probabilities of the sample values to sum to one. If
  739. ``bins`` is not specified, it defaults to ``freqdist.B()``.
  740. """
  741. LidstoneProbDist.__init__(self, freqdist, 1, bins)
  742. def __repr__(self):
  743. """
  744. :rtype: str
  745. :return: A string representation of this ``ProbDist``.
  746. """
  747. return '<LaplaceProbDist based on %d samples>' % self._freqdist.N()
  748. @compat.python_2_unicode_compatible
  749. class ELEProbDist(LidstoneProbDist):
  750. """
  751. The expected likelihood estimate for the probability distribution
  752. of the experiment used to generate a frequency distribution. The
  753. "expected likelihood estimate" approximates the probability of a
  754. sample with count *c* from an experiment with *N* outcomes and
  755. *B* bins as *(c+0.5)/(N+B/2)*. This is equivalent to adding 0.5
  756. to the count for each bin, and taking the maximum likelihood
  757. estimate of the resulting frequency distribution.
  758. """
  759. def __init__(self, freqdist, bins=None):
  760. """
  761. Use the expected likelihood estimate to create a probability
  762. distribution for the experiment used to generate ``freqdist``.
  763. :type freqdist: FreqDist
  764. :param freqdist: The frequency distribution that the
  765. probability estimates should be based on.
  766. :type bins: int
  767. :param bins: The number of sample values that can be generated
  768. by the experiment that is described by the probability
  769. distribution. This value must be correctly set for the
  770. probabilities of the sample values to sum to one. If
  771. ``bins`` is not specified, it defaults to ``freqdist.B()``.
  772. """
  773. LidstoneProbDist.__init__(self, freqdist, 0.5, bins)
  774. def __repr__(self):
  775. """
  776. Return a string representation of this ``ProbDist``.
  777. :rtype: str
  778. """
  779. return '<ELEProbDist based on %d samples>' % self._freqdist.N()
  780. @compat.python_2_unicode_compatible
  781. class HeldoutProbDist(ProbDistI):
  782. """
  783. The heldout estimate for the probability distribution of the
  784. experiment used to generate two frequency distributions. These
  785. two frequency distributions are called the "heldout frequency
  786. distribution" and the "base frequency distribution." The
  787. "heldout estimate" uses uses the "heldout frequency
  788. distribution" to predict the probability of each sample, given its
  789. frequency in the "base frequency distribution".
  790. In particular, the heldout estimate approximates the probability
  791. for a sample that occurs *r* times in the base distribution as
  792. the average frequency in the heldout distribution of all samples
  793. that occur *r* times in the base distribution.
  794. This average frequency is *Tr[r]/(Nr[r].N)*, where:
  795. - *Tr[r]* is the total count in the heldout distribution for
  796. all samples that occur *r* times in the base distribution.
  797. - *Nr[r]* is the number of samples that occur *r* times in
  798. the base distribution.
  799. - *N* is the number of outcomes recorded by the heldout
  800. frequency distribution.
  801. In order to increase the efficiency of the ``prob`` member
  802. function, *Tr[r]/(Nr[r].N)* is precomputed for each value of *r*
  803. when the ``HeldoutProbDist`` is created.
  804. :type _estimate: list(float)
  805. :ivar _estimate: A list mapping from *r*, the number of
  806. times that a sample occurs in the base distribution, to the
  807. probability estimate for that sample. ``_estimate[r]`` is
  808. calculated by finding the average frequency in the heldout
  809. distribution of all samples that occur *r* times in the base
  810. distribution. In particular, ``_estimate[r]`` =
  811. *Tr[r]/(Nr[r].N)*.
  812. :type _max_r: int
  813. :ivar _max_r: The maximum number of times that any sample occurs
  814. in the base distribution. ``_max_r`` is used to decide how
  815. large ``_estimate`` must be.
  816. """
  817. SUM_TO_ONE = False
  818. def __init__(self, base_fdist, heldout_fdist, bins=None):
  819. """
  820. Use the heldout estimate to create a probability distribution
  821. for the experiment used to generate ``base_fdist`` and
  822. ``heldout_fdist``.
  823. :type base_fdist: FreqDist
  824. :param base_fdist: The base frequency distribution.
  825. :type heldout_fdist: FreqDist
  826. :param heldout_fdist: The heldout frequency distribution.
  827. :type bins: int
  828. :param bins: The number of sample values that can be generated
  829. by the experiment that is described by the probability
  830. distribution. This value must be correctly set for the
  831. probabilities of the sample values to sum to one. If
  832. ``bins`` is not specified, it defaults to ``freqdist.B()``.
  833. """
  834. self._base_fdist = base_fdist
  835. self._heldout_fdist = heldout_fdist
  836. # The max number of times any sample occurs in base_fdist.
  837. self._max_r = base_fdist[base_fdist.max()]
  838. # Calculate Tr, Nr, and N.
  839. Tr = self._calculate_Tr()
  840. r_Nr = base_fdist.r_Nr(bins)
  841. Nr = [r_Nr[r] for r in range(self._max_r + 1)]
  842. N = heldout_fdist.N()
  843. # Use Tr, Nr, and N to compute the probability estimate for
  844. # each value of r.
  845. self._estimate = self._calculate_estimate(Tr, Nr, N)
  846. def _calculate_Tr(self):
  847. """
  848. Return the list *Tr*, where *Tr[r]* is the total count in
  849. ``heldout_fdist`` for all samples that occur *r*
  850. times in ``base_fdist``.
  851. :rtype: list(float)
  852. """
  853. Tr = [0.0] * (self._max_r + 1)
  854. for sample in self._heldout_fdist:
  855. r = self._base_fdist[sample]
  856. Tr[r] += self._heldout_fdist[sample]
  857. return Tr
  858. def _calculate_estimate(self, Tr, Nr, N):
  859. """
  860. Return the list *estimate*, where *estimate[r]* is the probability
  861. estimate for any sample that occurs *r* times in the base frequency
  862. distribution. In particular, *estimate[r]* is *Tr[r]/(N[r].N)*.
  863. In the special case that *N[r]=0*, *estimate[r]* will never be used;
  864. so we define *estimate[r]=None* for those cases.
  865. :rtype: list(float)
  866. :type Tr: list(float)
  867. :param Tr: the list *Tr*, where *Tr[r]* is the total count in
  868. the heldout distribution for all samples that occur *r*
  869. times in base distribution.
  870. :type Nr: list(float)
  871. :param Nr: The list *Nr*, where *Nr[r]* is the number of
  872. samples that occur *r* times in the base distribution.
  873. :type N: int
  874. :param N: The total number of outcomes recorded by the heldout
  875. frequency distribution.
  876. """
  877. estimate = []
  878. for r in range(self._max_r + 1):
  879. if Nr[r] == 0:
  880. estimate.append(None)
  881. else:
  882. estimate.append(Tr[r] / (Nr[r] * N))
  883. return estimate
  884. def base_fdist(self):
  885. """
  886. Return the base frequency distribution that this probability
  887. distribution is based on.
  888. :rtype: FreqDist
  889. """
  890. return self._base_fdist
  891. def heldout_fdist(self):
  892. """
  893. Return the heldout frequency distribution that this
  894. probability distribution is based on.
  895. :rtype: FreqDist
  896. """
  897. return self._heldout_fdist
  898. def samples(self):
  899. return self._base_fdist.keys()
  900. def prob(self, sample):
  901. # Use our precomputed probability estimate.
  902. r = self._base_fdist[sample]
  903. return self._estimate[r]
  904. def max(self):
  905. # Note: the Heldout estimation is *not* necessarily monotonic;
  906. # so this implementation is currently broken. However, it
  907. # should give the right answer *most* of the time. :)
  908. return self._base_fdist.max()
  909. def discount(self):
  910. raise NotImplementedError()
  911. def __repr__(self):
  912. """
  913. :rtype: str
  914. :return: A string representation of this ``ProbDist``.
  915. """
  916. s = '<HeldoutProbDist: %d base samples; %d heldout samples>'
  917. return s % (self._base_fdist.N(), self._heldout_fdist.N())
  918. @compat.python_2_unicode_compatible
  919. class CrossValidationProbDist(ProbDistI):
  920. """
  921. The cross-validation estimate for the probability distribution of
  922. the experiment used to generate a set of frequency distribution.
  923. The "cross-validation estimate" for the probability of a sample
  924. is found by averaging the held-out estimates for the sample in
  925. each pair of frequency distributions.
  926. """
  927. SUM_TO_ONE = False
  928. def __init__(self, freqdists, bins):
  929. """
  930. Use the cross-validation estimate to create a probability
  931. distribution for the experiment used to generate
  932. ``freqdists``.
  933. :type freqdists: list(FreqDist)
  934. :param freqdists: A list of the frequency distributions
  935. generated by the experiment.
  936. :type bins: int
  937. :param bins: The number of sample values that can be generated
  938. by the experiment that is described by the probability
  939. distribution. This value must be correctly set for the
  940. probabilities of the sample values to sum to one. If
  941. ``bins`` is not specified, it defaults to ``freqdist.B()``.
  942. """
  943. self._freqdists = freqdists
  944. # Create a heldout probability distribution for each pair of
  945. # frequency distributions in freqdists.
  946. self._heldout_probdists = []
  947. for fdist1 in freqdists:
  948. for fdist2 in freqdists:
  949. if fdist1 is not fdist2:
  950. probdist = HeldoutProbDist(fdist1, fdist2, bins)
  951. self._heldout_probdists.append(probdist)
  952. def freqdists(self):
  953. """
  954. Return the list of frequency distributions that this ``ProbDist`` is based on.
  955. :rtype: list(FreqDist)
  956. """
  957. return self._freqdists
  958. def samples(self):
  959. # [xx] nb: this is not too efficient
  960. return set(sum([list(fd) for fd in self._freqdists], []))
  961. def prob(self, sample):
  962. # Find the average probability estimate returned by each
  963. # heldout distribution.
  964. prob = 0.0
  965. for heldout_probdist in self._heldout_probdists:
  966. prob += heldout_probdist.prob(sample)
  967. return prob / len(self._heldout_probdists)
  968. def discount(self):
  969. raise NotImplementedError()
  970. def __repr__(self):
  971. """
  972. Return a string representation of this ``ProbDist``.
  973. :rtype: str
  974. """
  975. return '<CrossValidationProbDist: %d-way>' % len(self._freqdists)
  976. @compat.python_2_unicode_compatible
  977. class WittenBellProbDist(ProbDistI):
  978. """
  979. The Witten-Bell estimate of a probability distribution. This distribution
  980. allocates uniform probability mass to as yet unseen events by using the
  981. number of events that have only been seen once. The probability mass
  982. reserved for unseen events is equal to *T / (N + T)*
  983. where *T* is the number of observed event types and *N* is the total
  984. number of observed events. This equates to the maximum likelihood estimate
  985. of a new type event occurring. The remaining probability mass is discounted
  986. such that all probability estimates sum to one, yielding:
  987. - *p = T / Z (N + T)*, if count = 0
  988. - *p = c / (N + T)*, otherwise
  989. """
  990. def __init__(self, freqdist, bins=None):
  991. """
  992. Creates a distribution of Witten-Bell probability estimates. This
  993. distribution allocates uniform probability mass to as yet unseen
  994. events by using the number of events that have only been seen once. The
  995. probability mass reserved for unseen events is equal to *T / (N + T)*
  996. where *T* is the number of observed event types and *N* is the total
  997. number of observed events. This equates to the maximum likelihood
  998. estimate of a new type event occurring. The remaining probability mass
  999. is discounted such that all probability estimates sum to one,
  1000. yielding:
  1001. - *p = T / Z (N + T)*, if count = 0
  1002. - *p = c / (N + T)*, otherwise
  1003. The parameters *T* and *N* are taken from the ``freqdist`` parameter
  1004. (the ``B()`` and ``N()`` values). The normalizing factor *Z* is
  1005. calculated using these values along with the ``bins`` parameter.
  1006. :param freqdist: The frequency counts upon which to base the
  1007. estimation.
  1008. :type freqdist: FreqDist
  1009. :param bins: The number of possible event types. This must be at least
  1010. as large as the number of bins in the ``freqdist``. If None, then
  1011. it's assumed to be equal to that of the ``freqdist``
  1012. :type bins: int
  1013. """
  1014. assert bins is None or bins >= freqdist.B(), (
  1015. 'bins parameter must not be less than %d=freqdist.B()' % freqdist.B()
  1016. )
  1017. if bins is None:
  1018. bins = freqdist.B()
  1019. self._freqdist = freqdist
  1020. self._T = self._freqdist.B()
  1021. self._Z = bins - self._freqdist.B()
  1022. self._N = self._freqdist.N()
  1023. # self._P0 is P(0), precalculated for efficiency:
  1024. if self._N == 0:
  1025. # if freqdist is empty, we approximate P(0) by a UniformProbDist:
  1026. self._P0 = 1.0 / self._Z
  1027. else:
  1028. self._P0 = self._T / (self._Z * (self._N + self._T))
  1029. def prob(self, sample):
  1030. # inherit docs from ProbDistI
  1031. c = self._freqdist[sample]
  1032. return c / (self._N + self._T) if c != 0 else self._P0
  1033. def max(self):
  1034. return self._freqdist.max()
  1035. def samples(self):
  1036. return self._freqdist.keys()
  1037. def freqdist(self):
  1038. return self._freqdist
  1039. def discount(self):
  1040. raise NotImplementedError()
  1041. def __repr__(self):
  1042. """
  1043. Return a string representation of this ``ProbDist``.
  1044. :rtype: str
  1045. """
  1046. return '<WittenBellProbDist based on %d samples>' % self._freqdist.N()
  1047. ##//////////////////////////////////////////////////////
  1048. ## Good-Turing Probability Distributions
  1049. ##//////////////////////////////////////////////////////
  1050. # Good-Turing frequency estimation was contributed by Alan Turing and
  1051. # his statistical assistant I.J. Good, during their collaboration in
  1052. # the WWII. It is a statistical technique for predicting the
  1053. # probability of occurrence of objects belonging to an unknown number
  1054. # of species, given past observations of such objects and their
  1055. # species. (In drawing balls from an urn, the 'objects' would be balls
  1056. # and the 'species' would be the distinct colors of the balls (finite
  1057. # but unknown in number).
  1058. #
  1059. # Good-Turing method calculates the probability mass to assign to
  1060. # events with zero or low counts based on the number of events with
  1061. # higher counts. It does so by using the adjusted count *c\**:
  1062. #
  1063. # - *c\* = (c + 1) N(c + 1) / N(c)* for c >= 1
  1064. # - *things with frequency zero in training* = N(1) for c == 0
  1065. #
  1066. # where *c* is the original count, *N(i)* is the number of event types
  1067. # observed with count *i*. We can think the count of unseen as the count
  1068. # of frequency one (see Jurafsky & Martin 2nd Edition, p101).
  1069. #
  1070. # This method is problematic because the situation ``N(c+1) == 0``
  1071. # is quite common in the original Good-Turing estimation; smoothing or
  1072. # interpolation of *N(i)* values is essential in practice.
  1073. #
  1074. # Bill Gale and Geoffrey Sampson present a simple and effective approach,
  1075. # Simple Good-Turing. As a smoothing curve they simply use a power curve:
  1076. #
  1077. # Nr = a*r^b (with b < -1 to give the appropriate hyperbolic
  1078. # relationship)
  1079. #
  1080. # They estimate a and b by simple linear regression technique on the
  1081. # logarithmic form of the equation:
  1082. #
  1083. # log Nr = a + b*log(r)
  1084. #
  1085. # However, they suggest that such a simple curve is probably only
  1086. # appropriate for high values of r. For low values of r, they use the
  1087. # measured Nr directly. (see M&S, p.213)
  1088. #
  1089. # Gale and Sampson propose to use r while the difference between r and
  1090. # r* is 1.96 greater than the standard deviation, and switch to r* if
  1091. # it is less or equal:
  1092. #
  1093. # |r - r*| > 1.96 * sqrt((r + 1)^2 (Nr+1 / Nr^2) (1 + Nr+1 / Nr))
  1094. #
  1095. # The 1.96 coefficient correspond to a 0.05 significance criterion,
  1096. # some implementations can use a coefficient of 1.65 for a 0.1
  1097. # significance criterion.
  1098. #
  1099. ##//////////////////////////////////////////////////////
  1100. ## Simple Good-Turing Probablity Distributions
  1101. ##//////////////////////////////////////////////////////
  1102. @compat.python_2_unicode_compatible
  1103. class SimpleGoodTuringProbDist(ProbDistI):
  1104. """
  1105. SimpleGoodTuring ProbDist approximates from frequency to frequency of
  1106. frequency into a linear line under log space by linear regression.
  1107. Details of Simple Good-Turing algorithm can be found in:
  1108. - Good Turing smoothing without tears" (Gale & Sampson 1995),
  1109. Journal of Quantitative Linguistics, vol. 2 pp. 217-237.
  1110. - "Speech and Language Processing (Jurafsky & Martin),
  1111. 2nd Edition, Chapter 4.5 p103 (log(Nc) = a + b*log(c))
  1112. - http://www.grsampson.net/RGoodTur.html
  1113. Given a set of pair (xi, yi), where the xi denotes the frequency and
  1114. yi denotes the frequency of frequency, we want to minimize their
  1115. square variation. E(x) and E(y) represent the mean of xi and yi.
  1116. - slope: b = sigma ((xi-E(x)(yi-E(y))) / sigma ((xi-E(x))(xi-E(x)))
  1117. - intercept: a = E(y) - b.E(x)
  1118. """
  1119. SUM_TO_ONE = False
  1120. def __init__(self, freqdist, bins=None):
  1121. """
  1122. :param freqdist: The frequency counts upon which to base the
  1123. estimation.
  1124. :type freqdist: FreqDist
  1125. :param bins: The number of possible event types. This must be
  1126. larger than the number of bins in the ``freqdist``. If None,
  1127. then it's assumed to be equal to ``freqdist``.B() + 1
  1128. :type bins: int
  1129. """
  1130. assert (
  1131. bins is None or bins > freqdist.B()
  1132. ), 'bins parameter must not be less than %d=freqdist.B()+1' % (freqdist.B() + 1)
  1133. if bins is None:
  1134. bins = freqdist.B() + 1
  1135. self._freqdist = freqdist
  1136. self._bins = bins
  1137. r, nr = self._r_Nr()
  1138. self.find_best_fit(r, nr)
  1139. self._switch(r, nr)
  1140. self._renormalize(r, nr)
  1141. def _r_Nr_non_zero(self):
  1142. r_Nr = self._freqdist.r_Nr()
  1143. del r_Nr[0]
  1144. return r_Nr
  1145. def _r_Nr(self):
  1146. """
  1147. Split the frequency distribution in two list (r, Nr), where Nr(r) > 0
  1148. """
  1149. nonzero = self._r_Nr_non_zero()
  1150. if not nonzero:
  1151. return [], []
  1152. return zip(*sorted(nonzero.items()))
  1153. def find_best_fit(self, r, nr):
  1154. """
  1155. Use simple linear regression to tune parameters self._slope and
  1156. self._intercept in the log-log space based on count and Nr(count)
  1157. (Work in log space to avoid floating point underflow.)
  1158. """
  1159. # For higher sample frequencies the data points becomes horizontal
  1160. # along line Nr=1. To create a more evident linear model in log-log
  1161. # space, we average positive Nr values with the surrounding zero
  1162. # values. (Church and Gale, 1991)
  1163. if not r or not nr:
  1164. # Empty r or nr?
  1165. return
  1166. zr = []
  1167. for j in range(len(r)):
  1168. i = r[j - 1] if j > 0 else 0
  1169. k = 2 * r[j] - i if j == len(r) - 1 else r[j + 1]
  1170. zr_ = 2.0 * nr[j] / (k - i)
  1171. zr.append(zr_)
  1172. log_r = [math.log(i) for i in r]
  1173. log_zr = [math.log(i) for i in zr]
  1174. xy_cov = x_var = 0.0
  1175. x_mean = sum(log_r) / len(log_r)
  1176. y_mean = sum(log_zr) / len(log_zr)
  1177. for (x, y) in zip(log_r, log_zr):
  1178. xy_cov += (x - x_mean) * (y - y_mean)
  1179. x_var += (x - x_mean) ** 2
  1180. self._slope = xy_cov / x_var if x_var != 0 else 0.0
  1181. if self._slope >= -1:
  1182. warnings.warn(
  1183. 'SimpleGoodTuring did not find a proper best fit '
  1184. 'line for smoothing probabilities of occurrences. '
  1185. 'The probability estimates are likely to be '
  1186. 'unreliable.'
  1187. )
  1188. self._intercept = y_mean - self._slope * x_mean
  1189. def _switch(self, r, nr):
  1190. """
  1191. Calculate the r frontier where we must switch from Nr to Sr
  1192. when estimating E[Nr].
  1193. """
  1194. for i, r_ in enumerate(r):
  1195. if len(r) == i + 1 or r[i + 1] != r_ + 1:
  1196. # We are at the end of r, or there is a gap in r
  1197. self._switch_at = r_
  1198. break
  1199. Sr = self.smoothedNr
  1200. smooth_r_star = (r_ + 1) * Sr(r_ + 1) / Sr(r_)
  1201. unsmooth_r_star = (r_ + 1) * nr[i + 1] / nr[i]
  1202. std = math.sqrt(self._variance(r_, nr[i], nr[i + 1]))
  1203. if abs(unsmooth_r_star - smooth_r_star) <= 1.96 * std:
  1204. self._switch_at = r_
  1205. break
  1206. def _variance(self, r, nr, nr_1):
  1207. r = float(r)
  1208. nr = float(nr)
  1209. nr_1 = float(nr_1)
  1210. return (r + 1.0) ** 2 * (nr_1 / nr ** 2) * (1.0 + nr_1 / nr)
  1211. def _renormalize(self, r, nr):
  1212. """
  1213. It is necessary to renormalize all the probability estimates to
  1214. ensure a proper probability distribution results. This can be done
  1215. by keeping the estimate of the probability mass for unseen items as
  1216. N(1)/N and renormalizing all the estimates for previously seen items
  1217. (as Gale and Sampson (1995) propose). (See M&S P.213, 1999)
  1218. """
  1219. prob_cov = 0.0
  1220. for r_, nr_ in zip(r, nr):
  1221. prob_cov += nr_ * self._prob_measure(r_)
  1222. if prob_cov:
  1223. self._renormal = (1 - self._prob_measure(0)) / prob_cov
  1224. def smoothedNr(self, r):
  1225. """
  1226. Return the number of samples with count r.
  1227. :param r: The amount of frequency.
  1228. :type r: int
  1229. :rtype: float
  1230. """
  1231. # Nr = a*r^b (with b < -1 to give the appropriate hyperbolic
  1232. # relationship)
  1233. # Estimate a and b by simple linear regression technique on
  1234. # the logarithmic form of the equation: log Nr = a + b*log(r)
  1235. return math.exp(self._intercept + self._slope * math.log(r))
  1236. def prob(self, sample):
  1237. """
  1238. Return the sample's probability.
  1239. :param sample: sample of the event
  1240. :type sample: str
  1241. :rtype: float
  1242. """
  1243. count = self._freqdist[sample]
  1244. p = self._prob_measure(count)
  1245. if count == 0:
  1246. if self._bins == self._freqdist.B():
  1247. p = 0.0
  1248. else:
  1249. p = p / (self._bins - self._freqdist.B())
  1250. else:
  1251. p = p * self._renormal
  1252. return p
  1253. def _prob_measure(self, count):
  1254. if count == 0 and self._freqdist.N() == 0:
  1255. return 1.0
  1256. elif count == 0 and self._freqdist.N() != 0:
  1257. return self._freqdist.Nr(1) / self._freqdist.N()
  1258. if self._switch_at > count:
  1259. Er_1 = self._freqdist.Nr(count + 1)
  1260. Er = self._freqdist.Nr(count)
  1261. else:
  1262. Er_1 = self.smoothedNr(count + 1)
  1263. Er = self.smoothedNr(count)
  1264. r_star = (count + 1) * Er_1 / Er
  1265. return r_star / self._freqdist.N()
  1266. def check(self):
  1267. prob_sum = 0.0
  1268. for i in range(0, len(self._Nr)):
  1269. prob_sum += self._Nr[i] * self._prob_measure(i) / self._renormal
  1270. print("Probability Sum:", prob_sum)
  1271. # assert prob_sum != 1.0, "probability sum should be one!"
  1272. def discount(self):
  1273. """
  1274. This function returns the total mass of probability transfers from the
  1275. seen samples to the unseen samples.
  1276. """
  1277. return self.smoothedNr(1) / self._freqdist.N()
  1278. def max(self):
  1279. return self._freqdist.max()
  1280. def samples(self):
  1281. return self._freqdist.keys()
  1282. def freqdist(self):
  1283. return self._freqdist
  1284. def __repr__(self):
  1285. """
  1286. Return a string representation of this ``ProbDist``.
  1287. :rtype: str
  1288. """
  1289. return '<SimpleGoodTuringProbDist based on %d samples>' % self._freqdist.N()
  1290. class MutableProbDist(ProbDistI):
  1291. """
  1292. An mutable probdist where the probabilities may be easily modified. This
  1293. simply copies an existing probdist, storing the probability values in a
  1294. mutable dictionary and providing an update method.
  1295. """
  1296. def __init__(self, prob_dist, samples, store_logs=True):
  1297. """
  1298. Creates the mutable probdist based on the given prob_dist and using
  1299. the list of samples given. These values are stored as log
  1300. probabilities if the store_logs flag is set.
  1301. :param prob_dist: the distribution from which to garner the
  1302. probabilities
  1303. :type prob_dist: ProbDist
  1304. :param samples: the complete set of samples
  1305. :type samples: sequence of any
  1306. :param store_logs: whether to store the probabilities as logarithms
  1307. :type store_logs: bool
  1308. """
  1309. self._samples = samples
  1310. self._sample_dict = dict((samples[i], i) for i in range(len(samples)))
  1311. self._data = array.array(str("d"), [0.0]) * len(samples)
  1312. for i in range(len(samples)):
  1313. if store_logs:
  1314. self._data[i] = prob_dist.logprob(samples[i])
  1315. else:
  1316. self._data[i] = prob_dist.prob(samples[i])
  1317. self._logs = store_logs
  1318. def max(self):
  1319. # inherit documentation
  1320. return max((p, v) for (v, p) in self._sample_dict.items())[1]
  1321. def samples(self):
  1322. # inherit documentation
  1323. return self._samples
  1324. def prob(self, sample):
  1325. # inherit documentation
  1326. i = self._sample_dict.get(sample)
  1327. if i is None:
  1328. return 0.0
  1329. return 2 ** (self._data[i]) if self._logs else self._data[i]
  1330. def logprob(self, sample):
  1331. # inherit documentation
  1332. i = self._sample_dict.get(sample)
  1333. if i is None:
  1334. return float('-inf')
  1335. return self._data[i] if self._logs else math.log(self._data[i], 2)
  1336. def update(self, sample, prob, log=True):
  1337. """
  1338. Update the probability for the given sample. This may cause the object
  1339. to stop being the valid probability distribution - the user must
  1340. ensure that they update the sample probabilities such that all samples
  1341. have probabilities between 0 and 1 and that all probabilities sum to
  1342. one.
  1343. :param sample: the sample for which to update the probability
  1344. :type sample: any
  1345. :param prob: the new probability
  1346. :type prob: float
  1347. :param log: is the probability already logged
  1348. :type log: bool
  1349. """
  1350. i = self._sample_dict.get(sample)
  1351. assert i is not None
  1352. if self._logs:
  1353. self._data[i] = prob if log else math.log(prob, 2)
  1354. else:
  1355. self._data[i] = 2 ** (prob) if log else prob
  1356. ##/////////////////////////////////////////////////////
  1357. ## Kneser-Ney Probability Distribution
  1358. ##//////////////////////////////////////////////////////
  1359. # This method for calculating probabilities was introduced in 1995 by Reinhard
  1360. # Kneser and Hermann Ney. It was meant to improve the accuracy of language
  1361. # models that use backing-off to deal with sparse data. The authors propose two
  1362. # ways of doing so: a marginal distribution constraint on the back-off
  1363. # distribution and a leave-one-out distribution. For a start, the first one is
  1364. # implemented as a class below.
  1365. #
  1366. # The idea behind a back-off n-gram model is that we have a series of
  1367. # frequency distributions for our n-grams so that in case we have not seen a
  1368. # given n-gram during training (and as a result have a 0 probability for it) we
  1369. # can 'back off' (hence the name!) and try testing whether we've seen the
  1370. # n-1-gram part of the n-gram in training.
  1371. #
  1372. # The novelty of Kneser and Ney's approach was that they decided to fiddle
  1373. # around with the way this latter, backed off probability was being calculated
  1374. # whereas their peers seemed to focus on the primary probability.
  1375. #
  1376. # The implementation below uses one of the techniques described in their paper
  1377. # titled "Improved backing-off for n-gram language modeling." In the same paper
  1378. # another technique is introduced to attempt to smooth the back-off
  1379. # distribution as well as the primary one. There is also a much-cited
  1380. # modification of this method proposed by Chen and Goodman.
  1381. #
  1382. # In order for the implementation of Kneser-Ney to be more efficient, some
  1383. # changes have been made to the original algorithm. Namely, the calculation of
  1384. # the normalizing function gamma has been significantly simplified and
  1385. # combined slightly differently with beta. None of these changes affect the
  1386. # nature of the algorithm, but instead aim to cut out unnecessary calculations
  1387. # and take advantage of storing and retrieving information in dictionaries
  1388. # where possible.
  1389. @compat.python_2_unicode_compatible
  1390. class KneserNeyProbDist(ProbDistI):
  1391. """
  1392. Kneser-Ney estimate of a probability distribution. This is a version of
  1393. back-off that counts how likely an n-gram is provided the n-1-gram had
  1394. been seen in training. Extends the ProbDistI interface, requires a trigram
  1395. FreqDist instance to train on. Optionally, a different from default discount
  1396. value can be specified. The default discount is set to 0.75.
  1397. """
  1398. def __init__(self, freqdist, bins=None, discount=0.75):
  1399. """
  1400. :param freqdist: The trigram frequency distribution upon which to base
  1401. the estimation
  1402. :type freqdist: FreqDist
  1403. :param bins: Included for compatibility with nltk.tag.hmm
  1404. :type bins: int or float
  1405. :param discount: The discount applied when retrieving counts of
  1406. trigrams
  1407. :type discount: float (preferred, but can be set to int)
  1408. """
  1409. if not bins:
  1410. self._bins = freqdist.B()
  1411. else:
  1412. self._bins = bins
  1413. self._D = discount
  1414. # cache for probability calculation
  1415. self._cache = {}
  1416. # internal bigram and trigram frequency distributions
  1417. self._bigrams = defaultdict(int)
  1418. self._trigrams = freqdist
  1419. # helper dictionaries used to calculate probabilities
  1420. self._wordtypes_after = defaultdict(float)
  1421. self._trigrams_contain = defaultdict(float)
  1422. self._wordtypes_before = defaultdict(float)
  1423. for w0, w1, w2 in freqdist:
  1424. self._bigrams[(w0, w1)] += freqdist[(w0, w1, w2)]
  1425. self._wordtypes_after[(w0, w1)] += 1
  1426. self._trigrams_contain[w1] += 1
  1427. self._wordtypes_before[(w1, w2)] += 1
  1428. def prob(self, trigram):
  1429. # sample must be a triple
  1430. if len(trigram) != 3:
  1431. raise ValueError('Expected an iterable with 3 members.')
  1432. trigram = tuple(trigram)
  1433. w0, w1, w2 = trigram
  1434. if trigram in self._cache:
  1435. return self._cache[trigram]
  1436. else:
  1437. # if the sample trigram was seen during training
  1438. if trigram in self._trigrams:
  1439. prob = (self._trigrams[trigram] - self.discount()) / self._bigrams[
  1440. (w0, w1)
  1441. ]
  1442. # else if the 'rougher' environment was seen during training
  1443. elif (w0, w1) in self._bigrams and (w1, w2) in self._wordtypes_before:
  1444. aftr = self._wordtypes_after[(w0, w1)]
  1445. bfr = self._wordtypes_before[(w1, w2)]
  1446. # the probability left over from alphas
  1447. leftover_prob = (aftr * self.discount()) / self._bigrams[(w0, w1)]
  1448. # the beta (including normalization)
  1449. beta = bfr / (self._trigrams_contain[w1] - aftr)
  1450. prob = leftover_prob * beta
  1451. # else the sample was completely unseen during training
  1452. else:
  1453. prob = 0.0
  1454. self._cache[trigram] = prob
  1455. return prob
  1456. def discount(self):
  1457. """
  1458. Return the value by which counts are discounted. By default set to 0.75.
  1459. :rtype: float
  1460. """
  1461. return self._D
  1462. def set_discount(self, discount):
  1463. """
  1464. Set the value by which counts are discounted to the value of discount.
  1465. :param discount: the new value to discount counts by
  1466. :type discount: float (preferred, but int possible)
  1467. :rtype: None
  1468. """
  1469. self._D = discount
  1470. def samples(self):
  1471. return self._trigrams.keys()
  1472. def max(self):
  1473. return self._trigrams.max()
  1474. def __repr__(self):
  1475. '''
  1476. Return a string representation of this ProbDist
  1477. :rtype: str
  1478. '''
  1479. return '<KneserNeyProbDist based on {0} trigrams'.format(self._trigrams.N())
  1480. ##//////////////////////////////////////////////////////
  1481. ## Probability Distribution Operations
  1482. ##//////////////////////////////////////////////////////
  1483. def log_likelihood(test_pdist, actual_pdist):
  1484. if not isinstance(test_pdist, ProbDistI) or not isinstance(actual_pdist, ProbDistI):
  1485. raise ValueError('expected a ProbDist.')
  1486. # Is this right?
  1487. return sum(
  1488. actual_pdist.prob(s) * math.log(test_pdist.prob(s), 2) for s in actual_pdist
  1489. )
  1490. def entropy(pdist):
  1491. probs = (pdist.prob(s) for s in pdist.samples())
  1492. return -sum(p * math.log(p, 2) for p in probs)
  1493. ##//////////////////////////////////////////////////////
  1494. ## Conditional Distributions
  1495. ##//////////////////////////////////////////////////////
  1496. @compat.python_2_unicode_compatible
  1497. class ConditionalFreqDist(defaultdict):
  1498. """
  1499. A collection of frequency distributions for a single experiment
  1500. run under different conditions. Conditional frequency
  1501. distributions are used to record the number of times each sample
  1502. occurred, given the condition under which the experiment was run.
  1503. For example, a conditional frequency distribution could be used to
  1504. record the frequency of each word (type) in a document, given its
  1505. length. Formally, a conditional frequency distribution can be
  1506. defined as a function that maps from each condition to the
  1507. FreqDist for the experiment under that condition.
  1508. Conditional frequency distributions are typically constructed by
  1509. repeatedly running an experiment under a variety of conditions,
  1510. and incrementing the sample outcome counts for the appropriate
  1511. conditions. For example, the following code will produce a
  1512. conditional frequency distribution that encodes how often each
  1513. word type occurs, given the length of that word type:
  1514. >>> from nltk.probability import ConditionalFreqDist
  1515. >>> from nltk.tokenize import word_tokenize
  1516. >>> sent = "the the the dog dog some other words that we do not care about"
  1517. >>> cfdist = ConditionalFreqDist()
  1518. >>> for word in word_tokenize(sent):
  1519. ... condition = len(word)
  1520. ... cfdist[condition][word] += 1
  1521. An equivalent way to do this is with the initializer:
  1522. >>> cfdist = ConditionalFreqDist((len(word), word) for word in word_tokenize(sent))
  1523. The frequency distribution for each condition is accessed using
  1524. the indexing operator:
  1525. >>> cfdist[3]
  1526. FreqDist({'the': 3, 'dog': 2, 'not': 1})
  1527. >>> cfdist[3].freq('the')
  1528. 0.5
  1529. >>> cfdist[3]['dog']
  1530. 2
  1531. When the indexing operator is used to access the frequency
  1532. distribution for a condition that has not been accessed before,
  1533. ``ConditionalFreqDist`` creates a new empty FreqDist for that
  1534. condition.
  1535. """
  1536. def __init__(self, cond_samples=None):
  1537. """
  1538. Construct a new empty conditional frequency distribution. In
  1539. particular, the count for every sample, under every condition,
  1540. is zero.
  1541. :param cond_samples: The samples to initialize the conditional
  1542. frequency distribution with
  1543. :type cond_samples: Sequence of (condition, sample) tuples
  1544. """
  1545. defaultdict.__init__(self, FreqDist)
  1546. if cond_samples:
  1547. for (cond, sample) in cond_samples:
  1548. self[cond][sample] += 1
  1549. def __reduce__(self):
  1550. kv_pairs = ((cond, self[cond]) for cond in self.conditions())
  1551. return (self.__class__, (), None, None, kv_pairs)
  1552. def conditions(self):
  1553. """
  1554. Return a list of the conditions that have been accessed for
  1555. this ``ConditionalFreqDist``. Use the indexing operator to
  1556. access the frequency distribution for a given condition.
  1557. Note that the frequency distributions for some conditions
  1558. may contain zero sample outcomes.
  1559. :rtype: list
  1560. """
  1561. return list(self.keys())
  1562. def N(self):
  1563. """
  1564. Return the total number of sample outcomes that have been
  1565. recorded by this ``ConditionalFreqDist``.
  1566. :rtype: int
  1567. """
  1568. return sum(fdist.N() for fdist in itervalues(self))
  1569. def plot(self, *args, **kwargs):
  1570. """
  1571. Plot the given samples from the conditional frequency distribution.
  1572. For a cumulative plot, specify cumulative=True.
  1573. (Requires Matplotlib to be installed.)
  1574. :param samples: The samples to plot
  1575. :type samples: list
  1576. :param title: The title for the graph
  1577. :type title: str
  1578. :param conditions: The conditions to plot (default is all)
  1579. :type conditions: list
  1580. """
  1581. try:
  1582. import matplotlib.pyplot as plt #import statment fix
  1583. except ImportError:
  1584. raise ValueError(
  1585. 'The plot function requires matplotlib to be installed.'
  1586. 'See http://matplotlib.org/'
  1587. )
  1588. cumulative = _get_kwarg(kwargs, 'cumulative', False)
  1589. percents = _get_kwarg(kwargs, 'percents', False)
  1590. conditions = [c for c in _get_kwarg(kwargs, 'conditions', self.conditions()) if c in self] # conditions should be in self
  1591. title = _get_kwarg(kwargs, 'title', '')
  1592. samples = _get_kwarg(
  1593. kwargs, 'samples', sorted(set(v
  1594. for c in conditions
  1595. for v in self[c]))
  1596. ) # this computation could be wasted
  1597. if "linewidth" not in kwargs:
  1598. kwargs["linewidth"] = 2
  1599. ax = plt.gca()
  1600. if (len(conditions) != 0):
  1601. freqs = []
  1602. for condition in conditions:
  1603. if cumulative:
  1604. # freqs should be a list of list where each sub list will be a frequency of a condition
  1605. freqs.append(list(self[condition]._cumulative_frequencies(samples)))
  1606. ylabel = "Cumulative Counts"
  1607. legend_loc = 'lower right'
  1608. if percents:
  1609. freqs[-1] = [f / freqs[len(freqs) - 1] * 100 for f in freqs]
  1610. ylabel = "Cumulative Percents"
  1611. else:
  1612. freqs.append([self[condition][sample] for sample in samples])
  1613. ylabel = "Counts"
  1614. legend_loc = 'upper right'
  1615. # percents = [f * 100 for f in freqs] only in ConditionalProbDist?
  1616. i = 0
  1617. for freq in freqs:
  1618. kwargs['label'] = conditions[i] #label for each condition
  1619. i += 1
  1620. ax.plot(freq, *args, **kwargs)
  1621. ax.legend(loc=legend_loc)
  1622. ax.grid(True, color="silver")
  1623. ax.set_xticks(range(len(samples)))
  1624. ax.set_xticklabels([text_type(s) for s in samples], rotation=90)
  1625. if title:
  1626. ax.set_title(title)
  1627. ax.set_xlabel("Samples")
  1628. ax.set_ylabel(ylabel)
  1629. plt.show()
  1630. return ax
  1631. def tabulate(self, *args, **kwargs):
  1632. """
  1633. Tabulate the given samples from the conditional frequency distribution.
  1634. :param samples: The samples to plot
  1635. :type samples: list
  1636. :param conditions: The conditions to plot (default is all)
  1637. :type conditions: list
  1638. :param cumulative: A flag to specify whether the freqs are cumulative (default = False)
  1639. :type title: bool
  1640. """
  1641. cumulative = _get_kwarg(kwargs, 'cumulative', False)
  1642. conditions = _get_kwarg(kwargs, 'conditions', sorted(self.conditions()))
  1643. samples = _get_kwarg(
  1644. kwargs, 'samples', sorted(set(v for c in conditions
  1645. if c in self
  1646. for v in self[c]))) # this computation could be wasted
  1647. width = max(len("%s" % s) for s in samples)
  1648. freqs = dict()
  1649. for c in conditions:
  1650. if cumulative:
  1651. freqs[c] = list(self[c]._cumulative_frequencies(samples))
  1652. else:
  1653. freqs[c] = [self[c][sample] for sample in samples]
  1654. width = max(width, max(len("%d" % f) for f in freqs[c]))
  1655. condition_size = max(len("%s" % c) for c in conditions)
  1656. print(' ' * condition_size, end=' ')
  1657. for s in samples:
  1658. print("%*s" % (width, s), end=' ')
  1659. print()
  1660. for c in conditions:
  1661. print("%*s" % (condition_size, c), end=' ')
  1662. for f in freqs[c]:
  1663. print("%*d" % (width, f), end=' ')
  1664. print()
  1665. # Mathematical operators
  1666. def __add__(self, other):
  1667. """
  1668. Add counts from two ConditionalFreqDists.
  1669. """
  1670. if not isinstance(other, ConditionalFreqDist):
  1671. return NotImplemented
  1672. result = ConditionalFreqDist()
  1673. for cond in self.conditions():
  1674. newfreqdist = self[cond] + other[cond]
  1675. if newfreqdist:
  1676. result[cond] = newfreqdist
  1677. for cond in other.conditions():
  1678. if cond not in self.conditions():
  1679. for elem, count in other[cond].items():
  1680. if count > 0:
  1681. result[cond][elem] = count
  1682. return result
  1683. def __sub__(self, other):
  1684. """
  1685. Subtract count, but keep only results with positive counts.
  1686. """
  1687. if not isinstance(other, ConditionalFreqDist):
  1688. return NotImplemented
  1689. result = ConditionalFreqDist()
  1690. for cond in self.conditions():
  1691. newfreqdist = self[cond] - other[cond]
  1692. if newfreqdist:
  1693. result[cond] = newfreqdist
  1694. for cond in other.conditions():
  1695. if cond not in self.conditions():
  1696. for elem, count in other[cond].items():
  1697. if count < 0:
  1698. result[cond][elem] = 0 - count
  1699. return result
  1700. def __or__(self, other):
  1701. """
  1702. Union is the maximum of value in either of the input counters.
  1703. """
  1704. if not isinstance(other, ConditionalFreqDist):
  1705. return NotImplemented
  1706. result = ConditionalFreqDist()
  1707. for cond in self.conditions():
  1708. newfreqdist = self[cond] | other[cond]
  1709. if newfreqdist:
  1710. result[cond] = newfreqdist
  1711. for cond in other.conditions():
  1712. if cond not in self.conditions():
  1713. for elem, count in other[cond].items():
  1714. if count > 0:
  1715. result[cond][elem] = count
  1716. return result
  1717. def __and__(self, other):
  1718. """
  1719. Intersection is the minimum of corresponding counts.
  1720. """
  1721. if not isinstance(other, ConditionalFreqDist):
  1722. return NotImplemented
  1723. result = ConditionalFreqDist()
  1724. for cond in self.conditions():
  1725. newfreqdist = self[cond] & other[cond]
  1726. if newfreqdist:
  1727. result[cond] = newfreqdist
  1728. return result
  1729. # @total_ordering doesn't work here, since the class inherits from a builtin class
  1730. def __le__(self, other):
  1731. if not isinstance(other, ConditionalFreqDist):
  1732. raise_unorderable_types("<=", self, other)
  1733. return set(self.conditions()).issubset(other.conditions()) and all(
  1734. self[c] <= other[c] for c in self.conditions()
  1735. )
  1736. def __lt__(self, other):
  1737. if not isinstance(other, ConditionalFreqDist):
  1738. raise_unorderable_types("<", self, other)
  1739. return self <= other and self != other
  1740. def __ge__(self, other):
  1741. if not isinstance(other, ConditionalFreqDist):
  1742. raise_unorderable_types(">=", self, other)
  1743. return other <= self
  1744. def __gt__(self, other):
  1745. if not isinstance(other, ConditionalFreqDist):
  1746. raise_unorderable_types(">", self, other)
  1747. return other < self
  1748. def __repr__(self):
  1749. """
  1750. Return a string representation of this ``ConditionalFreqDist``.
  1751. :rtype: str
  1752. """
  1753. return '<ConditionalFreqDist with %d conditions>' % len(self)
  1754. @compat.python_2_unicode_compatible
  1755. @add_metaclass(ABCMeta)
  1756. class ConditionalProbDistI(dict):
  1757. """
  1758. A collection of probability distributions for a single experiment
  1759. run under different conditions. Conditional probability
  1760. distributions are used to estimate the likelihood of each sample,
  1761. given the condition under which the experiment was run. For
  1762. example, a conditional probability distribution could be used to
  1763. estimate the probability of each word type in a document, given
  1764. the length of the word type. Formally, a conditional probability
  1765. distribution can be defined as a function that maps from each
  1766. condition to the ``ProbDist`` for the experiment under that
  1767. condition.
  1768. """
  1769. @abstractmethod
  1770. def __init__(self):
  1771. """
  1772. Classes inheriting from ConditionalProbDistI should implement __init__.
  1773. """
  1774. def conditions(self):
  1775. """
  1776. Return a list of the conditions that are represented by
  1777. this ``ConditionalProbDist``. Use the indexing operator to
  1778. access the probability distribution for a given condition.
  1779. :rtype: list
  1780. """
  1781. return list(self.keys())
  1782. def __repr__(self):
  1783. """
  1784. Return a string representation of this ``ConditionalProbDist``.
  1785. :rtype: str
  1786. """
  1787. return '<%s with %d conditions>' % (type(self).__name__, len(self))
  1788. class ConditionalProbDist(ConditionalProbDistI):
  1789. """
  1790. A conditional probability distribution modeling the experiments
  1791. that were used to generate a conditional frequency distribution.
  1792. A ConditionalProbDist is constructed from a
  1793. ``ConditionalFreqDist`` and a ``ProbDist`` factory:
  1794. - The ``ConditionalFreqDist`` specifies the frequency
  1795. distribution for each condition.
  1796. - The ``ProbDist`` factory is a function that takes a
  1797. condition's frequency distribution, and returns its
  1798. probability distribution. A ``ProbDist`` class's name (such as
  1799. ``MLEProbDist`` or ``HeldoutProbDist``) can be used to specify
  1800. that class's constructor.
  1801. The first argument to the ``ProbDist`` factory is the frequency
  1802. distribution that it should model; and the remaining arguments are
  1803. specified by the ``factory_args`` parameter to the
  1804. ``ConditionalProbDist`` constructor. For example, the following
  1805. code constructs a ``ConditionalProbDist``, where the probability
  1806. distribution for each condition is an ``ELEProbDist`` with 10 bins:
  1807. >>> from nltk.corpus import brown
  1808. >>> from nltk.probability import ConditionalFreqDist
  1809. >>> from nltk.probability import ConditionalProbDist, ELEProbDist
  1810. >>> cfdist = ConditionalFreqDist(brown.tagged_words()[:5000])
  1811. >>> cpdist = ConditionalProbDist(cfdist, ELEProbDist, 10)
  1812. >>> cpdist['passed'].max()
  1813. 'VBD'
  1814. >>> cpdist['passed'].prob('VBD')
  1815. 0.423...
  1816. """
  1817. def __init__(self, cfdist, probdist_factory, *factory_args, **factory_kw_args):
  1818. """
  1819. Construct a new conditional probability distribution, based on
  1820. the given conditional frequency distribution and ``ProbDist``
  1821. factory.
  1822. :type cfdist: ConditionalFreqDist
  1823. :param cfdist: The ``ConditionalFreqDist`` specifying the
  1824. frequency distribution for each condition.
  1825. :type probdist_factory: class or function
  1826. :param probdist_factory: The function or class that maps
  1827. a condition's frequency distribution to its probability
  1828. distribution. The function is called with the frequency
  1829. distribution as its first argument,
  1830. ``factory_args`` as its remaining arguments, and
  1831. ``factory_kw_args`` as keyword arguments.
  1832. :type factory_args: (any)
  1833. :param factory_args: Extra arguments for ``probdist_factory``.
  1834. These arguments are usually used to specify extra
  1835. properties for the probability distributions of individual
  1836. conditions, such as the number of bins they contain.
  1837. :type factory_kw_args: (any)
  1838. :param factory_kw_args: Extra keyword arguments for ``probdist_factory``.
  1839. """
  1840. self._probdist_factory = probdist_factory
  1841. self._factory_args = factory_args
  1842. self._factory_kw_args = factory_kw_args
  1843. for condition in cfdist:
  1844. self[condition] = probdist_factory(
  1845. cfdist[condition], *factory_args, **factory_kw_args
  1846. )
  1847. def __missing__(self, key):
  1848. self[key] = self._probdist_factory(
  1849. FreqDist(), *self._factory_args, **self._factory_kw_args
  1850. )
  1851. return self[key]
  1852. class DictionaryConditionalProbDist(ConditionalProbDistI):
  1853. """
  1854. An alternative ConditionalProbDist that simply wraps a dictionary of
  1855. ProbDists rather than creating these from FreqDists.
  1856. """
  1857. def __init__(self, probdist_dict):
  1858. """
  1859. :param probdist_dict: a dictionary containing the probdists indexed
  1860. by the conditions
  1861. :type probdist_dict: dict any -> probdist
  1862. """
  1863. self.update(probdist_dict)
  1864. def __missing__(self, key):
  1865. self[key] = DictionaryProbDist()
  1866. return self[key]
  1867. ##//////////////////////////////////////////////////////
  1868. ## Adding in log-space.
  1869. ##//////////////////////////////////////////////////////
  1870. # If the difference is bigger than this, then just take the bigger one:
  1871. _ADD_LOGS_MAX_DIFF = math.log(1e-30, 2)
  1872. def add_logs(logx, logy):
  1873. """
  1874. Given two numbers ``logx`` = *log(x)* and ``logy`` = *log(y)*, return
  1875. *log(x+y)*. Conceptually, this is the same as returning
  1876. ``log(2**(logx)+2**(logy))``, but the actual implementation
  1877. avoids overflow errors that could result from direct computation.
  1878. """
  1879. if logx < logy + _ADD_LOGS_MAX_DIFF:
  1880. return logy
  1881. if logy < logx + _ADD_LOGS_MAX_DIFF:
  1882. return logx
  1883. base = min(logx, logy)
  1884. return base + math.log(2 ** (logx - base) + 2 ** (logy - base), 2)
  1885. def sum_logs(logs):
  1886. return reduce(add_logs, logs[1:], logs[0]) if len(logs) != 0 else _NINF
  1887. ##//////////////////////////////////////////////////////
  1888. ## Probabilistic Mix-in
  1889. ##//////////////////////////////////////////////////////
  1890. class ProbabilisticMixIn(object):
  1891. """
  1892. A mix-in class to associate probabilities with other classes
  1893. (trees, rules, etc.). To use the ``ProbabilisticMixIn`` class,
  1894. define a new class that derives from an existing class and from
  1895. ProbabilisticMixIn. You will need to define a new constructor for
  1896. the new class, which explicitly calls the constructors of both its
  1897. parent classes. For example:
  1898. >>> from nltk.probability import ProbabilisticMixIn
  1899. >>> class A:
  1900. ... def __init__(self, x, y): self.data = (x,y)
  1901. ...
  1902. >>> class ProbabilisticA(A, ProbabilisticMixIn):
  1903. ... def __init__(self, x, y, **prob_kwarg):
  1904. ... A.__init__(self, x, y)
  1905. ... ProbabilisticMixIn.__init__(self, **prob_kwarg)
  1906. See the documentation for the ProbabilisticMixIn
  1907. ``constructor<__init__>`` for information about the arguments it
  1908. expects.
  1909. You should generally also redefine the string representation
  1910. methods, the comparison methods, and the hashing method.
  1911. """
  1912. def __init__(self, **kwargs):
  1913. """
  1914. Initialize this object's probability. This initializer should
  1915. be called by subclass constructors. ``prob`` should generally be
  1916. the first argument for those constructors.
  1917. :param prob: The probability associated with the object.
  1918. :type prob: float
  1919. :param logprob: The log of the probability associated with
  1920. the object.
  1921. :type logprob: float
  1922. """
  1923. if 'prob' in kwargs:
  1924. if 'logprob' in kwargs:
  1925. raise TypeError('Must specify either prob or logprob ' '(not both)')
  1926. else:
  1927. ProbabilisticMixIn.set_prob(self, kwargs['prob'])
  1928. elif 'logprob' in kwargs:
  1929. ProbabilisticMixIn.set_logprob(self, kwargs['logprob'])
  1930. else:
  1931. self.__prob = self.__logprob = None
  1932. def set_prob(self, prob):
  1933. """
  1934. Set the probability associated with this object to ``prob``.
  1935. :param prob: The new probability
  1936. :type prob: float
  1937. """
  1938. self.__prob = prob
  1939. self.__logprob = None
  1940. def set_logprob(self, logprob):
  1941. """
  1942. Set the log probability associated with this object to
  1943. ``logprob``. I.e., set the probability associated with this
  1944. object to ``2**(logprob)``.
  1945. :param logprob: The new log probability
  1946. :type logprob: float
  1947. """
  1948. self.__logprob = logprob
  1949. self.__prob = None
  1950. def prob(self):
  1951. """
  1952. Return the probability associated with this object.
  1953. :rtype: float
  1954. """
  1955. if self.__prob is None:
  1956. if self.__logprob is None:
  1957. return None
  1958. self.__prob = 2 ** (self.__logprob)
  1959. return self.__prob
  1960. def logprob(self):
  1961. """
  1962. Return ``log(p)``, where ``p`` is the probability associated
  1963. with this object.
  1964. :rtype: float
  1965. """
  1966. if self.__logprob is None:
  1967. if self.__prob is None:
  1968. return None
  1969. self.__logprob = math.log(self.__prob, 2)
  1970. return self.__logprob
  1971. class ImmutableProbabilisticMixIn(ProbabilisticMixIn):
  1972. def set_prob(self, prob):
  1973. raise ValueError('%s is immutable' % self.__class__.__name__)
  1974. def set_logprob(self, prob):
  1975. raise ValueError('%s is immutable' % self.__class__.__name__)
  1976. ## Helper function for processing keyword arguments
  1977. def _get_kwarg(kwargs, key, default):
  1978. if key in kwargs:
  1979. arg = kwargs[key]
  1980. del kwargs[key]
  1981. else:
  1982. arg = default
  1983. return arg
  1984. ##//////////////////////////////////////////////////////
  1985. ## Demonstration
  1986. ##//////////////////////////////////////////////////////
  1987. def _create_rand_fdist(numsamples, numoutcomes):
  1988. """
  1989. Create a new frequency distribution, with random samples. The
  1990. samples are numbers from 1 to ``numsamples``, and are generated by
  1991. summing two numbers, each of which has a uniform distribution.
  1992. """
  1993. fdist = FreqDist()
  1994. for x in range(numoutcomes):
  1995. y = random.randint(1, (1 + numsamples) // 2) + random.randint(
  1996. 0, numsamples // 2
  1997. )
  1998. fdist[y] += 1
  1999. return fdist
  2000. def _create_sum_pdist(numsamples):
  2001. """
  2002. Return the true probability distribution for the experiment
  2003. ``_create_rand_fdist(numsamples, x)``.
  2004. """
  2005. fdist = FreqDist()
  2006. for x in range(1, (1 + numsamples) // 2 + 1):
  2007. for y in range(0, numsamples // 2 + 1):
  2008. fdist[x + y] += 1
  2009. return MLEProbDist(fdist)
  2010. def demo(numsamples=6, numoutcomes=500):
  2011. """
  2012. A demonstration of frequency distributions and probability
  2013. distributions. This demonstration creates three frequency
  2014. distributions with, and uses them to sample a random process with
  2015. ``numsamples`` samples. Each frequency distribution is sampled
  2016. ``numoutcomes`` times. These three frequency distributions are
  2017. then used to build six probability distributions. Finally, the
  2018. probability estimates of these distributions are compared to the
  2019. actual probability of each sample.
  2020. :type numsamples: int
  2021. :param numsamples: The number of samples to use in each demo
  2022. frequency distributions.
  2023. :type numoutcomes: int
  2024. :param numoutcomes: The total number of outcomes for each
  2025. demo frequency distribution. These outcomes are divided into
  2026. ``numsamples`` bins.
  2027. :rtype: None
  2028. """
  2029. # Randomly sample a stochastic process three times.
  2030. fdist1 = _create_rand_fdist(numsamples, numoutcomes)
  2031. fdist2 = _create_rand_fdist(numsamples, numoutcomes)
  2032. fdist3 = _create_rand_fdist(numsamples, numoutcomes)
  2033. # Use our samples to create probability distributions.
  2034. pdists = [
  2035. MLEProbDist(fdist1),
  2036. LidstoneProbDist(fdist1, 0.5, numsamples),
  2037. HeldoutProbDist(fdist1, fdist2, numsamples),
  2038. HeldoutProbDist(fdist2, fdist1, numsamples),
  2039. CrossValidationProbDist([fdist1, fdist2, fdist3], numsamples),
  2040. SimpleGoodTuringProbDist(fdist1),
  2041. SimpleGoodTuringProbDist(fdist1, 7),
  2042. _create_sum_pdist(numsamples),
  2043. ]
  2044. # Find the probability of each sample.
  2045. vals = []
  2046. for n in range(1, numsamples + 1):
  2047. vals.append(tuple([n, fdist1.freq(n)] + [pdist.prob(n) for pdist in pdists]))
  2048. # Print the results in a formatted table.
  2049. print(
  2050. (
  2051. '%d samples (1-%d); %d outcomes were sampled for each FreqDist'
  2052. % (numsamples, numsamples, numoutcomes)
  2053. )
  2054. )
  2055. print('=' * 9 * (len(pdists) + 2))
  2056. FORMATSTR = ' FreqDist ' + '%8s ' * (len(pdists) - 1) + '| Actual'
  2057. print(FORMATSTR % tuple(repr(pdist)[1:9] for pdist in pdists[:-1]))
  2058. print('-' * 9 * (len(pdists) + 2))
  2059. FORMATSTR = '%3d %8.6f ' + '%8.6f ' * (len(pdists) - 1) + '| %8.6f'
  2060. for val in vals:
  2061. print(FORMATSTR % val)
  2062. # Print the totals for each column (should all be 1.0)
  2063. zvals = list(zip(*vals))
  2064. sums = [sum(val) for val in zvals[1:]]
  2065. print('-' * 9 * (len(pdists) + 2))
  2066. FORMATSTR = 'Total ' + '%8.6f ' * (len(pdists)) + '| %8.6f'
  2067. print(FORMATSTR % tuple(sums))
  2068. print('=' * 9 * (len(pdists) + 2))
  2069. # Display the distributions themselves, if they're short enough.
  2070. if len("%s" % fdist1) < 70:
  2071. print(' fdist1: %s' % fdist1)
  2072. print(' fdist2: %s' % fdist2)
  2073. print(' fdist3: %s' % fdist3)
  2074. print()
  2075. print('Generating:')
  2076. for pdist in pdists:
  2077. fdist = FreqDist(pdist.generate() for i in range(5000))
  2078. print('%20s %s' % (pdist.__class__.__name__[:20], ("%s" % fdist)[:55]))
  2079. print()
  2080. def gt_demo():
  2081. from nltk import corpus
  2082. emma_words = corpus.gutenberg.words('austen-emma.txt')
  2083. fd = FreqDist(emma_words)
  2084. sgt = SimpleGoodTuringProbDist(fd)
  2085. print('%18s %8s %14s' % ("word", "freqency", "SimpleGoodTuring"))
  2086. fd_keys_sorted = (
  2087. key for key, value in sorted(fd.items(), key=lambda item: item[1], reverse=True)
  2088. )
  2089. for key in fd_keys_sorted:
  2090. print('%18s %8d %14e' % (key, fd[key], sgt.prob(key)))
  2091. if __name__ == '__main__':
  2092. demo(6, 10)
  2093. demo(5, 5000)
  2094. gt_demo()
  2095. __all__ = [
  2096. 'ConditionalFreqDist',
  2097. 'ConditionalProbDist',
  2098. 'ConditionalProbDistI',
  2099. 'CrossValidationProbDist',
  2100. 'DictionaryConditionalProbDist',
  2101. 'DictionaryProbDist',
  2102. 'ELEProbDist',
  2103. 'FreqDist',
  2104. 'SimpleGoodTuringProbDist',
  2105. 'HeldoutProbDist',
  2106. 'ImmutableProbabilisticMixIn',
  2107. 'LaplaceProbDist',
  2108. 'LidstoneProbDist',
  2109. 'MLEProbDist',
  2110. 'MutableProbDist',
  2111. 'KneserNeyProbDist',
  2112. 'ProbDistI',
  2113. 'ProbabilisticMixIn',
  2114. 'UniformProbDist',
  2115. 'WittenBellProbDist',
  2116. 'add_logs',
  2117. 'log_likelihood',
  2118. 'sum_logs',
  2119. 'entropy',
  2120. ]