Ein einfacher nichtparametrischer Signifikanztest für beliebige Verteilungen Ingenieurwissenschaften

Es gibt zwei Datensätze, die man vergleichen möchte. Beide weisen eine Variation der Werte auf. Es stellt sich die Frage, ob sich die Datensätze signifikant oder nur zufällig unterscheiden? Sinnvollerweise geht es bei einer Datensichtung zunächst um diese Ja- oder Nein-Frage, bevor man nach der Größe des eventuellen Unterschieds fragt und versucht herauszufinden, ob es diesen tatsächlich gibt und was die Ursache dafür sein könnte. Eins sollte klar sein. Ein Signifikanztest ist – auch wenn es das Wort Test vielleicht suggeriert – kein Beweis für einen tatsächlichen oder ursächlichen Unterschied. Er liefert nur ein Indiz. Und wie das bei einem Indizienprozess der Fall ist, braucht man viele Indizien. Allerdings lässt sich bei wirklich vielen Indizien meist auch ein Beweis finden, der in jedem Fall zu bevorzugen ist.

Damit überhaupt von einem Test gesprochen werden kann, wird zunächst etwas verklausuliert eine Hypothese über die Datensätze aufgestellt. Das macht es allerdings nicht besser, denn der Signifikanztest offenbart nicht, ob und wie wahrscheinlich diese Hypothese zutrifft. J. Cohen formuliert es so:

What’s wrong with significance testing? Well, among many other things, it does not tell us what we want to know, and we so much want to know what we want know that, out of desperation, we nevertheless believe in that it does! What we want to know is „Given these data, what is the probability that the hypothesis is true?“ But as most of us know, what it tells us is „Given that hypothesis is true, what is the probability of these (or more extreme) data?“

Die Begriffe Test und Hypothese suggerieren irreführend eine mathematische Beweisführung. Statistische Signifikanz ist aber so zu verstehen, wie sie vom Erfinder des Signifikanztests Ronald Fischer gedacht war. Als erste Kennzahl, um nicht der Heuristik zu erliegen, im Rauschen Muster erkennen zu wollen und zur Abschätzung einer sinnvollen Stichprobengröße.

Personen, die nicht mit Statistik vertraut sind, stellen häufig die Frage, warum man nicht einfach die Mittelwerte der jeweiligen Datensätze berechnen und vergleichen kann, um einen statistischen Unterschied festzustellen?

Nun, das kann man durchaus machen, wenn sich die Mittelwerte um einen Faktor 2 oder größer unterscheiden und man Grund zur Annahme hat, dass die Daten nicht stark streuen. Für solche Fälle braucht man keinen Signifikanztest. Allerdings geht es häufig nicht um offensichtliche, sondern um geringe Unterschiede, sagen wir einmal 10%. Hier wird es schwierig. Es gilt die Regel, um so kleiner der Unterschied, den man signifikant feststellen möchte und umso größer die Streuung der Werte, desto mehr Daten braucht man.

Üblicherweise muss man aber bei einem Sigifikanztest wie dem Student-t-Test davon ausgehen, dass man die theoretisch exakte Verteilung der Daten (in diesem Fall eine Normalverteilung) irgendwie schon kennt. Es muss zwar nicht unbedingt eine Normalverteilung sein, aber die Verteilung muss durch eine mathematische Funktion mit Parametern beschrieben werden  können. Es dürfen also nur die Parameter unbekannt sein, nicht aber die Art der Verteilung. Wikipedia formuliert es so.

Überprüft wird statistische Signifikanz durch statistische Tests, die so gewählt werden müssen, dass sie dem Datenmaterial und den zu testenden Parametern bezüglich der Wahrscheinlichkeitsfunktion entsprechen. Nur dann ist es möglich, aus der Wahrscheinlichkeitsverteilung für Zufallsvariablen mathematisch korrekt den jeweiligen p-Wert zu errechnen als die Wahrscheinlichkeit, ein Stichprobenergebnis wie das beobachtete oder ein extremeres zu erhalten.

Ich will zeigen, wie die Voraussetzung einer parametrischen Funktion umgangen werden kann. Das löst natürlich nicht das grundsätzliche Problem, das Signifikanz nicht mit der Wahrscheinlichkeit einer Tatsache gleichzusetzen ist. Ich werde in dem Zuge allerdings noch einmal anhand einer visuellen Aufbereitung der Daten klären, worin das Problem mit der Signifikanz liegt.

Es sollte uns doch egal sein können, wie die Daten verteilt sind. Sie sind eben so verteilt, wie sie verteilt sind und nicht so, wie wir sie gerne hätten, dass sie theoretisch in der Grundgesamtheit mathematisch schön verteilt wären. Wieso sollte man überhaupt zusätzlich irgendwelche hypothetischen Annahmen über die Art der Verteilung treffen müssen? Die ganze Information, die zur Verfügung steht, steckt doch bereits zu 100% in den verfügbaren Daten. Mehr als diese Information hat man nicht.

Was macht man, wenn man feststellt, dass man alle statistischen Verteil-Funktionen und Tests im Lehrbuch durchprobiert und keine so recht passen mag? Es reicht ja schon, wenn eine doppelte Glockenkurve vorliegt oder die Form der Glocke willkürlich schief verzerrt ist? Kann man in solchen Fällen dann nicht herausfinden, ob eine signifikante Abweichung vorliegt? Oder ist es ok, wenn man trotzdem von einer Normalverteilung ausgeht? Muss man dann einen zusätzlichen Fehler in Kauf nehmen oder sogar berechnen?

Es gibt also gleich zwei Probleme bei der Durchführung eines Signifikanztest. Erstens, man erhält ganz grundsätzlich nur ein Indiz statt eines Beweises. Zweitens, man muss sicherstellen, dass sich die Daten mit einer parametrischen Wahrscheinlichkeitsfunktion fitten lassen und den entsprechenden Test wählen, sonst erhält man irgendeinen nicht zu interpretierenden Wert.

Im folgenden will ich aufzeigen, dass die Aussage der Wikipedia falsch ist. Man kann eine signifikante Abweichung mit einem allgemeinen Test ohne die Annahme einer zugrunde liegenden parametrischen Verteilung feststellen. Einziger Wermutstropfen des Verfahrens ist, dass man auf den Computer zurückgreifen muss, um sehr viele Zufallszahlen zu erzeugen und zu sortieren. Das hat zur Folge, dass man die Signifikanz nicht mehr mit der Hand oder per Excel anhand einer einfachen Formel bestimmen kann. Andererseits sind die Daten heute sowieso im Computer, also warum ihn dann nicht auch für weiterführende Berechnungen nutzen? Der Programmcode dafür umfasst in Python ganze 11 Zeilen (siehe unten) und sollte sich in jeder Programmiersprache, selbst in Excel VBA realisieren lassen.

Man kann von einem nichtparametrischen Verfahren sprechen, da die Berechnung nicht auf den Parametern einer statistischen Verteilfunktionen basiert, sondern auf den Daten selbst operiert. Dieses Verfahren funktioniert deshalb für beliebige Verteilungen. Es ist somit egal, was für eine Verteilung vorliegt.

Der Trick besteht in der Verwendung eines Zufallszahlengenerators.

Ein Verfahren, bei dem man unterschiedliche Variationen per Zufallszahlen berechnet, nennt man Monte-Carlo Simulation. Diese Verfahren wurde ausnahmsweise nicht von einem Graf von Monte-Carlo siebzehnhundert irgendwann erfunden, sondern ist ein Codename aus der amerikanischen Atombomben Entwicklung. Mit dem Verfahren kann man Integrationsprobleme lösen, zu denen man keine Stammfunktion bilden kann – die sogenannte Monte-Carlo-Integration.

Im Allgemeinen sprechen wir vom Monte-Carlo-Verfahren wenn wir einen Zufallszahlen-Generator verwenden, um einzelne Ergebnisse zu berechnen und dann aus der Summe der vielen Einzelergebnisse wiederum Informationen gewinnen können. Man kann die Strategie vielleicht treffender als Massiv-Random-Sampling (massive Zufalls-Probe) bezeichnen. Man erkauft sich die Einfachheit des Verfahrens durch pure Rechenpower, denn die Lösung wird um so genauer, desto mehr Zufallsfälle hinzugezogen werden. Die Lösung konvergiert irgendwo zwischen 1 Mio und unendlich vielen Zufallszahlen zur exakten Lösung, für unsere Beispiel bereits ab ca. 10.000. Durch die Konvergenz kann man für praktische Berechnungen solange Zufallszahlen hinzuzählen, bis eine gewünschte Güte erreicht ist.

Der Trick besteht allerdings nicht in der Erzeugung irgendwelcher Zufallszahlen, sondern Zahlenwerten, die sich nahtlos in die vorliegenden Daten einfügen, die für sich genommen also die gleiche Glockenkurve ergeben wie die eigentlichen Daten.

Schauen wir uns einmal folgende Datensätze zweier Produktionszellen (aus einer Kaggle Challenge) in der üblichen Glockenkurven-Darstellung an.

1_24_802_dist 0_17_431_dist

set_0 (schwarz) repräsentiert jeweils die Daten einer Produktionszelle ohne Ausschuss. set_1 (rot) repräsentiert die Daten der jeweiligen Anlage mit Anschluss. Der Ausschuss wurde erst am Ende der Produktionskette von mehreren hundert Fertigungsschritten festgestellt. Kann man in den Daten einen signifikanten Unterschied feststellen, der darauf hindeutet, dass diese beiden Produktionszellen eine erhöhte Wahrscheinlichkeit für einen Ausschuss bedingen – also im Endeffekt irgendetwas damit zu tun haben könnten?

Mit einer Glockenkurvendarstellung kann man allerdings wenig anfangen, da ein Histogramm sehr stark von der Klassifizierung abhängt. Wir sortieren deshalb die Werte der Größe nach, weisen jedem Wert einen Anteil zu und tragen sie dann der Reihe nach auf. Wir erhalten dadurch die kumulierte Wahrscheinlichkeit.

1_24_802_cumdist 0_17_431_cumdist

Bei 1. Mio Daten wird dem kleinsten Wert die Wahrscheinlichkeit ein Millionstel, dem zweit kleinsten Wert die Wahrscheinlichkeit zwei Millionstel zugewiesen usw. (Es ist nicht genau ein Millionstel, da die Wahrscheinlichkeit nie von 0 bis 1 definiert ist, sondern für 99 Werte z.B nur von 0.01 – 0.99).

Anschließend verzerren wir die Darstellung der y-Achse noch so, dass Werte, die perfekt normalverteilt wären, auf einer Geraden in dem Diagramm liegen. Ein solches Diagramm nennt sich Wahrscheinlichkeits-Diagramm.

1_24_802_prob 0_17_431_prob

Die Vorteile eines Wahrscheinlichkeitsdiagramms liegen auf der Hand.

  • Man sieht sofort, wie sehr die Werte auf einer Geraden liegen und dementsprechend normalverteilt sind oder eben nicht.
  • Es sind alle Werte numerisch exakt dargestellt und müssen nicht wie in einem Histogramm zunächst klassifiziert werden.
  • Die Steigung der Geraden entspricht der Standardabweichung, also dem Maß für die Streuung der Ergebnisse.
  • Der Median und Quantile können entsprechend den  Probability-Werten direkt abgelesen werden.
  • Man kann übersichtlich mehrere Verteilungen in einem Diagramm darstellen.

Mit etwas Eingewöhnung ist man nie wieder auf ein Histogramm angewiesen.

Wie man an den Daten sieht, sind in beiden Fällen die Daten nicht normalverteilt. Den ersten Datensatz könnte man vielleicht mit Augen zu noch als annähernd normalverteilt betrachten, den zweiten allerdings mit seinen einzelnen Häufungen keineswegs.

Es ist zwar ein Unterschied deutlich sichtbar, aber könnte es nicht sein, dass set_1 nur durch Zufall von set_0 abweicht, also vom Prinzip derselben Verteilung wie set_0 entspricht?

An den Diagrammen verdeutlicht sich visuell das Problem des Signifikanztests. Signifikanz könnte alleine dadurch entstehen, dass set_0 zufällig etwas mehr nach links und set_1 zufällig etwas nach rechts driftet. Dieser Drift könnte bei weiteren (zufällig anderen) Werten verschwinden und liegt an der einmaligen Stichprobennahme.

Die Idee besteht nun darin, den zufallsbedingten Drift der beiden Datensätze zu schätzen (also zu berechnen) und zu prüfen, wie stark sich durch einen solchen Drift die Kurven annähern könnten.

Aus den Datensätzen werden dazu Konfidenzbänder mit zwei Hüllkurven berechnet. Die Kurve hätte dann zufällig auch irgendwie anders in dem Konfidenzband verlaufen können, nur nicht außerhalb des Bandes, zumindest nicht sehr wahrscheinlich. Und über genau diese Wahrscheinlichkeit, die man zulassen möchte, definiert man das Konfidenzband. Man bekommt also ein gutes visuelles Gefühl dafür, wie groß der Drift ausfallen könnte, der bereits in den Daten steckt. Denn – um es noch einmal zu wiederholen – es ist ja nicht davon auszugehen, dass die Kurve im Mittel bereits richtig liegt. Es ist sogar sehr unwahrscheinlich, dass der Kurvenverlauf bei Erhöhung der Daten so bleibt. Wegen dieser Einschränkung sprechen wir von Signifikanz und nicht von der Wahrscheinlichkeit eines statistischen Unterschieds.

Allerdings steigt die Wahrscheinlichkeit, dass die Daten im Mittel stimmen, mit zunehmender Datenanzahl bzw. Stichprobengröße. Aber auch hier ist Vorsicht geboten. Es könnte ja sein, dass es einen Drift über die Zeit der Datenerhebung gibt. Wie gesagt, man darf nur in der juristischen Dimension eines Indizes über Signifikanz denken. Eine Gegenpartei könnte z.B. den zeitlichen Drift beweisen oder zumindest ein ebenso starkes Indiz dafür vorliegen. Das hat auf die ermittelte Signifikanz keinen Einfluss, doch der Widerspruch mindert die Bedeutung des Unterschieds und damit auch, dass es wahr sein könnte – denn Wahrheit ist immer widerspruchsfrei.

Mit ca. 500.000 Datenpunkten für set_0 kommt in diesen Beispielen die Stichprobe allerdings einer Grundgesamtheit schon recht nahe. Das Konfidenzband fällt hier bis auf Strichstärke mit der Kurve zusammen. Wir können dadurch (in diesem Beispiel) vernachlässigen, dass set_0 einen zufälligen Drift aufweist. Ist das nicht der Fall, könnte man z.b. die Konfidenzbänder beider Kurven berechnen und deren Überlappung begutachten. Das ist allerdings durchaus komplizierter.

Zur Berechnung der Konfidenzbänder verwenden wir einen Zufallszahlengenerator, der Datenpunkte geniert, die der Verteilung von set_0 entsprechen. Dazu wird ein herkömmlicher Zufallszahlengenerator verwendet, der einen Wert zwischen 0 und 1 generiert. Mit diesem Wert geht man im kumulierten Diagramm nun von der y-Achse aus bis zur schwarzen Linie und liest den zugehörigen Wert auf der x-Achse ab. Dazu dient ein herkömmlicher Interpolations-Algorithmus, am besten ein Spline und mit Extrapolationsbedingungen.

Es werden so viele Zufallszahlen auf diese Weise erzeugt, wie set_1 Datenpunkte aufweist. Dieses neue Set nennen wir set_R. Wir erstellen nun noch viele weitere Male ein zufälliges set_R und sortieren in jedem set_R die Werte der Größe nach. Man stelle sich dazu zu jedem roten Punkt in den folgenden Diagrammen 10.000 grüne Punkte vor, die rechts und links daneben auf einer Horizontalen verlaufen. Die Quantile aller Datenpunkte bilden die beiseitigen Konfidenzbändern. Für einen Konfidenzbereich von 0,95 werden z.B. die 0,025 und 0,975 Quantile ermittelt, also bei 10.000 Werten jeweils 250 links und rechts abgeschnitten. Mit 95% Wahrscheinlichkeit liegt eine erneute Stichprobe im Umfang von set_1  in diesem Bereich.

 

1_24_802_probconf0_17_431_probconf

Ein Drift innerhalb des Konfidenzbands ist folglich alleine durch die Größe der Stichprobe möglich. Man kann durch die Konfidenzbänder rein visuell prüfen, ob der Kurvenverlauf von set_1 innerhalb des Konfidenzbands von set_0 verläuft – für die gleiche Stichprobengröße von set_1 wohlgemerkt, denn das Konfidenzband von set_0 (mit der Stichprobengöße von set_0) liegt auf Strichstärke.

Weiterhin kann man nun einen ooc-Wert (out of confidence) als Kennzahl für die Signifikanz angeben, der quantifiziert, wie groß der Anteil der Werte außerhalb des Konfidenzbandes ist, und einen dev-Wert (Deviation) bestimmen, der angibt, wie groß die relative Abweichung der Kurven zueinander in x-Richtung ist.

Sie liefern für die dargestellten Beispiel mit 100.000 Samples und 95% Konfidenzbereich

ooc        /   dev

57,7%  /  3,8%

79,2%  / 3,3%

In beiden Fällen ist somit eine Abweichung  zwischen 3% und 4% mit einer Signifikanz von 60% bzw 80% ooc festzustellen. Was sagt uns das nun?

Es macht aus der visuellen Vermutung, dass es einen Unterschied gibt, eine signifikante Bestätigung. Mehr nicht. Ob es den Einfluss tatsächlich gibt, lässt sich mit diesem Verfahren, wie mit jedem Signifikanztest, nicht sagen. Aber das Verfahren ermöglicht einfache Kennzahlen, mit denen automatisch große Datenmengen bewertet werden können, ohne sich Gedanken um einen geeigneten Test zu machen. Obendrein lassen sich die Kennzahlen visuell interpretieren.

Der ooc-Wert gibt den Anteil der Datenpunkte an, die außerhalb des Konfidenzbandes gleicher Stichprobe liegen. Der dev-Wert gibt den relativen Unterschied an und zwar nicht als Unterschied der Mittelwerte, sondern als über die Länge der Kurve und über die Spanne der Werte normierter Betrag des Abstandes der Kurven in x-Richtung. Es folgt der Programmcode für die Berechnung der beiden Werte in Phyton und Matlab.


Pyhton Code mit Kommentaren

# Einbindung der Bibliotheken
import numpy as np
import numpy.random as rnd

# Funktionsdefinition
def nonparstat(set_0, set_1, confidencelevel=0.95, samples = 10000):
# Nichtparametrische Berechnung der statistischen Signifikanz
# [ ooc, dev ] = nonparstat(set_0,set_1,confidencelevel,samples)
#
#   Ausgangsgrößen
#   ooc:                Relative Anteil der signifikanten Abweichung der
#                       Daten bezogen auf die Gesamtmenge der Daten in set_1.
#                       Wert zwischen 0 und 1 für 0% - 100% Signifikanz.
#                       Die signifikante Abweichung kann somit auch nur für bestimmte
#                       Wertebereiche vorliegen, in anderen aber nicht.
#                       Beträgt ooc z.B. 0.1, weichen 10% der Daten
#                       signifikant ab.
#   dev:                Relative Abweichung bezogen auf den gesamten Wertebereich.
#                       Beträgt dev z.B. 0.1, weicht set_1 um 10% von set_0 ab.
#
#   Eingangsgrößen
#   set_0:              Datensatz, gegen den geprüft werden soll.
#   set_1:              Datensatz, dessen signifikante Abweichung gegenüber set_0 zu 
#                       berechnen ist.
#   confidencelevel:    Der Vertrauensbereich. Üblichweise 0.8,0.9,0.95,0.99, 
#                       Default 0.95, falls nicht spezifiziert.
#                       Der Vertrauensbereich berechnet sich aus der zulässigen 
#                       Irrtumswahrscheinlichkeit 
#                       (0.99 bedeutet z.B., dass mit 1-0.99=0.01=1% Wahrscheinlichkeit
#                       das Ergebnis doch nicht signifkant ist.)
#   samples:            Anzahl der Zufallszahlen zur Berechnung des Vertrauensbereichs. 
#                       Desto größer der Vertrauensbereich gewählt ist, desto mehr 
#                       Samples müssen verwendet werden. 
#                       Default 10000, falls nicht spezifiziert

# Erstellung eines Wahrscheinlichkeitsvektors passend zur Länge von set_0
prob_0 = np.linspace(1/len(set_0), 1-1/(len(set_0)),len(set_0))

# Sortieren der Datensets der Werte nach.
set_0_sort = np.sort(set_0, axis=None)
set_1_sort = np.sort(set_1, axis=None)

# Berechnung der Zufallszahlen aus der Verteilung von prob und set_0 per Interpolation 
# und Sortierung
set_R_sort = np.sort(np.interp(rnd.random((samples,len(set_1))),prob_0,set_0_sort,left=set_0_sort[1], right=set_0_sort[-1]), axis=1)

# Ermittlung der Quantile aus den sortierten Reihen von set_R
set_R_sort_mean = np.nanpercentile(set_R_sort,50,axis=0)
set_R_sort_qmin = np.nanpercentile(set_R_sort,(1-confidencelevel)/2*100,axis=0)
set_R_sort_qmax = np.nanpercentile(set_R_sort,(confidencelevel+(1-confidencelevel)/2)*100,axis=0)

# Berechnung der Wertes aus set_1, die ausserhalb des Konfidenzintervalls 
# (out of confidence) liegen.
ooc = np.sum(np.logical_or(set_1_sort < set_R_sort_qmin , set_1_sort > set_R_sort_qmax))/len(set_1)

# Berechnung der mittleren, relative Abweichung bezogen auf den gesamten Wertebereich 
# berechnet.
dev = np.double(np.sum(np.abs(set_1_sort-set_R_sort_mean)) /np.abs(max(set_0)-min(set_0))/len(set_1))

return (ooc , dev)

Matlab Code ohne Kommentare

function [ ooc, dev ] = nonparstat(set_0,set_1,confidencelevel,samples)

if nargin < 4
    samples=10000;
end
if nargin < 3
    confidencelevel=0.95;
end

prob=[1/length(set_0) : 1/(length(set_0)) : 1-1/(length(set_0)),1]';

set_0_sort=sort(set_0);
set_1_sort=sort(set_1);
set_R_sort=sort(reshape(interp1(prob,set_0_sort,random('uniform',0,1,samples*length(set_1),1),'spline'),[length(set_1),samples]));
    
set_R_sort_mean=median(set_R_sort')';
set_R_sort_qmin=quantile(set_R_sort',(1-confidencelevel)/2)';
set_R_sort_qmax=quantile(set_R_sort',confidencelevel+(1-confidencelevel)/2)';

ooc = sum(or((set_1_sort < set_R_sort_qmin),(set_1_sort > set_R_sort_qmax)))/length(set_1);      
dev = sum(abs(set_1_sort-set_R_sort_mean)) /abs(max(set_0)-min(set_0))/length(set_1);

end

Anmerkung zum Code. Zur Ausgabe der Konfidenzbänder muss der Code lediglich so umgewandelt werden, dass set_R_sort_qmin und set_R_sort_qmax als Vektorarrays ausgegeben werden.

Ein einfacher nichtparametrischer Signifikanztest für beliebige Verteilungen
5 (100%) 1 vote

Schreibe einen Kommentar

Deine E-Mail-Adresse wird nicht veröffentlicht. Erforderliche Felder sind mit * markiert.