2010-07-09 12 views
51

Wie kann ich die empirische CDF eines Arrays von Zahlen in Matplotlib in Python plotten? Ich suche nach dem cdf-Analog von pylabs "hist" -Funktion.Wie man empirische CDF in Matplotlib in Python plotten?

Eines, was ich denken kann ist:

from scipy.stats import cumfreq 
a = array([...]) # my array of numbers 
num_bins = 20 
b = cumfreq(a, num_bins) 
plt.plot(b) 

Ist allerdings das richtig? Gibt es einen leichteren/besseren Weg?

danke.

Antwort

15

Das sieht (fast) genau das sein, was Sie wollen. Zwei Dinge:

Zunächst werden die Ergebnisse sind ein Tupel von vier Elementen. Der dritte ist die Größe der Behälter. Die zweite ist der Startpunkt der kleinsten Bin. Die erste ist die Anzahl der Punkte in der oder unter jeder Tonne. (Die letzte ist die Anzahl der Punkte außerhalb der Grenzen, aber da Sie keine festgelegt haben, werden alle Punkte in Binned.)

Zweitens wollen Sie die Ergebnisse so skalieren, dass der endgültige Wert 1 ist, um folgen Sie den üblichen Konventionen eines CDF, aber ansonsten ist es richtig.

Hier ist, was es tut, unter der Haube:

def cumfreq(a, numbins=10, defaultreallimits=None): 
    # docstring omitted 
    h,l,b,e = histogram(a,numbins,defaultreallimits) 
    cumhist = np.cumsum(h*1, axis=0) 
    return cumhist,l,b,e 

Es macht die Histogrammierung, erzeugt dann eine kumulative Summe der Zählungen in jedem Fach. Der i-te Wert des Ergebnisses ist also die Anzahl der Array-Werte, die kleiner oder gleich dem Maximum des i-ten Bin ist. Der endgültige Wert ist also nur die Größe des ursprünglichen Arrays.

Schließlich ist es zu zeichnen, werden Sie den Anfangswert des Faches verwenden müssen, und die Binabmessung zu bestimmen, welche x-Achsen-Werte Sie benötigen.

Eine weitere Option ist numpy.histogram zu verwenden, die die Normalisierung tun können, und gibt die bin Kanten. Sie müssen die kumulative Summe der resultierenden Zählungen selbst durchführen.

a = array([...]) # your array of numbers 
num_bins = 20 
counts, bin_edges = numpy.histogram(a, bins=num_bins, normed=True) 
cdf = numpy.cumsum(counts) 
pylab.plot(bin_edges[1:], cdf) 

(bin_edges[1:] ist der obere Rand jedes Fachs.)

+17

Nur eine kurze Notiz: dieser Code tatsächlich Ihnen die empirischen CDF nicht geben (eine Funktion Schritt Erhöhung um 1/n bei jedem der n Datenpunkte). Stattdessen gibt dieser Code eine Schätzung der CDF basierend auf einer histogrammbasierten Schätzung der PDF. Diese histogrammbasierte Schätzung kann durch sorgfältige/unpassende Auswahl der Bins manipuliert/vorgespannt werden, so dass es eine nicht so gute Charakterisierung der echten CDF wie der tatsächlichen ECDF ist. –

+2

Ich mag auch nicht den Punkt, dass dies Binning auferlegt; siehe Dave's kurze Antwort, die einfach 'numpy.sort' verwendet, um die CDF ohne Binning zu plotten. –

3

Was tun Sie mit dem CDF tun wollen? Um es zu plotten, das ist ein Anfang. Sie könnten ein paar verschiedene Werte versuchen, wie folgt aus:

from __future__ import division 
import numpy as np 
from scipy.stats import cumfreq 
import pylab as plt 

hi = 100. 
a = np.arange(hi) ** 2 
for nbins in (2, 20, 100): 
    cf = cumfreq(a, nbins) # bin values, lowerlimit, binsize, extrapoints 
    w = hi/nbins 
    x = np.linspace(w/2, hi - w/2, nbins) # care 
    # print x, cf 
    plt.plot(x, cf[0], label=str(nbins)) 

plt.legend() 
plt.show() 

Histogram listet verschiedene Regeln für die Anzahl von Behältern, zum Beispiel num_bins ~ sqrt(len(a)).

(Fine Druck: zwei ganz verschiedene Dinge laufen hier auf,

  • Binning/Histogramm die Rohdaten
  • plot eine glatte Kurve durch die klassierte Werte sagen 20 interpoliert

. Jede diesen gehen kann weg auf Daten, die „klumpig“ oder haben langen Schwanz sind, auch für 1d Daten - 2d, werden 3D-Daten immer schwieriger
Siehe auch 012.351.Density_estimation und using scipy gaussian kernel density estimation ).

65

können Sie die ECDF Funktion aus der scikits.statsmodels Bibliothek verwenden:

import numpy as np 
import scikits.statsmodels as sm 
import matplotlib.pyplot as plt 

sample = np.random.uniform(0, 1, 50) 
ecdf = sm.tools.ECDF(sample) 

x = np.linspace(min(sample), max(sample)) 
y = ecdf(x) 
plt.step(x, y) 

Mit Version 0.4 scicits.statsmodels wurde statsmodels umbenannt. ECDF befindet sich jetzt im Modul distributions (während statsmodels.tools.tools.ECDF abgeschrieben wird).

import numpy as np 
import statsmodels.api as sm # recommended import according to the docs 
import matplotlib.pyplot as plt 

sample = np.random.uniform(0, 1, 50) 
ecdf = sm.distributions.ECDF(sample) 

x = np.linspace(min(sample), max(sample)) 
y = ecdf(x) 
plt.step(x, y) 
plt.show() 
+2

@bmu (und @Luca): genial; Danke, dass du den Code mit dem aktuellen statsmodel aktuell gemacht hast! – ars

+0

Für scikits.statsmodels v0.3.1 musste 'scikits.statsmodels.tools als smtools' importieren und' ecdf = smtools.tools.EDCF (...) ' – alexei

3

Ich habe eine triviale neben AFoglia Methode, die CDF

n_counts,bin_edges = np.histogram(myarray,bins=11,normed=True) 
cdf = np.cumsum(n_counts) # cdf not normalized, despite above 
scale = 1.0/cdf[-1] 
ncdf = scale * cdf 

Normalisieren der histo macht seine Integral Einheit zu normalisieren, die die CDF bedeutet nicht normalisiert werden. Du musst es selbst skalieren.

13

Haben Sie das kumulative = True-Argument für pyplot.hist versucht?

+1

Sehr gute Bemerkung. Dies erzwingt jedoch Binning; Siehe Daves Antwort mit np.sort. –

+0

Schöne und einfache Option, aber der Nachteil ist begrenzte Anpassung der resultierenden Liniendiagramm, z. konnte nicht herausfinden, wie man Marker hinzufügt. Ging für 'scikits.statsmodels' Antwort. – alexei

62

Wenn Sie linspace mögen und bevorzugen Einzeiler, können Sie tun:

plt.plot(np.sort(a), np.linspace(0, 1, len(a), endpoint=False)) 

meinen Geschmack gegeben, ich habe fast immer:

# a is the data array 
sorted_ = np.sort(a) 
yvals = np.arange(len(sorted_))/float(len(sorted_)) 
plt.plot(sorted_, yvals) 

was für mich funktioniert, auch wenn es >O(1e6) Datenwerte. Wenn Sie wirklich unten Probe benötigen würde ich

sorted_ = np.sort(a)[::down_sampling_step] 

bearbeiten reagieren gesetzt/bearbeiten zu kommentieren, warum ich endpoint=False oder die yvals wie oben definiert. Im Folgenden sind einige technische Details.

Die empirische CDF ist in der Regel formal als

CDF(x) = "number of samples <= x"/"number of samples" 

definiert, um genau diese formale Definition passen Sie yvals = np.arange(1,len(sorted_)+1)/float(len(sorted_)) verwenden brauchen würde, damit wir yvals = [1/N, 2/N ... 1] bekommen. Dieser Schätzer ist ein unverzerrter Schätzer, der in der Grenze der unendlichen Abtastwerte Wikipedia ref. mit der wahren CDF konvergiert.

Ich neige dazu, yvals = [0, 1/N, 2/N ... (N-1)/N] zu verwenden, da (a) es einfacher, Code/mehr idomatic ist, (b), aber immer noch formal gerechtfertigt, da man immer CDF(x) mit 1-CDF(x) im Konvergenz Beweis austauschen kann, und (c) arbeitet mit der (einfaches) Downsampling-Verfahren wie oben beschrieben.

In einigen besonderen Fällen ist es sinnvoll

yvals = (arange(len(sorted_))+0.5)/len(sorted_) 

, die zwischen diesen beiden Zwischen Konventionen zu definieren ist. In der Tat sagt es "gibt es eine 1/(2N) Chance von einem Wert weniger als die niedrigste, die ich in meiner Probe gesehen habe, und eine 1/(2N) Chance von einem Wert größer als die größte, die ich bisher gesehen habe.

Für große Stichproben und vernünftige Verteilungen ist die im Hauptteil der Antwort gegebene Konvention jedoch leicht zu schreiben, ist ein unvoreingenommener Schätzer der wahren CDF und arbeitet mit der Downsampling-Methode.

+3

Diese Antwort sollte mehr Upvotes erhalten, da es bisher die einzige ist, die kein Binning erzwingt. Ich habe den Code mit linspace nur ein wenig vereinfacht. –

+1

@hans_meine Ihre Bearbeitung, d.h. 'yvals = linspace (0,1, len (sortierte))', erzeugt "yvals", die keine unvoreingenommene Schätzfunktion der wahren CDF sind. – Dave

+0

Dann hätten wir linspace mit 'endpoint = False' verwenden sollen, oder? –

3

Wenn Sie die tatsächliche wahre ECDF anzeigen möchten (was David B bemerkte, ist eine Schrittfunktion, die 1/n an jedem der n Datenpunkte erhöht), ist mein Vorschlag, Code zu schreiben, um zwei "Plot" -Punkte zu generieren Datenpunkt:

a = array([...]) # your array of numbers 
sorted=np.sort(a) 
x2 = [] 
y2 = [] 
y = 0 
for x in sorted: 
    x2.extend([x,x]) 
    y2.append(y) 
    y += 1.0/len(a) 
    y2.append(y) 
plt.plot(x2,y2) 

auf diese Weise erhalten Sie eine grafische Darstellung mit den n Stufen erhalten, die von einem ECDF charakteristisch sind, die für die Datensätze schön ist besonders, die klein genug sind für die Stufen sichtbar zu sein. Außerdem ist es nicht notwendig, ein Binning mit Histogrammen durchzuführen (was die Einführung einer Verzerrung in die gezogene ECDF riskiert).

2

können wir einfach verwenden, um die Funktion von stepmatplotlib, die ein stufenweise Stück macht, die die Definition der empirischen CDF ist:

import numpy as np 
from matplotlib import pyplot as plt 

data = np.random.randn(11) 

levels = np.linspace(0, 1, len(data) + 1) # endpoint 1 is included by default 
plt.step(sorted(list(data) + [max(data)]), levels) 

Die letzte vertikale Linie bei max(data) wurde manuell hinzugefügt. Ansonsten stoppt die Handlung einfach bei Level 1 - 1/len(data).

Alternativ können wir die where='post' Option step()

levels = np.linspace(1./len(data), 1, len(data)) 
plt.step(sorted(data), levels, where='post') 

wobei in diesem Fall die anfängliche vertikale Linie von Null ist nicht aufgetragen verwenden.

1

(Dies ist eine Kopie meiner Antwort auf die Frage: Plotting CDF of a pandas series in python)

A CDF oder kumulative Verteilungsfunktion Grundstück ist im Grunde ein Diagramm mit auf der X-Achse der sortierten Werte und auf der Y-Achse der kumulativen Verteilung. Also würde ich eine neue Serie mit den sortierten Werten als Index und der kumulativen Verteilung als Werte erstellen.

Zuerst ein Beispiel Serie erstellen:

import pandas as pd 
import numpy as np 
ser = pd.Series(np.random.normal(size=100)) 

Sortieren der Serie:

ser = ser.order() 

Nun, bevor Sie fortfahren, fügen Sie wieder die letzte (und größte) Wert. Dieser Schritt ist wichtig, vor allem für kleine Probengrößen, um eine unvoreingenommene CDF zu erhalten:

ser[len(ser)] = ser.iloc[-1] 

Erstellen Sie eine neue Serie mit den sortierten Werten als Index und der kumulativen Verteilung als

cum_dist = np.linspace(0.,1.,len(ser)) 
ser_cdf = pd.Series(cum_dist, index=ser) 

Schließlich schätzt, Grundstück die Funktion als Schritte:

ser_cdf.plot(drawstyle='steps') 
5

One-Liner basiert auf Dave Antwort:

plt.plot(np.sort(arr), np.linspace(0, 1, len(arr), endpoint=False)) 

Edit: Dies wurde auch von Hans_meine in den Kommentaren vorgeschlagen.

+1

Dies ist die einfachste Antwort, das Problem elegant zu lösen.Dies sollte die akzeptierte Antwort sein! – Alex

1

Dies ist mit Bokeh

`` `

from bokeh.plotting import figure, show 
from statsmodels.distributions.empirical_distribution import ECDF 
ecdf = ECDF(pd_series) 
p = figure(title="tests", tools="save", background_fill_color="#E8DDCB") 
p.line(ecdf.x,ecdf.y) 
show(p) 

` ``

1

Unter der Annahme, dass vals Ihre Werte hält, dann können Sie einfach die CDF plotten wie folgt:

y = numpy.arange(0, 101) 
x = numpy.percentile(vals, y) 
plot(x, y) 

Um es zwischen 0 und 1 zu skalieren, teilen Sie einfach y durch 100.

0

Es ist ein Einstrich im Seebornen mit dem kumulativen = True-Parameter. Hier gehen Sie,

import seaborn as sns 
sns.kdeplot(a, cumulative=True) 
0

Keine der Antworten bisher bedeckt, was ich wollte, als ich hier gelandet, das ist:

def empirical_cdf(x, data): 
    "evaluate ecdf of data at points x" 
    return np.mean(data[None, :] <= x[:, None], axis=1) 

Er wertet die empirische CDF eines bestimmten Datensatzes an einer Reihe von Punkten x, die nicht sortiert werden müssen. Es gibt kein intermediäres Binning und keine externen Bibliotheken.

Eine äquivalente Methode, die besser für große x skaliert ist, die Daten zu sortieren und verwenden np.searchsorted:

def empirical_cdf(x, data): 
    "evaluate ecdf of data at points x" 
    data = np.sort(data) 
    return np.searchsorted(data, x)/float(data.size)