Forschungsmethoden und Analyseverfahren in der Biomechanik
Datum
15. Juni 2023
1 Datensatz: palmerpenguins
Um das Procedere einer Clusteranalyse zu illustrieren, beginnen wir mit einem Datensatz, der zu diesem Zweck seit neuestem auch in einigen Web-Tutorials zum Thema Clusteranalyse genutzt wird: dem Datensatz palmerpenguins. Er ist durch Installation des Pakets palmerpenguins verfügbar und enthält biometrische Daten zu unterschiedlichen Pinguin-Arten aus einer Forschungsstation (Palmer Station) in der Antarktis.
# A tibble: 6 × 8
species island bill_length_mm bill_depth_mm flipper_length_mm body_mass_g
<fct> <fct> <dbl> <dbl> <int> <int>
1 Adelie Torgersen 39.1 18.7 181 3750
2 Adelie Torgersen 39.5 17.4 186 3800
3 Adelie Torgersen 40.3 18 195 3250
4 Adelie Torgersen NA NA NA NA
5 Adelie Torgersen 36.7 19.3 193 3450
6 Adelie Torgersen 39.3 20.6 190 3650
# ℹ 2 more variables: sex <fct>, year <int>
Wir bilden zur Vereinfachung zunächst einen Subdatensatz, der nur die Variablen species und die biometrischen Variablen zu den Flossen und Schnäbeln und zudem nur vollständige Datensätze enthält.
Folgende Grafik illustriert die Bedeutung der Variablen bill_length_mm und bill_depth_mm im Datensatz:
1.1 Datentransformation
In den meisten Problemstellungen, die man mit Clusteranalysen untersucht, wird man teilweise sehr unterschiedliche Variablen zur Identifikation der Cluster einbinden. Diese haben zumeist auch eine sehr unterschiedliche Ausprägung, Streuung oder auch eine komplett andere Skalierung und Metrik. Im Beispiel der Pinguine finden wir Flossenlängen im Bereich von ca. 170 bis 230 mm, Schnabellängen hingegen nur im Bereich zwischen 30 und 60 mm, beim Körpergewicht wären es Zahlenwerte von mehreren Tausend Gramm. Um nicht jene Variablen mit der größten Ausprägung oder Streuung stärker zu gewichten, ist es ratsam, die Variablen vorher in z-Werte umzurechnen, und alle weiteren Prozeduren dann mit z-standardisierten Werten zu unternehmen. Dies stellt sicher, dass wir alle Variablen auf eine einheitliche und vergleichbare Metrik bringen, bevor wir die Gruppen anhand der Variablenausprägung identifizieren.
1.2 Prüfung der Korrelation der beteiligten Variablen
Bevor man mit einer Clusteranalyse beginnt, ist es sinnvoll, die in Frage kommenden Variablen auf mögliche Korrelation untereinander zu prüfen. Variablen, die hoch miteinander korrelieren, bringen keine neuen Gruppierungsinformationen mit und sollten nur in Form einer der beiden Variablen im Datensatz zur Clusteranalyse repräsentiert sein.
Da alle drei biometrischen Variablen metrischer Natur sind, können wir sie direkt gegeneinander in Form einer Scatterplotmatrix ansehen:
plot(PENG[,c("FLz", "BLz", "BDz", "BMz")])
Wir sehen, dass bestenfalls eine der Paarungen der Variablen auf “problemlose” Art miteinander korreliert, sondern je nach Variablenpaarung eine Gruppenbildung erkennbar wird. Einzig die Paarung FLz und BMz zeigt eine Korrelation, die weniger deutlich auf eine Subgruppenbildung hinweist. Wir werden daher im ersten Durchgang die Variable des Körpergewichts der Pinguine nicht berücksichtigen, da wir einen Großteil derselben Info durch die Flossenlänge abbilden können.
Im realen Datensatz sind de facto drei Gruppen (bzw. Arten von Pinguinen) enthalten, je nach Variablenpaarung erkennen wir eine Gruppenbildung zwischen 2 und 4 Clustern. Inwiefern diese jeweils aus zwei Variablen im Scatterplot erkennbaren Gruppen irgendetwas mit den tatsächlich vorhandenen Arten der Pinguine zu tun haben, können wir an dieser Stelle noch nicht beurteilen. Wir können jedoch erkennen, dass wir alle drei Variablen in die Clusteranalyse integrieren können, da keine der Variablen stark mit einer anderen korreliert und diese damit aus der Analyse verdrängen würde.
Wir werden versuchen, die in den Daten enthaltenen drei Arten, die in der Variable species codiert sind, möglichst korrekt mittels Clusteranalyse anhand ihrer biometrischen Merkmale zu identifizieren.
1.3 Hierarchische Clusterlösung
Angenommen, man wüsste nicht, wie viele verschiedene Arten von Pinguinen in diesen Daten enthalten sind, würde man sich möglicherweise mit einer hierarchischen Clusteranalyse dem Problem nähern. Diese definiert sich zunächst durch die Wahl eines Proximitätsmaßes, welches eine Matrix aus Proximitätswerten erzeugt ( Distanzmatrix), und in Folge durch die Wahl eines Fusionierungsalgorithmus.
Die so entstandene Hierarchie unterschiedlicher Clusterlösungen muss zunächst auf Plausibilität und ihre Passung zur Datenlage und eventuellen Erwartungen geprüft werden. Dazu eignet sich bei hierarchischen Methoden vorzugsweise das Dendrogramm:
# Dendrogramm, die Option `labels=FALSE` verhindert die Ausgabe der für uns unnützen Datenlabelsplot(hCluster, labels=FALSE)
Im Dendrogramm ist die Länge der vertikalen Linien zwischen den horizontalen Verbindungen entscheidend. Sie repräsentiert die Heterogenität der Gruppen, welche in dem jeweiligen Schritt vereinigt werden. Je höher die Querverbindung gegenüber dem darunter liegenden Vereinigungsschritt liegt und je länger die darunterliegenden vertikalen Linien ausfallen, desto heterogenere Gruppen werden in diesem Schritt vereinigt. Trennen wir das Dendrogramm horizontal auf Höhe ca. 25, erhalten wir eine 2-Gruppen-Lösung, da wir zwei verbleibende vertikale Linien durchtrennen. Trennen wir bei ca. 15, erhalten wir eine 3-Gruppen-Lösung. Etwas darunter dann jeweils eine 4-, 5- und 6-Gruppen-Lösungen. Alle weiteren, feineren Unterteilungen können wir ohnehin ignorieren.
Wir erkennen, dass die Vereinigung der zwei großen letzten Gruppen die höchste Heterogenität erzeugt. Dies ist keine Seltenheit und auch in anderen Datensätzen oftmals der Fall. Allerdings erzeugt auch die Vereinigung von drei auf zwei Gruppen (im rechten Teil des Dendrogramms auf Höhe von ca. 18) eine große Heterogenität in dieser Gruppe. Unterhalb der 3-Gruppen-Lösung erscheint jedoch keine klar erkennbare höhere Anzahl an Gruppen als bessere Lösung. Wir sehen uns daher die 2- und 3-Gruppen-Lösungen etwas näher an.
# Clusterlösungen mit 2 bis 4 Gruppen speichernPENG$hC2 <-factor(cutree(hCluster, k=2), labels =c("C2.1", "C2.2"))PENG$hC3 <-factor(cutree(hCluster, k=3), labels =c("C3.1", "C3.2", "C3.3"))
In grafischer Form eines Scatterplots lassen sich die Lösungen wie folgt illustrieren (hierzu können wir auch die Daten in Originalmetrik verwenden):
Wir erkennen, dass wir mit einer hierarchischen 3-Cluster-Lösung die “wahren” Gruppierungen schon rein visuell sehr gut abbilden können. Wie gut die Passung unserer so ermittelten 3-Gruppen Lösung mit den tatsächlichen Arten ist, lässt sich anhand einer Konfusionsmatrix beurteilen:
Wir können erkennen, dass wir – bis auf wenige Ausnahmen – die wahren Arten durch eine Einteilung auf Basis der durch biometrische Merkmale gebildeter Cluster sehr gut vorhersagen können. Einzig die Art Chinstrap erscheint noch problematisch, da sie mit 7 Tieren im “falschen” Cluster auftaucht, bzw. das Cluster C3.1 sieben Chinstrap Tiere fälschlicherweise mit einschließt. Gentoo ist zu 100% korrekt durch das Cluster C3.3 vorhergesagt, und Adelie findet sich fast ausnahmslos in C3.1 wieder.
Eventuell ließe sich diese schon exzellente Klassifikation durch Hinzunahme weiterer Variablen verbessern. Da wir nur noch das Körpergewicht zur Verfügung und bisher nicht berücksichtigt haben, ziehen wir diese Variable hinzu:
An dieser Konfusionsmatrix erkennt man zwei wichtige Dinge:
Nicht immer wird durch Hinzunahme von Variablen ein besseres Klassifikationsergebnis erzielt. Im Fall des Körpergewichts der Pinguine führt dies dazu, dass sich die Art Chinstrap ungünstiger auf die drei CLuster verteilt als zuvor.
Nicht immer ist die Reihenfolge der Cluster gleichbedeutend. Im Fall oben war Gentoo in C3.3 repräsentiert, nun ist dies in C3.2 der Fall. Hier muss man immer inhaltlich zuordnen, welches der Cluster inhaltlich welche Eigenschaft abbildet. Hierbei kann ein Profillinienplot sehr hilfreich sein, da eine reale Clusteranalyse ja womöglich ohne die Kenntnis der wahren Pinguinarten interpretiert werden müsste.
Die anhand der Konfusionsmatrix bessere Klassifikation erscheint mit der 3-Clusterlösung aus der vorherigen Analyse erreicht worden zu sein. Wir verwenden daher für die folgenden Schritte auch diese Lösung in Form der Variable hC3.
1.3.1 Profilplot über alle Variableneigenschaften
Je nach Anzanhl der zur Clusterbildung verwendeten Variablen zeigen Abbildungen wie jene oben ja nur noch einen Ausschnitt der genutzten Daten an. Was im zweidimensionalen Raum als nicht eindeutige Trennung (bzw. große Überlappung) zweier Gruppen erscheint, kann im höherdimensionalen Raum mit z. B. mehr als 3 Variablen eindeutiger ausfallen. Dazu kann man sich anstelle von Scatterplots einen Profillinienplot anfertigen.
Um dies mit ggplot2 ökonomisch zu erledigen, müssen wir die Daten zuvor aggregieren. Wir verwenden dazu wiederum die z-standardisierten Merkmale und diesmal das Paket plyr:
Man sieht, dass die drei identifizierten Cluster jeweils eigene “Muster” erzeugen.C3.1 enthält die Tiere mit unterdurchschnittlich langen Flossen und Schnäbeln, die dafür aber “dickere” (bzw. höhere) Schnäbel zeigen. C3.2enthält Tiere wie C3.1, allerdings mit überduchschnittlich langen Schnäbeln. C3.3 ist diametral zu C3.1, d. h. überdurchschnittlich lange Flossen und Schnäbel, dafür flachere (weniger hohe) Schnäbel.
Offene Punkte einer hierarchischen Clusterlösung bestehen also vorwiegend in der Optimierung der erhaltenen Lösung bzw. in der korrekten Zuordnung der wenigen falsch klassifizierten Tiere. Dies werden wir mit der hierarchischen Methode alleine nicht erreichen können.
1.4 Partitionierendes Verfahren mit k-means
Wenn wir also eine Vorstellung davon haben, dass wir es mit drei Punguinarten zu tun haben, können wir entweder direkt mit einem partitionierenden Verfahren einsteigen oder aber dieses dazu nutzen, die bestehende Lösung aus einem hierarchischen Vorgehen zu optimieren.
Die Methode nach k-means erlaubt es, auch ohne Angabe von Clusterzentren direkt eine 3-Cluster-Lösung zu generieren:
Wir sehen, dass diese Lösung ähnlich gute Ergebnisse wie die hierarchische Lösung erzeugt. Wiederum ist die Trennung der Arten Adelie und Chinstrap das “Problem”, denn einige der Tiere dieser Arten tauchen im jeweils falschen Cluster auf.
Man könnte nun versuchen, die Leistung von k-means dahgingehend zu optimieren, in dem man die aus einer hierarchischen 3-Cluster-Lösung bereits errechneten Zentren als Startwerte für das k-means Verfahren vorgibt:
# Matrix mit Startwerten für Zentren erhalten wir in Form einer Aggregierung aus dem bereits erzeugten `MW`c.init <- MW[,2:4]# kmeans rechnenkCluster3.opt <-kmeans(PENG[,c("FLz", "BLz", "BDz")], centers = c.init)PENG$C3.opt <-factor(kCluster3.opt$cluster, labels=c("C1.opt", "C2.opt", "C3.opt"))table(PENG$species, PENG$C3.opt)
Leider ändert sich durch die Optimierung die Fehlklassifikation nicht bedeutsam, so dass man mit einem Verbleibenden Anteil von insgesamt ca. 3% (10 von 342) falsch klassifizierter Tiere leben muss. Sieht man sich die Lage der Cluster in der Abbildung realer Pinguinarten genauer an, so erkennt man, dass einige der Tiere aus Adelie schlichtweg die Eigenschaften von Chinstrap Tieren zeigen und umgekehrt. Anhand ihrer wenigen biometrischen Merkmale wird man sie daher mit Clusterbildenden verfahren kaum korrekt zuordnen können, da sie durch ihre “Ausreißerrolle” stets der anderen Gruppe näher sein werden als ihrer eigenen Art. Abhilfe könnte hierbei womöglich kein anderes Clusterbildungsverfahren schaffen. Dies geht allerdings über den Scope dieser Übung hinaus. Interessierte sind auf das Tutorial von Mouselimis verwiesen.
Wertvoller wäre, um die Arten besser trennen zu können, zusätzliche Variablen zu erheben, welche eine bessere Trennung von Adelie und Chinstrap ermöglichen. Die wenigen hier genutzten biometrischen Variablen sind dazu offensichtlich nicht in der Lage.
