#### Stochastik IV, Übungsblatt 6

### Aufgabe 4

#### Achtung: Bei Ausführung des kompletten Skripts beachte die rechnerintensiven 
#### Berechnungen der minimalen Stressfunktionswerte in Abhängigkeit der Dimension, 
#### die bei der hier im Skript notierten 'Isotonic Kruskal-Shepard' 
#### multidimensionalen Skalierung (MDS) nötig sind!

library(MASS)

# Einlesen der Daten und Einbinden in Suchpfad

data.crabs <- 
read.table("D:/Univ Augsburg/Lehre/SS07/VO_MultiVar/Datensaetze/crabs.txt",
head = TRUE, sep = "\t", quote = "")

attach(data.crabs)

crabs <- data.crabs[, 4:8]  # die fünf morphologischen Messvariablen

## Teilaufgabe (a)

# i. Klassische multidimensionale Skalierung (MDS) 

crabs.klassisch <- cmdscale(dist(crabs))

# ii. 'Sammon Mapping' MDS 

crabs.sammon <- sammon(dist(crabs))

# iii. 'isotonic Kruskal-Shepard' MDS 

crabs.isotonicKS <- isoMDS(dist(crabs))

# Graphische Darstellung der 2D Repräsentationen der Daten in einer Figur,
# mit Farb- und Formunterscheidung der vier Gruppen

par(mfrow = c(2, 2))

plot(crabs.klassisch, main = "Klassische multidimensionale Skalierung",
     xlab = "Dimension 1", ylab = "Dimension 2",
     cex = 1.5, col = (sp - 1) + (2 * (sex - 1)) + 1, 
     pch = (sp - 1) + (2 * (sex - 1)) + 1)

# NOTIZ: "Blue Male": col=1, pch=1; "Blue Female": col=3, pch=3; 
#        "Orange Male": col=2, pch=2; "Orange Female": col=4, pch=4

plot(crabs.sammon$points, main = "'Sammon Mapping' multidimensionale Skalierung",
     xlab = "Dimension 1", ylab = "Dimension 2",
     cex = 1.5, col = (sp - 1) + (2 * (sex - 1)) + 1, 
     pch = (sp - 1) + (2 * (sex - 1)) + 1)

plot(crabs.isotonicKS$points, main = "'Isotonic Kruskal-Shepard' multidimensionale Skalierung",
     xlab = "Dimension 1", ylab = "Dimension 2",
     cex = 1.5, col = (sp - 1) + (2 * (sex - 1)) + 1, 
     pch = (sp - 1) + (2 * (sex - 1)) + 1)

frame()

legend(0.03, 0.86, pch = c(1, 3, 2, 4), col = c(1, 3, 2, 4),
       c("Blue Male", "Blue Female", "Orange Male", "Orange Female"))

# Diagnostische Shepard-Diagramme

crabs.klassisch.sh <- Shepard(dist(crabs), crabs.klassisch)
crabs.sammon.sh <- Shepard(dist(crabs), crabs.sammon$points)
crabs.isotonicKS.sh <- Shepard(dist(crabs), crabs.isotonicKS$points)

# Klassische MDS

windows()

plot(crabs.klassisch.sh, pch = ".", 
     main = "Shepard-Diagramm klassische multidimensionale Skalierung",
     xlab = "Ursprüngliche Distanzen", ylab = "Gefittete Distanzen", 
     col = "green", cex = 2)
lines(crabs.klassisch.sh$x, crabs.klassisch.sh$yf, type ="S", col = "red", lwd = 2)
abline(0,1, col = "blue", lwd = 2, lty = 2)
legend(30, 20, lty = c(1, 3, 2), col = c(2, 3, 4), lwd = 2,
       c("Treppenfunktion", "Distanzpaare", "Gerade y = x"))

# 'Sammon Mapping' MDS

windows()

plot(crabs.sammon.sh, pch = ".", 
     main = "Shepard-Diagramm 'Sammon Mapping' multidimensionale Skalierung",
     xlab = "Ursprüngliche Distanzen", ylab = "Gefittete Distanzen", 
     col = "green", cex = 2)
lines(crabs.sammon.sh$x, crabs.sammon.sh$yf, type ="S", col = "red", lwd = 2)
abline(0,1, col = "blue", lwd = 2, lty = 2)
legend(30, 20, lty = c(1, 3, 2), col = c(2, 3, 4), lwd = 2,
       c("Treppenfunktion", "Distanzpaare", "Gerade y = x"))

# 'Isotonic Kruskal-Shepard' MDS

windows()

plot(crabs.isotonicKS.sh, pch = ".", 
     main = "Shepard-Diagramm 'Isotonic Kruskal-Shepard' multidimensionale Skalierung",
     xlab = "Ursprüngliche Distanzen", ylab = "Gefittete Distanzen", 
     col = "green", cex = 2)
lines(crabs.isotonicKS.sh$x, crabs.isotonicKS.sh$yf, type ="S", col = "red", lwd = 2)
abline(0,1, col = "blue", lwd = 2, lty = 2)
legend(30, 20, lty = c(1, 3, 2), col = c(2, 3, 4), lwd = 2,
       c("Treppenfunktion", "Distanzpaare", "Gerade y = x"))

## Teilaufgabe (b)

# 'Sammon Mapping' MDS

stress.sammon <- 1:7
for ( i in 1:7 )
    stress.sammon[i] <- sammon(dist(crabs), k = i)$stress

plot(stress.sammon, 
     main = "Stressfunktion 'Sammon Mapping' multidimensionale Skalierung",
     xlab = "Anzahl Dimensionen", 
     ylab = "Minimaler Stressfunktionswert", type ="o", col = "lightblue", lwd = 2)
abline(h = 0, lty = 2)

# Ellenbogen-Kriterium: 'k = 2 oder 3'; im Fall k = 2 ist die graphische 
# Veranschaulichung natürlich leichter als im Fall k = 3.

# 'Isotonic Kruskal-Shepard' MDS 
# Achtung: Berechnungen dauern 'eine Weile'!

stress.isotonicKS <- 1:7
for ( i in 1:7 )
    stress.isotonicKS[i] <- isoMDS(dist(crabs), k = i)$stress

plot(stress.isotonicKS, 
     main = "Stressfunktion 'Isotonic Kruskal-Shepard' multidimensionale Skalierung",
     xlab = "Anzahl Dimensionen", 
     ylab = "Minimaler Stressfunktionswert", type ="o", col = "lightblue", lwd = 2)
abline(h = 0, lty = 2)

# Ellenbogen-Kriterium: k = 3 oder 4; im Gegensatz zu 'k = 1 und 2' ergibt sich hier 
# die Notwendigkeit für 'gute' graphische Veranschaulichungen.