#### Data laden #### Z = read.csv("Uebung3ss.csv", row.names=1) colnames(Z) summary(Z) #### Explorative Datenanalyse #### par(mfrow=c(2,4)) for(i in 1:8){ hist(Z[,i], main=colnames(Z)[i]) } for(i in 1:8){ hist(log(Z[,i]), main=colnames(Z)[i]) } par(mfrow=c(1,1)) pc = princomp(Z, cor = TRUE) biplot(pc) ### x, A, D, E gehen zusammen; y, B1 (inverse), B2, B3 gehen zusammen; die 2 gruppen sind orthogonal # Y={A, B1, B2, B3, D, E } sind eine Zusammensetzung # library("compositions") # Pakett um Zusammensetzungen zu analysieren Y = acomp(Z[, -(1:2)]) # definiere eine relative Komposition ('A'itchison 'COMP'osition) colnames(Y) ## Dreieckdiagramme ## plot(acomp(Y[,c("A","B1","B2")])) plot(acomp(Y[,c("A","B1","B2")]), center=TRUE) #### analysieren Sie ihre Zusammenhänge mit PCA #### pc_comp = princomp(Y, cor = TRUE) # nonsense! fÜr Zusammensetzungen cor=FALSE biplot(pc_comp) pc_comp = princomp(Y, cor = FALSE) biplot(pc_comp) ## gibt es bessere Diagramme? JA: ## ## Dreieckdiagramme ## plot(acomp(Y[,c("B1","B2","B3")])) plot(acomp(Y[,c("A","D","E")])) plot(acomp(Y[,c("A","D","E")]), center=TRUE) plot(acomp(Y[,c("A","D","B1")])) ## conclusion: (B1,B2,B3) und (A,D,E) sind jeweils gut als 1-dimensionale Subkompositionen zu verstehen ## A vs. (B1,B2,B3) weisst darauf hin, dass die unhabgängig sein könnten #### Regression von Y vs. X={x,y} #### ## logratio transformation ## # alr, ilr, clr?? ## erst mal keine, sondern adhoc: (mdB12 = lm(log(B1/B2)~x+y, data=Z)) (mdB23 = lm(log(B2/B3)~x+y, data=Z)) (mdB13 = lm(log(B1/B3)~x+y, data=Z)) # FRAGE: Kann man eine der drei Modelle als lineare Kombination der anderen 2 bauen? coef(mdB12) coef(mdB23) coef(mdB13) ## wie? # FRAGE: sind alle Koeffitiente significant? summary(mdB12) summary(mdB23) summary(mdB13) ## testen was geschieht mit der Subkomposition (A,D,E) ## und mit einem "gemsichten" log-Verhältnis? (mdB1A = lm(log(B1/A)~x+y, data=Z)) summary(mdB1A) (mdB2A = lm(log(B2/A)~x+y, data=Z)) summary(mdB2A) coef(mdB1A) - coef(mdB2A) coef(mdB12) ## Gleichungen in einem Papier schreiben! ## gesamte Linearmodell ## (mdAll = lm(alr(Y)~data.matrix(Z[,1:2]))) # kann man die alr-Koeffiziente hier aus den Koeffizienten für Modellen mdB12, B13, B23, AD, AE, ED, B1A vorhersagen? # JA!: coef(mdAll)[,1] coef(mdAE) coef(mdAll)[,2] coef(mdB1A)-coef(mdAE) # usw ## Annahmen sind OK? Diagnostics ## par(mfrow=c(2,4)) plot(mdB12) (lmdB12 = lm(log(B1/B2)~log(x)+log(y), data=Z)) plot(lmdB12) # Modelle für alle andere log-Verhältnisse abgleichen ## Tests fuer Koeffiziente ## (lmdAll = lm(alr(Y)~log(data.matrix(Z[,1:2])))) summary(lmdAll) ## Kann man das Modell vereinfachen? wie? ## # Ja! # wir haben gesehen, dass (B1,B2,B3) = f(x) und (A,D,E) = f(y) # das bedeutet dass log(B1/E) = f(x,y), und das gilt auch für alle "Kreuzverhältnisse" # Lösung: adhoc logratios bauen, erst innerhalb die Subkompositionen, danach ein einzigen gekreuzten Bilanz: V = balanceBase(Y, expr=~( B1/(B2/B3) )/( (A/D)/E ) ) print(V) lrs = ilr(Y, V) (lmdAll2 = lm(lrs~log(data.matrix(Z[,1:2])))) summary(lmdAll2) ## kann man die Koeffizienten von lmdAll und lmdAll2 voneinander berechnen? coefficients(lmdAll2) coefficients(lmdAll)