"There are no routine statistical questions, only questionable statistical routines." - Sir David R. Cox

“Viel besser eine ungefähre Antwort auf die richtige Frage, die oft vage ist, als eine exakte Antwort auf eine falsche Frage, die immer genauer formuliert werden kann.” — John Tukey

Wir verwenden die explorative Datenanalyse um, das Verständnis unserer Daten zu entwickeln. Die einfachste Methode ist Fragen als Werkzeug der Untersuchung zu verwenden. Wenn wir eine Frage stellen, konzentrieren wir uns auf einem bestimmten Datenteil. Dies hilft uns zu entscheiden, welche Charts, Modelle oder Transformationen wir verwenden sollten.

EDA ist ein kreativer Prozess. Wie in den meisten kreativen Prozessen ist der Schlüssel zu einer guten Frage, noch mehr Fragen zu stellen. Es ist schwierig solche Fragen am Anfang der Analyse zu stellen, weil wir noch nicht wissen, welche Schlussfolgerungen im Datensatz enthalten sind. Von der anderen Seite hebt jede neue Frage einen anderen Aspekt hervor und erhöht die Chance, eine Entdeckung zu machen. Wir können schnell zu dem interessantesten Teil der Stichprobe übergehen und durch konsequente Fragen die Situation klären.

Es gibt keine Regeln, welche Fragen wir stellen müssen, um die Untersuchung durchzuführen. Aber zwei Typen von Fragen werden immer hilfreich sein:

• Welcher Typ haben die Änderungen in meinen Variablen?
  
• Welchen Typ haben die Kovariationen zwischen den Variablen?

Definieren wir den Hauptbegriff.

Variationen stellen die Variabilität der Werte einer Variablen von einer Dimension in die andere dar. Es gibt eine Vielzahl an Variationen im alltäglichen Leben. Wenn man eine kontinuierliche Variable sieben Mal hintereinander misst, erhält man sieben verschieden Ergebnisse. Das gilt auch wenn Sie konstante Größen messen (z.B. Lichtgeschwindigkeit). Jede der Dimensionen wird einen Fehler enthalten, der sich immer ändert. Die Variablen eines Typs können sich auch ändern, wenn man verschiedene konkrete Gegenstände (z.B. die Augenfarbe verschiedener Menschen) oder ein Gegenstand zu verschiedenen Zeitpunkten misst (z.B. das Energieniveau eines Elektrons zu unterschiedlichen Zeitpunkten). Jede Variable hat einen eigenen Charakter der Variationen, die eine interessante Information aufdecken können. Die beste Weise, diese Information zu verstehen, ist die Werteverteilung der Variablen zu visualisieren. Das ist gerade der Fall, wenn ein Bild mehr als tausend Worte sagt.

2.1. Gesamtstatistik

Die Gesamtstatistik unserer Zeitreihe kann man mithilfe der Funktion table.Stats()::PerformenceAnalitics verfolgen.

> table.Stats(env$dataSet %>% tk_xts())
Using column `Data` for date_var.
                     ftlm      stlm      rbci      pcci
Observations    7955.0000 7908.0000 7934.0000 7960.0000
NAs               42.0000   89.0000   63.0000   37.0000
Minimum           -0.7597   -1.0213   -0.9523   -0.5517
Quartile 1        -0.0556   -0.1602   -0.0636   -0.0245
Median            -0.0001    0.0062   -0.0016   -0.0001
Arithmetic Mean    0.0007    0.0025    0.0007    0.0001
Geometric Mean    -0.0062       NaN   -0.0084   -0.0011
Quartile 3         0.0562    0.1539    0.0675    0.0241
Maximum            2.7505    3.0407    2.3872    1.8859
SE Mean            0.0014    0.0033    0.0015    0.0006
LCL Mean (0.95)   -0.0020   -0.0040   -0.0022   -0.0010
UCL Mean (0.95)    0.0034    0.0090    0.0035    0.0012
Variance           0.0152    0.0858    0.0172    0.0026
Stdev              0.1231    0.2929    0.1311    0.0506
Skewness           4.2129    1.7842    2.3037    6.4718
Kurtosis          84.6116   16.7471   45.0133  247.4208
                   v.fatl    v.satl    v.rftl    v.rstl
Observations    7959.0000 7933.0000 7954.0000 7907.0000
NAs               38.0000   64.0000   43.0000   90.0000
Minimum           -0.3967   -0.0871   -0.1882   -0.4719
Quartile 1        -0.0225   -0.0111   -0.0142   -0.0759
Median            -0.0006    0.0003    0.0000    0.0024
Arithmetic Mean    0.0002    0.0002    0.0002    0.0011
Geometric Mean    -0.0009    0.0000   -0.0003   -0.0078
Quartile 3         0.0220    0.0110    0.0138    0.0751
Maximum            1.4832    0.3579    0.6513    1.3093
SE Mean            0.0005    0.0002    0.0003    0.0015
LCL Mean (0.95)   -0.0009   -0.0003   -0.0005   -0.0020
UCL Mean (0.95)    0.0012    0.0007    0.0009    0.0041
Variance           0.0023    0.0005    0.0009    0.0188
Stdev              0.0483    0.0219    0.0308    0.1372
Skewness           5.2643    2.6705    3.9472    1.5682
Kurtosis         145.8441   36.9378   74.4182   13.5724
                   v.ftlm    v.stlm    v.rbci    v.pcci
Observations    7954.0000 7907.0000 7933.0000 7959.0000
NAs               43.0000   90.0000   64.0000   38.0000
Minimum           -0.9500   -0.2055   -0.6361   -1.4732
Quartile 1        -0.0280   -0.0136   -0.0209   -0.0277
Median            -0.0002   -0.0001   -0.0004   -0.0002
Arithmetic Mean    0.0000    0.0001    0.0000    0.0000
Geometric Mean    -0.0018   -0.0003   -0.0009       NaN
Quartile 3         0.0273    0.0143    0.0207    0.0278
Maximum            1.4536    0.3852    1.1254    1.9978
SE Mean            0.0006    0.0003    0.0005    0.0006
LCL Mean (0.95)   -0.0012   -0.0005   -0.0009   -0.0013
UCL Mean (0.95)    0.0013    0.0007    0.0009    0.0013
Variance           0.0032    0.0007    0.0018    0.0034
Stdev              0.0561    0.0264    0.0427    0.0579
Skewness           1.2051    0.8513    2.0643    3.0207
Kurtosis          86.2425   23.0651   86.3768  233.1964

Was wir anhand dieser Tabelle sehen wir:

• Alle Prädiktoren haben eine (relativ) kleine Anzahl undefinierter Variablen NA.
• Alle Prädiktoren haben eine stark ausgeprägte Rechtsschiefe (skewness).
• Alle Prädiktoren haben eine hohe Kurtosis (kurtosis).

2.2. Visualisierung der Gesamtstatistik

“Der größte Wert eines Bildes besteht darin, dass es uns darauf aufmerksam macht, was wir nie erwartet haben zu sehen”. — John Tukey

Schauen wir uns die Variation und die Kovariation unserer Variablen im Datensatz dataSet an. Da die Anzahl der Variablen (14) es nicht erlaubt, diese in einem Chart zu platzieren, teilen wir sie in drei Gruppen auf.

require(GGally)
evalq(ggpairs(dataSet, columns = 2:6, 
              mapping = aes(color = Class),
              title = "DigFilter1"), 
      env)

Abb. 6. Die erste Gruppe der Prädiktoren

evalq(ggpairs(dataSet, columns = 7:10, 
              mapping = aes(color = Class),
              title = "DigFilter2"), 
      env)

Abb. 7. Die zweite Gruppe der Prädiktoren

evalq(ggpairs(dataSet, columns = 11:14, 
              mapping = aes(color = Class),
              title = "DigFilter3"), 
      env)

Abb. 8. Die dritte Gruppe der Prädiktoren

Was wir in den Charts sehen:

• alle Prädiktoren haben eine Verteilungsform, die nahe an der normalen ist, aber mit einer stark ausgeprägten Rechtsschiefe;
• alle Prädikatoren haben einen kleinen Interquartilsabstand (IQR);
• alle Prädiktoren haben erhebliche Ausreißer;
• die Anzahl der Beispiele auf den zwei Ebenen der Zielvariablen "Class" unterscheidet sich ein wenig von einander.

3. Datenaufbereitung

3.1. Bereinigen von Daten

Die erste Phase der Aufbereitung von Ausgangsdaten stellt das Entfernen oder die Imputation undefinierter Werte und Lücken in den Daten dar. Obwohl viele Modelle die Verwendung undefinierter Daten (NA) und fehlender Daten in Datensätzen zulassen, ist es besser, diese vor den Manipulationen zu entfernen. Diese Operation wird für den vollen Datensatz unabhängig vom verwendeten Modell durchgeführt.

Die Gesamtstatistik für unsere Rohdaten weist darauf hin, dass NA im Datensatz vorhanden sind. Das sind künstliche A, die sich bei der Berechnung digitaler Filter gebildet haben. Ihre Anzahl ist relativ klein, deswegen kann man sie entfernen. Wir haben bereits den Datensatz dataSet erhalten, der für die weitere Verarbeitung bereit ist. Bereinigen wir die Daten.

Im Allgemeinen versteht man unter einer "Bereinigung" folgendes:

• Entfernen von Prädiktoren mit einer Null- oder Nahe-Null-Varianz (method = c(“zv”, “nzv”));
• Entfernen stark korrelierter Variablen. Der Grenzwert der Korrelation wird von Nutzer gesetzt (method = “corr”), standardmäßig ist er gleich 0.9. Diese Phase ist nicht immer nötig, das hängt von den weiteren Methoden der Umwandlung ab;
• Entfernen von Prädiktoren, die nur einen einmaligen Wert in einer beliebigen Klasse haben (method = “conditionalX”).

Alle diese Operationen sind in der Funktion preProcess()::caret durch die oben angeführten Methoden implementiert und werden für den vollen Datensatz vor der Aufteilung in Einlern- und Testdatensätze durchgeführt.

require(caret)
evalq({preProClean <- preProcess(x = dataSet,method = c("zv", "nzv", "conditionalX", "corr"))
      dataSetClean <- predict(preProClean, dataSet %>% na.omit)},
env)

Schauen wir, ob es entfernte Prädikatoren gibt, und was nach dem Bereinigen übrig geblieben ist:

> env$preProClean$method$remove
#[1] "v.rbci" 
> dim(env$dataSetClean)
[1] 7906   13
> colnames(env$dataSetClean)
 [1] "Data"   "ftlm"   "stlm"   "rbci"   "pcci"  
 [6] "v.fatl" "v.satl" "v.rftl" "v.rstl" "v.ftlm"
[11] "v.stlm" "v.pcci" "Class"

3.2. Erkennung und Verarbeitung von Ausreißern (outliers)

Probleme mit der Datenqualität, solche wie Schiefe oder Ausreißer, hängen oft voneinander ab. Dies macht eine vorläufige Verarbeitung aufwändig und erschwert die Aufgabe, Zusammenhänge und Tendenzen in der Stichprobe festzustellen.

Was ist ein Ausreißer?

Definieren wir den Begriff: unter einem Ausreißer verstehen wir eine Beobachtung, die zu groß oder zu klein im Vergleich zu den meisten anderen Beobachtungen ist. Eine gründliche Klassifikation von Ausreißern sowie Methoden deren Erkennung und Verarbeitung ist in [2] angeführt.

Typen der Ausreißer

Ausreißer verzerren die Verteilung von Variablen und die anhand solcher Daten trainierten Modelle. Es gibt eine Vielzahl von Methoden für die Erkennung und Verarbeitung von Ausreißern. Sie hängen davon ab, ob wir den Ausreißer global oder lokal identifizieren. Lokale Ausreißer stellen Ausreißer einer Variablen dar. Globale Ausreißer sind Ausreißer in einem mehrdimensionalen Raum, der durch eine Matrix oder einen Datenrahmen definiert wird.

Was verursacht Ausreißer?

Alle Ausreißer kann man nach ihrem Ursprung aufteilen:

Künstliche

• Eingabefehler — Fehler, die während des Sammelns, Schreibens oder der Verarbeitung von Daten auftreten;
• experimentelle Fehler;
• Stichprobenfehler.

Natürliche — Ausreißer, die durch die Natur der Variablen bedingt sind.

Welche Auswirkungen haben Ausreißer?

Die Ausreißer können die Ergebnisse der Datenanalyse oder einer statistischen Modellierung erheblich beeinträchtigen: sie erhöhen die Fehlervarianz und reduzieren die Leistung der statistischen Tests; wenn die Ausreißer nicht zufällig verteilt wurden, können sie auch die Normalität reduzieren; darüber hinaus können sie die grundlegenden Annahmen der Regressionsanalyse, Varianzanalyse und anderer statistischer Annahmen des Modells beeinflussen.

Wie kann man lokale Ausreißer erkennen?

Am häufigsten werden die Ausreißer durch das Visualisieren von Daten ermittelt. Eine einfache und breit verwendete Methode ist boxplot. Schauen wir uns zum Beispiel den Prädiktor ftlm an:

evalq(ggplot(dataSetClean, aes(x = factor(0), 
                               y = ftlm,
                               color = 'red')) + 
        geom_boxplot() + xlab("") + 
        scale_x_discrete(breaks = NULL) + 
        coord_flip(),
      env)

Abb.9. Boxplot ftlm

Erläuterungen zum Bild:

IQR — Interquartilsabstand, oder der Abstand zwischen dem ersten und dritten Quartil.

In diesem Zusammenhang gibt es mehrere Definitionen von Ausreißern:

• Jeder Wert, der kleiner als -1.5*IQR und größer +1.5*IQR ist, ist ein Ausreißer. Manchmal wird der Koeffizient auf 2 oder 3 gesetzt. Alle Werte zwischen 1.5*IQR und 3*IQR werden als "milde" Ausreißer bezeichnet, alle Werte, die oberhalb von 3*IQR liegen, sind extreme Ausreißer.
• Jeder Wert außerhalb des Bereichs vom 5. und 95. Perzentil kann als Ausreißer bezeichnet werden.
• Punkte mit einem Abstand von drei Standardabweichungen und mehr sind auch Ausreißer.

Weiterhin werden wir die erste Definition des Ausreißers verwenden (über IQR).

Wie werden Ausreißer verarbeitet?

Die meisten Methoden für die Verarbeitung von Ausreißern sind ähnlich den Methoden für die Verarbeitung undefinierter Werte NA: Entfernen, Umwandeln, Aufteilen, Imputation und andere statistische Methoden.

• Entfernen von Ausreißern: wir entfernen die Werte von Ausreißern, wenn sie infolge eines Eingabefehlers oder eines Verarbeitungsfehlers auftreten oder wenn es zu wenige Ausreißer gibt. Darüber hinaus können wir die Verteilung an beiden Enden schneiden, um Ausreißer zu entfernen. Zum Beispiel, jeweils 1% von oben und von unten zu entfernen.


• Transformation und Binning (Gruppieren von Werten in Intervalle):
  
  • eine Transformation der Variablen kann Ausreißer ausschließen (dies wird im nächsten Abschnitt des Artikels betrachtet);
  • natürlicher Logarithmus reduziert die Änderungen, die die Extremwerte verursachen (ebenso im nächsten Abschnitt betrachtet);
  • Diskretisierung ist eine weitere Art, eine Variable umzuwandeln (s. den nächsten Abschnitt);
  • wir können auch die Zuweisung von Gewichten für verschiedene Beobachtungen verwenden (dieses Thema wird nicht betrachtet).



• Imputation: die gleichen Methoden, die wir für die Imputation undefinierter Werte verwenden, können auf Ausreißer angewandt werden. Man kann Mittelwert, Medianwert oder Modalwert verwenden. Vor der Imputation eines Wertes, müssen wir feststellen, ob wir mit einem natürlichen oder künstlichen Ausreißer zu tun haben. Wenn es um einen künstlichen Ausreißer geht, kann er imputiert werden.

Entfernen von Ausreißern

Wir entfernen die Werte von Ausreißern, wenn sie infolge eines Eingabefehlers oder eines Verarbeitungsfehlers aufgetreten sind oder wenn es zu wenige Ausreißer gibt (nur beim Bestimmen statistischer Messwerte der Variablen).

Die Daten einer Variablen (z.B. ftlm) können ohne Ausreißer wie folgt extrahiert werden:

evalq({dataSetClean$ftlm -> x  
  out.ftlm <- x[!x %in% boxplot.stats(x)$out]}, 
  env)

Oder so:

evalq({dataSetClean$ftlm -> x 
  out.ftlm1 <- x[x > quantile(x, .25) - 1.5*IQR(x) & 
          x < quantile(x, .75) + 1.5*IQR(x)]},
  env) 

> evalq(all.equal(out.ftlm, out.ftlm1), env)
[1] TRUE

> nrow(env$dataSetClean) - length(env$out.ftlm)
[1] 402 

boxplot(env$out.ftlm, main = "ftlm  without outliers", 
        boxwex = 0.5)

Die oben beschriebene Methode kann nicht für Matrizes und Dataframes verwendet werden, denn jede Variable in einem Dataframe kann verschiedene Anzahl von Ausreißern haben. Für solche Stichproben eignet sich die Methode der Ersetzung lokaler Ausreißer durch NA mit der weiteren Verwendung von Standardmethoden für die Verarbeitung von NA. Schreiben wir eine Funktion, die lokale Ausreißer durch NA ersetzt:

#-------remove_outliers-------------------------------
remove_outliers <- function(x, na.rm = TRUE, ...) {
  qnt <- quantile(x, probs = c(.25, .75), 
                  na.rm = na.rm, ...)
  H <- 1.5 * IQR(x, na.rm = na.rm)
  y <- x
  y[x < (qnt[1] - H)] <- NA
  y[x > (qnt[2] + H)] <- NA
  y
}

In unserem Datensatz dataSetClean ersetzen wir Ausreißer durch NA in allen Variablen außer c(Data, Class) und schauen, wie sich die Verteilung des erhaltenen Satzes x.out ändert:

evalq({
  dataSetClean %>% select(-c(Data,Class)) %>% as.data.frame() -> x 
  foreach(i = 1:ncol(x), .combine = "cbind") %do% {
    remove_outliers(x[ ,i])
  } -> x.out
  colnames(x.out) <- colnames(x)
  },  
env)
par(mfrow = c(1, 1))
chart.Boxplot(env$x, 
              main = "x.out with outliers",
              xlab = "")

Abb. 11. Daten mit Ausreißern

chart.Boxplot(env$x.out, 
              main = "x.out without outliers",
              xlab = "")

Abb.12. Daten ohne Ausreißer

Imputation von NA, die die Ausreißer ersetzten

Unter Imputation versteht man die Ersetzung fehlender, falscher oder ungültiger Werte durch andere Werte. Die Input-Daten für das Trainieren des Modells dürfen nur gültige Werte beinhalten. Mögliche Varianten:

• wir ersetzen NA durch mean, mediana, mod (dabei ändern sich die statistischen Charakteristiken des Datensatzes nicht)
• wir ersetzen die Ausreißer, die größer als 1.5*IQR sind, durch 0.95 Perzentil, und die Ausreißer, die kleiner als - 1.5*IQR sind, durch 0.05 Perzentil.

Schreiben wir eine Funktion, die die letzte Variante umsetzt, und schauen wir uns die Verteilung nach der Umwandlung an:

#-------capping_outliers-------------------------------
capping_outliers <- function(x, na.rm = TRUE, ...) {
  qnt <- quantile(x, probs = c(.25, .75), 
                  na.rm = na.rm, ...)
  caps <- quantile(x, probs = c(.05, .95), 
                   na.rm = na.rm, ...)
  H <- 1.5 * IQR(x, na.rm = na.rm)
  y <- x
  y[x < (qnt[1] - H)] <- caps[1] 
  y[x > (qnt[2] + H)] <- caps[2] 
  y
}

evalq({dataSetClean %>% select(-c(Data,Class)) %>%
    as.data.frame() -> x 
    foreach(i = 1:ncol(x), .combine = "cbind") %do% {
      capping_outliers(x[ ,i])
    } -> x.cap
    colnames(x.cap) <- colnames(x)
   },  
env)
chart.Boxplot(env$x.cap, 
              main = "x.cap with capping outliers",
              xlab = "")

Abb.13. Datensatz mit imputierten Ausreißern

Schauen wir auf die Variation und die Kovariation unserer Variablen im Datensatz dataSetCap (das ist der gleiche Datensatz dataSet, aber bereinigt und mit imputierten lokalen Ausreißern). Da die Anzahl der Variablen (13) es nicht erlaubt, diese in einem Chart zu platzieren, teilen wir diese in drei Gruppen auf.

evalq(x.cap %>% tbl_df() %>% 
        cbind(Data = dataSetClean$Data, .,
              Class = dataSetClean$Class) -> 
        dataSetCap, 
      env)
require(GGally)
evalq(ggpairs(dataSetCap, columns = 2:7, 

              mapping = aes(color = Class),
              title = "PredCap1"), 
      env)

Abb.14. Variation und Kovariation des ersten Teils des Datensatzes mit imputierten Ausreißern.

 Der zweite Teil des Datensatzes:

evalq(ggpairs(dataSetCap, columns = 8:13, 
              mapping = aes(color = Class),
              title = "PredCap2"), 
      env)

Wie kann man globale Ausreißer erkennen?

Zwei- und mehrdimensionale Ausreißer werden gewöhnlich mithilfe des Einfluss- oder Distanz-Index. Um globale Ausreißer festzustellen, werden verschiedene Abstände verwendet. Wir können dafür mehrere Pakete, z.B. DMwR, mvoutliers, Rlof, verwenden. Globale Ausreißer werden mithilfe von LOF (local outlier factor) bewertet. Berechnen und vergleichen wir LOF für den Datensatz mit Ausreißern x und für den Datensatz mit imputierten Ausreißern x.cap.

##------DMwR2-------------------
require(DMwR2)
evalq(lof.x <- lofactor(x,10), env)
evalq(lof.x.cap <- lofactor(x.cap,10), env)
par(mfrow = c(1, 3))
boxplot(env$lof.x, main = "lof.x", 
        boxwex = 0.5)
boxplot(env$lof.x.cap, main = "lof.x.cap", 
        boxwex = 0.5)
hist(env$lof.x.cap, breaks = 20)
par(mfrow = c(1, 1))

Abb.16. Globaler Ausreißerfaktor für einen Datensatz mit Ausreißern und einen Datensatz mit imputierten Ausreißern

Im Paket Rlof ist die Funktion lof() impelentiert. Sie findet den lokalen Faktor outlier [3] der Daten unter Verwendung von k Nachbarn. Lokaler Ausreißerfaktor (LOF) ist ein Wahrscheinlichkeitswert der Zugehörigkeit zu Ausreißern und wird für jede Beobachtung separat berechnet. Basierend auf diesem Wert trifft der Nutzer die Entscheidung, ob eine Beobachtung als Ausreißer betrachtet wird oder nicht.

Der LOF berücksichtigt die lokale Dichte, um festzustellen, ob eine Beobachtung ein Ausreißer ist. Das ist eine effizientere Implementierung des LOF unter Verwendung einer anderen Datenstruktur und Funktion zur Berechnung des Abstands im Vergleich zur Funktion lofactor(), die im Paket “dprep” verfügbar ist. Er unterstützt mehrere k Werte, die parallel berechnet werden, sowie verschiedene Abstandswerte, außer dem Euklidischen Abstand. Die Berechnungen werden parallel auf mehreren Kernen des Prozessors durchgeführt. Berechnen wir lofactor für dieselben Datensätze (x und x.cap) für 5, 6, 7, 8, 9 und 10 Nachbarn unter Verwendung der Minkowski Methode für die Berechnung des Abstands. Zeichnen wir die Histogramme dieser lofactor.

require(Rlof)
evalq(Rlof.x <- lof(x, c(5:10), cores = 2,
                       method = 'minkowski'),
        env)
  evalq(Rlof.x.cap <- lof(x.cap, c(5:10), 
                          cores = 2, 
                          method = 'minkowski'),
        env)
par(mfrow = c(2, 3))  
hist(env$Rlof.x.cap[ ,6], breaks = 20)
hist(env$Rlof.x.cap[ ,5], breaks = 20)
hist(env$Rlof.x.cap[ ,4], breaks = 20)
hist(env$Rlof.x.cap[ ,3], breaks = 20)
hist(env$Rlof.x.cap[ ,2], breaks = 20)
hist(env$Rlof.x.cap[ ,1], breaks = 20)
par(mfrow = c(1, 1))

Abb.17. Globaler Ausreißerfaktor für k Nachbarn

 Fast alle Punkte liegen im Bereich lofactor =1.6. Außerhalb des Bereichs:

> sum(env$Rlof.x.cap[ ,6] >= 1.6)
[1] 32

Die unbedeutende Anzahl der Ausreißer ist akzeptabel für diese Datensatzgröße.

Schreiben wir eine Funktion, die Zwischenberechnungen durchführt:

#-----prep.outlier--------------
prep.outlier <- function(x, na.rm = TRUE, ...) {
  qnt <- quantile(x, probs = c(.25, .75), 
                  na.rm = na.rm, ...)
  H <- 1.5 * IQR(x, na.rm = na.rm)
  caps <- quantile(x, probs = c(.05, .95), 
                   na.rm = na.rm, ...)
  list(lower = qnt[1] - H, upper = qnt[2] + H, 
       med = median(x), 
       cap1 = caps[1], cap2 = caps[2])
}

Berechnen wir die Parameter, die für die Erkennung und Imputation von Ausreißern notwendig sind. Legen wir die vorläufige Länge des Trainingsdatensatzes als 4000 Balken am Anfang fest und die weiteren 2000 Balken werden wir als Testdatensatz verwenden.

evalq(
  {train <- x[1:4000, ] 
  foreach(i = 1:ncol(train), .combine = "cbind") %do% {
    prep.outlier(train[ ,i]) %>% unlist()
  } -> pre.outl
  colnames(pre.outl) <- colnames(x)
  #pre.outl %<>% t()
  },  
  env)

Schauen wir uns das Ergebnis an:

> env$pre.outl
                   ftlm        stlm         rbci          pcci
lower.25% -0.2224942912 -0.59629203 -0.253231002 -9.902232e-02
upper.75%  0.2214486206  0.59242529  0.253529797  9.826936e-02
med       -0.0001534451  0.00282525 -0.001184966  8.417127e-05
cap1.5%   -0.1700418145 -0.40370452 -0.181326658 -6.892085e-02
cap2.95%   0.1676526431  0.39842675  0.183671973  6.853935e-02
                 v.fatl        v.satl        v.rftl        v.rstl
lower.25% -0.0900973332 -4.259328e-02 -0.0558921804 -0.2858430788
upper.75%  0.0888110249  4.178418e-02  0.0555115004  0.2889057397
med       -0.0008581219 -2.130064e-05 -0.0001707447 -0.0001721546
cap1.5%   -0.0658731640 -2.929586e-02 -0.0427927888 -0.1951978435
cap2.95%   0.0662353821  3.089833e-02  0.0411091859  0.1820803387
                 v.ftlm        v.stlm        v.pcci
lower.25% -0.1115823754 -5.366875e-02 -0.1115905239
upper.75%  0.1108670403  5.367466e-02  0.1119495436
med       -0.0003560178 -6.370034e-05 -0.0003173464
cap1.5%   -0.0765431363 -3.686945e-02 -0.0765950814
cap2.95%   0.0789209957  3.614423e-02  0.0770439553

Wie wir sehen, werden für jede Variable im Datensatz das erste und das dritte Quartil und der Medianwert sowie das 5. und 95. Perzentil definiert. Das ist alles, was man für die Erkennung und Verarbeitung von Ausreißern benötigt.

Wir brauchen eine Funktion, die die Ausreißer jeden Datensatzes nach vordefinierten Parametern verarbeiten wird. Mögliche Varianten der Verarbeitung: Ersetzen von Ausreißern durch NA,Ersetzen von Ausreißern durch den Medianwert, Ersetzen von Ausreißern durch 5./95. Perzentil.

#---------treatOutlier---------------------------------
  treatOutlier <- function(x, impute = TRUE, fill = FALSE,
                         lower, upper, med, cap1, cap2){ 
  if (impute) {
    x[x < lower] <- cap1 
    x[x > upper] <- cap2 
    return(x)
  }
  if (!fill) {
    x[x < lower | x > upper] <- NA 
    return(x)  
  } else {
    x[x < lower | x > upper] <- med
    return(x)
  }
} 

Da wir die notwendigen Parameter bereits im Trainigsdatensatz definiert haben, verarbeiten wir die Ausreißer im Trainingsdatensatz, indem wir diese durch das 5./95. Perzentil ersetzen. Danach verarbeiten wir auf die gleiche Weise Ausreißer im Testdatensatz und vergleichen die Verteilungen in den erhaltenen Datensätzen, indem wir drei Charts zeichnen.

#------------
evalq(
  {
  foreach(i = 1:ncol(train), .combine = 'cbind') %do% {
    stopifnot(exists("pre.outl", envir = env))
    lower = pre.outl['lower.25%', i] 
    upper = pre.outl['upper.75%', i]
    med = pre.outl['med', i]
    cap1 = pre.outl['cap1.5%', i] 
    cap2 = pre.outl['cap2.95%', i] 
    treatOutlier(x = train[ ,i], impute = T, fill = T, lower = lower,
                 upper = upper, med = med, cap1 = cap1, cap2 = cap2) 
  } -> train.out
  colnames(train.out) <- colnames(train)
  },
  env
) 
#-------------
evalq(
  {test <- x[4001:6000, ] 
  foreach(i = 1:ncol(test), .combine = 'cbind') %do% {
    stopifnot(exists("pre.outl", envir = env))
    lower = pre.outl['lower.25%', i] 
    upper = pre.outl['upper.75%', i]
    med = pre.outl['med', i]
    cap1 = pre.outl['cap1.5%', i] 
    cap2 = pre.outl['cap2.95%', i] 
    treatOutlier(x = test[ ,i], impute = T, fill = T, lower = lower,
                 upper = upper, med = med, cap1 = cap1, cap2 = cap2) 
  } -> test.out
  colnames(test.out) <- colnames(test)
  },
  env
)
#---------------
evalq(boxplot(train, main = "train  with outliers"), env)
evalq(boxplot(train.out, main = "train.out  without outliers"), env)
evalq(boxplot(test.out, main = "test.out  without outliers"), env)
#------------------------

Abb.18. Trainingsdatensatz mit Ausreißern

Abb.19. Trainingsdatensatz mit imputierten Ausreißern

Abb.20. Testdatensatz mit imputierten Ausreißern

Nicht alle Modelle sind empfindlich gegenüber Ausreißern. Unempfindlich sind zum Beispiel die Modelle der "Entscheidungsbäume"(DT) und der "Random Forest" (RF). 

Beim Erkennen und Verarbeiten von Ausreißern können auch weitere Pakete hilfreich sein, dazu gehören “univOutl”, “mvoutlier”, “outlier”, funModeling::prep.outlier().

3.3. Beseitigen der Schiefe (skewness)

evalq({
  sk <- skewness(x) 
  sk.out <- skewness(x.out) 
  sk.cap <- skewness(x.cap)
  }, 
  env)
> env$sk
             ftlm     stlm     rbci     pcci   v.fatl
Skewness 4.219857 1.785286 2.304655 6.491546 5.274871
           v.satl   v.rftl   v.rstl   v.ftlm    v.stlm
Skewness 2.677162 3.954098 1.568675 1.207227 0.8516043
           v.pcci
Skewness 3.031012
> env$sk.out
                ftlm        stlm        rbci       pcci
Skewness -0.04272076 -0.07893945 -0.02460354 0.01485785
             v.fatl      v.satl      v.rftl      v.rstl
Skewness 0.00780424 -0.02640635 -0.04663711 -0.04290957
                v.ftlm     v.stlm       v.pcci
Skewness -0.0009597876 0.01997082 0.0007462494
> env$sk.cap
                ftlm        stlm        rbci       pcci
Skewness -0.03329392 -0.07911245 -0.02847851 0.01915228
             v.fatl      v.satl      v.rftl      v.rstl
Skewness 0.01412182 -0.02617518 -0.03412228 -0.04596505
              v.ftlm      v.stlm      v.pcci
Skewness 0.008181183 0.009661169 0.002252508

require(PerformanceAnalytics)
evalq({
  k <- kurtosis(x) 
  k.out <- kurtosis(x.out) 
  k.cap <- kurtosis(x.cap)
}, 
env)
> env$k
                    ftlm     stlm     rbci     pcci
Excess Kurtosis 84.61177 16.77141 45.01858 247.9795
                  v.fatl   v.satl  v.rftl   v.rstl
Excess Kurtosis 145.9547 36.99944 74.4307 13.57613
                  v.ftlm   v.stlm   v.pcci
Excess Kurtosis 86.36448 23.06635 233.5408
> env$k.out
                        ftlm       stlm       rbci
Excess Kurtosis -0.003083449 -0.1668102 -0.1197043
                       pcci      v.fatl      v.satl
Excess Kurtosis -0.05113439 -0.02738558 -0.04341552
                     v.rftl     v.rstl     v.ftlm
Excess Kurtosis -0.01219999 -0.1316499 -0.0287925
                    v.stlm      v.pcci
Excess Kurtosis -0.1530424 -0.09950709
> env$k.cap
                      ftlm       stlm       rbci
Excess Kurtosis -0.2314336 -0.3075185 -0.2982044
                      pcci     v.fatl     v.satl
Excess Kurtosis -0.2452504 -0.2389486 -0.2331203
                    v.rftl     v.rstl     v.ftlm
Excess Kurtosis -0.2438431 -0.2673441 -0.2180059
                    v.stlm     v.pcci
Excess Kurtosis -0.2763058 -0.2698028

Im Ausgangsdatensatz x sind die Spitzen der Verteilung der Variablen sehr steil (Kurtosis viel größer als 0). Im Datensatz x.out mit entfernten Ausreißern sind die Hochs nahe an einer "normalen" Steilheit, im Datensatz mit imputierten Ausreißern sind die Wölbungen flacher. Beide Datensätze erfordern keine Korrekturen.

Anhang

Links

[3] Breuning, M., Kriegel, H., Ng, R.T, and Sander. J. (2000). LOF: Identifying density-based local outliers. In Proceedings of the ACM SIGMOD International Conference on Management of Data.