# Stochastik IV, Übungsblatt 4

# Aufgabe 4

library(ellipse)

# Einlesen der Daten und Einbinden in Suchpfad

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

attach(data.sweat)

# Scatterplot Matrix (ohne "Diagonale")

par(mfrow=c(3,3))

frame()

# (1) SWEAT vs. SODIUM

plot(SODIUM, SWEAT, pch = 19)
points(mean(SODIUM), mean(SWEAT), pch = 24)
points(50, 4, pch = 22)

# Konfidenzellipse SWEAT vs. SODIUM

lines(ellipse(cov(cbind(SODIUM, SWEAT)), 
              t=sqrt(((3 * (20 - 1)) / (20 * (20 - 3))) * qf(0.95, 3, 20 - 3)),
   	         centre=c(mean(SODIUM), mean(SWEAT))), col = 2)

# Simultane Konfidenzintervalle SWEAT vs. SODIUM

abline(v = (mean(SODIUM) + (sqrt(((3 * (20 - 1)) / (20 - 3)) 
                                   * qf(0.95, 3, 20 - 3)) *
                            sqrt(var(SODIUM) / 20))), col = 3)
abline(v = (mean(SODIUM) - (sqrt(((3 * (20 - 1)) / (20 - 3)) 
                                   * qf(0.95, 3, 20 - 3)) *
                            sqrt(var(SODIUM) / 20))), col = 3)
abline(h = (mean(SWEAT) + (sqrt(((3 * (20 - 1)) / (20 - 3)) 
                                   * qf(0.95, 3, 20 - 3)) *
                         sqrt(var(SWEAT) / 20))), col = 3)
abline(h = (mean(SWEAT) - (sqrt(((3 * (20 - 1)) / (20 - 3)) 
                                   * qf(0.95, 3, 20 - 3)) *
                         sqrt(var(SWEAT) / 20))), col = 3)

# Bonferroni korrigierte Konfidenzintervalle SWEAT vs. SODIUM

abline(v = (mean(SODIUM) + (qt(1 - (0.05 / (2 * 3)), 20 - 1) 
            * sqrt(var(SODIUM) / 20))), col = 4)
abline(v = (mean(SODIUM) - (qt(1 - (0.05 / (2 * 3)), 20 - 1) 
            * sqrt(var(SODIUM) / 20))), col = 4)
abline(h = (mean(SWEAT) + (qt(1 - (0.05 / (2 * 3)), 20 - 1) 
            * sqrt(var(SWEAT) / 20))), col = 4)
abline(h = (mean(SWEAT) - (qt(1 - (0.05 / (2 * 3)), 20 - 1) 
            * sqrt(var(SWEAT) / 20))), col = 4)

# (2) SWEAT vs. POTASSIUM

plot(POTASS, SWEAT, pch = 19)
points(mean(POTASS), mean(SWEAT), pch = 24)
points(10, 4, pch = 22)

# Konfidenzellipse SWEAT vs. POTASSIUM

lines(ellipse(cov(cbind(POTASS, SWEAT)), 
              t=sqrt(((3 * (20 - 1)) / (20 * (20 - 3))) * qf(0.95, 3, 20 - 3)),
   	         centre=c(mean(POTASS), mean(SWEAT))), col = 2)

# Simultane Konfidenzintervalle SWEAT vs. POTASSIUM

abline(v = (mean(POTASS) + (sqrt(((3 * (20 - 1)) / (20 - 3)) 
                                   * qf(0.95, 3, 20 - 3)) *
                            sqrt(var(POTASS) / 20))), col = 3)
abline(v = (mean(POTASS) - (sqrt(((3 * (20 - 1)) / (20 - 3)) 
                                   * qf(0.95, 3, 20 - 3)) *
                            sqrt(var(POTASS) / 20))), col = 3)
abline(h = (mean(SWEAT) + (sqrt(((3 * (20 - 1)) / (20 - 3)) 
                                   * qf(0.95, 3, 20 - 3)) *
                         sqrt(var(SWEAT) / 20))), col = 3)
abline(h = (mean(SWEAT) - (sqrt(((3 * (20 - 1)) / (20 - 3)) 
                                   * qf(0.95, 3, 20 - 3)) *
                         sqrt(var(SWEAT) / 20))), col = 3)

# Bonferroni korrigierte Konfidenzintervalle SWEAT vs. POTASSIUM

abline(v = (mean(POTASS) + (qt(1 - (0.05 / (2 * 3)), 20 - 1) 
            * sqrt(var(POTASS) / 20))), col = 4)
abline(v = (mean(POTASS) - (qt(1 - (0.05 / (2 * 3)), 20 - 1) 
            * sqrt(var(POTASS) / 20))), col = 4)
abline(h = (mean(SWEAT) + (qt(1 - (0.05 / (2 * 3)), 20 - 1) 
            * sqrt(var(SWEAT) / 20))), col = 4)
abline(h = (mean(SWEAT) - (qt(1 - (0.05 / (2 * 3)), 20 - 1) 
            * sqrt(var(SWEAT) / 20))), col = 4)

# (3) SODIUM vs. SWEAT

plot(SWEAT, SODIUM, pch = 19)
points(mean(SWEAT), mean(SODIUM), pch = 24)
points(4, 50, pch = 22)

# Konfidenzellipse SODIUM vs. SWEAT

lines(ellipse(cov(cbind(SWEAT, SODIUM)), 
              t=sqrt(((3 * (20 - 1)) / (20 * (20 - 3))) * qf(0.95, 3, 20 - 3)),
   	         centre=c(mean(SWEAT), mean(SODIUM))), col = 2)

# Simultane Konfidenzintervalle SODIUM vs. SWEAT

abline(v = (mean(SWEAT) + (sqrt(((3 * (20 - 1)) / (20 - 3)) 
                                   * qf(0.95, 3, 20 - 3)) *
                            sqrt(var(SWEAT) / 20))), col = 3)
abline(v = (mean(SWEAT) - (sqrt(((3 * (20 - 1)) / (20 - 3)) 
                                   * qf(0.95, 3, 20 - 3)) *
                            sqrt(var(SWEAT) / 20))), col = 3)
abline(h = (mean(SODIUM) + (sqrt(((3 * (20 - 1)) / (20 - 3)) 
                                   * qf(0.95, 3, 20 - 3)) *
                         sqrt(var(SODIUM) / 20))), col = 3)
abline(h = (mean(SODIUM) - (sqrt(((3 * (20 - 1)) / (20 - 3)) 
                                   * qf(0.95, 3, 20 - 3)) *
                         sqrt(var(SODIUM) / 20))), col = 3)

# Bonferroni korrigierte Konfidenzintervalle SODIUM vs. SWEAT

abline(v = (mean(SWEAT) + (qt(1 - (0.05 / (2 * 3)), 20 - 1) 
            * sqrt(var(SWEAT) / 20))), col = 4)
abline(v = (mean(SWEAT) - (qt(1 - (0.05 / (2 * 3)), 20 - 1) 
            * sqrt(var(SWEAT) / 20))), col = 4)
abline(h = (mean(SODIUM) + (qt(1 - (0.05 / (2 * 3)), 20 - 1) 
            * sqrt(var(SODIUM) / 20))), col = 4)
abline(h = (mean(SODIUM) - (qt(1 - (0.05 / (2 * 3)), 20 - 1) 
            * sqrt(var(SODIUM) / 20))), col = 4)

frame()

# (4) SODIUM vs. POTASSIUM

plot(POTASS, SODIUM, pch = 19)
points(mean(POTASS), mean(SODIUM), pch = 24)
points(10, 50, pch = 22)

# Konfidenzellipse SODIUM vs. POTASSIUM

lines(ellipse(cov(cbind(POTASS, SODIUM)), 
              t=sqrt(((3 * (20 - 1)) / (20 * (20 - 3))) * qf(0.95, 3, 20 - 3)),
   	         centre=c(mean(POTASS), mean(SODIUM))), col = 2)

# Simultane Konfidenzintervalle SODIUM vs. POTASSIUM

abline(v = (mean(POTASS) + (sqrt(((3 * (20 - 1)) / (20 - 3)) 
                                   * qf(0.95, 3, 20 - 3)) *
                            sqrt(var(POTASS) / 20))), col = 3)
abline(v = (mean(POTASS) - (sqrt(((3 * (20 - 1)) / (20 - 3)) 
                                   * qf(0.95, 3, 20 - 3)) *
                            sqrt(var(POTASS) / 20))), col = 3)
abline(h = (mean(SODIUM) + (sqrt(((3 * (20 - 1)) / (20 - 3)) 
                                   * qf(0.95, 3, 20 - 3)) *
                            sqrt(var(SODIUM) / 20))), col = 3)
abline(h = (mean(SODIUM) - (sqrt(((3 * (20 - 1)) / (20 - 3)) 
                                   * qf(0.95, 3, 20 - 3)) *
                           sqrt(var(SODIUM) / 20))), col = 3)

# Bonferroni korrigierte Konfidenzintervalle SODIUM vs. POTASSIUM

abline(v = (mean(POTASS) + (qt(1 - (0.05 / (2 * 3)), 20 - 1) 
            * sqrt(var(POTASS) / 20))), col = 4)
abline(v = (mean(POTASS) - (qt(1 - (0.05 / (2 * 3)), 20 - 1) 
            * sqrt(var(POTASS) / 20))), col = 4)
abline(h = (mean(SODIUM) + (qt(1 - (0.05 / (2 * 3)), 20 - 1) 
            * sqrt(var(SODIUM) / 20))), col = 4)
abline(h = (mean(SODIUM) - (qt(1 - (0.05 / (2 * 3)), 20 - 1) 
            * sqrt(var(SODIUM) / 20))), col = 4)

# (5) POTASSIUM vs. SWEAT

plot(SWEAT, POTASS, pch = 19)
points(mean(SWEAT), mean(POTASS), pch = 24)
points(4, 10, pch = 22)

# Konfidenzellipse POTASSIUM vs. SWEAT

lines(ellipse(cov(cbind(SWEAT, POTASS)), 
              t=sqrt(((3 * (20 - 1)) / (20 * (20 - 3))) * qf(0.95, 3, 20 - 3)),
   	         centre=c(mean(SWEAT), mean(POTASS))), col = 2)

# Simultane Konfidenzintervalle POTASSIUM vs. SWEAT

abline(v = (mean(SWEAT) + (sqrt(((3 * (20 - 1)) / (20 - 3)) 
                                   * qf(0.95, 3, 20 - 3)) *
                            sqrt(var(SWEAT) / 20))), col = 3)
abline(v = (mean(SWEAT) - (sqrt(((3 * (20 - 1)) / (20 - 3)) 
                                   * qf(0.95, 3, 20 - 3)) *
                            sqrt(var(SWEAT) / 20))), col = 3)
abline(h = (mean(POTASS) + (sqrt(((3 * (20 - 1)) / (20 - 3)) 
                                   * qf(0.95, 3, 20 - 3)) *
                         sqrt(var(POTASS) / 20))), col = 3)
abline(h = (mean(POTASS) - (sqrt(((3 * (20 - 1)) / (20 - 3)) 
                                   * qf(0.95, 3, 20 - 3)) *
                         sqrt(var(POTASS) / 20))), col = 3)

# Bonferroni korrigierte Konfidenzintervalle POTASSIUM vs. SWEAT

abline(v = (mean(SWEAT) + (qt(1 - (0.05 / (2 * 3)), 20 - 1) 
            * sqrt(var(SWEAT) / 20))), col = 4)
abline(v = (mean(SWEAT) - (qt(1 - (0.05 / (2 * 3)), 20 - 1) 
            * sqrt(var(SWEAT) / 20))), col = 4)
abline(h = (mean(POTASS) + (qt(1 - (0.05 / (2 * 3)), 20 - 1) 
            * sqrt(var(POTASS) / 20))), col = 4)
abline(h = (mean(POTASS) - (qt(1 - (0.05 / (2 * 3)), 20 - 1) 
            * sqrt(var(POTASS) / 20))), col = 4)

# (6) POTASSIUM vs. SODIUM

plot(SODIUM, POTASS, pch = 19)
points(mean(SODIUM), mean(POTASS), pch = 24)
points(50, 10, pch = 22)

# Konfidenzellipse POTASSIUM vs. SODIUM

lines(ellipse(cov(cbind(SODIUM, POTASS)), 
              t=sqrt(((3 * (20 - 1)) / (20 * (20 - 3))) * qf(0.95, 3, 20 - 3)),
   	         centre=c(mean(SODIUM), mean(POTASS))), col = 2)

# Simultane Konfidenzintervalle POTASSIUM vs. SODIUM

abline(v = (mean(SODIUM) + (sqrt(((3 * (20 - 1)) / (20 - 3)) 
                                   * qf(0.95, 3, 20 - 3)) *
                            sqrt(var(SODIUM) / 20))), col = 3)
abline(v = (mean(SODIUM) - (sqrt(((3 * (20 - 1)) / (20 - 3)) 
                                   * qf(0.95, 3, 20 - 3)) *
                            sqrt(var(SODIUM) / 20))), col = 3)
abline(h = (mean(POTASS) + (sqrt(((3 * (20 - 1)) / (20 - 3)) 
                                   * qf(0.95, 3, 20 - 3)) *
                            sqrt(var(POTASS) / 20))), col = 3)
abline(h = (mean(POTASS) - (sqrt(((3 * (20 - 1)) / (20 - 3)) 
                                   * qf(0.95, 3, 20 - 3)) *
                           sqrt(var(POTASS) / 20))), col = 3)

# Bonferroni korrigierte Konfidenzintervalle POTASSIUM vs. SODIUM

abline(v = (mean(SODIUM) + (qt(1 - (0.05 / (2 * 3)), 20 - 1) 
            * sqrt(var(SODIUM) / 20))), col = 4)
abline(v = (mean(SODIUM) - (qt(1 - (0.05 / (2 * 3)), 20 - 1) 
            * sqrt(var(SODIUM) / 20))), col = 4)
abline(h = (mean(POTASS) + (qt(1 - (0.05 / (2 * 3)), 20 - 1) 
            * sqrt(var(POTASS) / 20))), col = 4)
abline(h = (mean(POTASS) - (qt(1 - (0.05 / (2 * 3)), 20 - 1) 
            * sqrt(var(POTASS) / 20))), col = 4)

frame()

legend(0.03, 0.86, c("Ellipse", "Simultane", "Bonferroni"), bty = "n",
       col = c(2,3,4), lty = 1)
points(0.20, 0.21, pch = 24)
text(0.41, 0.21, "X quer ")
points(0.20, 0.10, pch = 22)
text(0.41, 0.11, "mu null ")