Dzisiaj kontynuujemy rozpisywanie krok po kroku wizualizacji cen mieszkań. Ten wpis poświęcony jest wykresowi pudełkowemu. Kory w programie R użyte poniżej można znaleźć na stronie tutaj.
Podobnie jak poprzedni pierwsza część to przygotowanie danych, druga to ich wizualizacja.
Częśc 1.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33
|
library(lattice) library(RColorBrewer) library(MASS) # wczytujemy dane o mieszkaniach mieszkaniaKWW2011 <- read.table("http://tofesi.mimuw.edu.pl/~cogito/smarterpoland/mieszkaniaKWW2011/mieszkaniaKWW2011.csv", row.names=NULL, sep=";", header=TRUE, colClasses=c("factor", "factor", "numeric", "numeric", "factor", "numeric", "numeric", "factor", "Date")) # wybieramy tylko mieszkania z Warszawy o powierzchni do 300m2 mieszkaniaKWW2011Warszawa <- mieszkaniaKWW2011[mieszkaniaKWW2011$miasto =="Warszawa", ] mieszkaniaKWW2011Warszawa$dzielnica <- factor(mieszkaniaKWW2011Warszawa$dzielnica) mieszkaniaKWW2011Warszawa <- mieszkaniaKWW2011Warszawa[mieszkaniaKWW2011Warszawa$powierzchnia<300, ] # do analiz trafią tylko oferty złożone po 13 września 2007 dataStart <- as.Date("13-09-2007","%d-%m-%Y") mieszkaniaKWW2011Warszawa2 <- mieszkaniaKWW2011Warszawa mieszkaniaKWW2011Warszawa2 <- mieszkaniaKWW2011Warszawa2[mieszkaniaKWW2011Warszawa2$data >= dataStart,] # # dodajemy zmienną określającą miesiąc w którym złożona została oferta mieszkaniaKWW2011Warszawa2$dataF <- factor(as.character(mieszkaniaKWW2011Warszawa2$data, "%Y-%m")) # # z wszystkich meiszkań do obiektu mieszkaniaKWW2011Warszawa3 wybieramy tylko te o średnim rozmiarze, tz o powierzchni d 49 do 68m2 mieszkaniaKWW2011Warszawa2$rozmiar <- cut(mieszkaniaKWW2011Warszawa2$powierzchnia, c(0,49,68,200), labels=c("do 49 m2", "od 49 do 68 m2", "powyzej 68 m2")) mieszkaniaKWW2011Warszawa3 <- mieszkaniaKWW2011Warszawa2[which(mieszkaniaKWW2011Warszawa2$rozmiar == "od 49 do 68 m2"),] # # na potrzeby tego wykresu wybieramy tylko mieszkania z Zoliborza mieszkaniaKWW2011Warszawa3Zoliborz <- mieszkaniaKWW2011Warszawa3[which(mieszkaniaKWW2011Warszawa3$dzielnica == "Zoliborz"),] at = seq(1000,24000,1000) |
Dane są wczytane, czas na wykres. Wykorzystamy funkcję bwplot() z pakietu lattice.
Narysujemy jak zmieniają się ceny m2 średniej wielkości mieszkań w kolejnych miesiącach.
|
bwplot(cenam2~dataF, data=mieszkaniaKWW2011Warszawa3Zoliborz, main="Zoliborz") |

Nie wygląda to najlepiej, zajmijmy się na początek osiami. Ponieważ cena jest zmienną silnie prawo skośną, przedstawimy ją na osi logarytmicznej. Ponieważ etykiety na osi OX zachodzą na siebie to zmniejszymy je i pokażemy pionowo.
|
bwplot(cenam2~dataF, data=mieszkaniaKWW2011Warszawa3Zoliborz, scales=list(y=list(log=T, at=at), x=list(rot=90, cex=0.6)), main="Zoliborz") |

Przy takiej rozpiętości na osi OY trudno analizować delikatne zmiany w cenach mieszkań, więc w dalszej części zawęzimy zainteresowania do przedziału osi OY od 7 do 14k.
|
bwplot(cenam2~dataF, data=mieszkaniaKWW2011Warszawa3Zoliborz, scales=list(y=list(log=T, at=at), x=list(rot=90, cex=0.6)), ylim=c(7000, 14000), main="Zoliborz") |

Trudno zauważyć jakiś trend. Dorysujmy więc krzywą trendu liniowego przedefiniowując funkcję rysującą panel. Nowa funkcja rysująca panel rysuje linie pomocnicze siatki, rysuje wykresy pudełkowe używając panel.bwplot() oraz dorysowuje linię odpornej regresji liniowej.
|
bwplot(cenam2~dataF, data=mieszkaniaKWW2011Warszawa3Zoliborz, scales=list(y=list(log=T, at=at), x=list(rot=90, cex=0.6)), ylim=c(7000, 14000), main="Zoliborz", panel = function(...) { panel.abline(h=log(at,10), col="grey85") panel.abline(v=4.5 + c(0,12,24,36), col="grey85") panel.bwplot(..., col="grey20") args = list(...) mod = rlm(args$y~as.numeric(args$x)) panel.abline(mod) } ) |

Trend liniowy liczony na wszystkich punktach to tylko jedno z możliwych podejść do zagadnienia oceny trendu. Dodajemy trend lokalnie ważony wielomianami stopnia pierwszego oraz trend liniowy wyznaczony tylko na podstawie median cen w kolejnych miesiącach.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23
|
bwplot(cenam2~dataF, data=mieszkaniaKWW2011Warszawa3Zoliborz, scales=list(y=list(log=T, at=at), x=list(rot=90, cex=0.6)), ylim=c(7000, 14000), main="Zoliborz", panel = function(...) { panel.abline(h=log(at,10), col="grey85") panel.abline(v=4.5 + c(0,12,24,36), col="grey85") panel.bwplot(..., col="grey20") args <- list(...) mod <- rlm(args$y~as.numeric(args$x)) panel.abline(mod) mediany <- sapply(miesiace, function(x) median(args$y[args$x == x], na.rm=T)) mod2 <- rlm(mediany~seq_along(miesiace)) panel.abline(mod2) indtmp <- c(1, length(mediany)) llines(indtmp, mediany[indtmp]) panel.loess(...) } ) |

Mało czytelne są te linie trendu. Narysujemy je grubszą kreską i dodatkowo użyjemy kolorów z pakietu RColorBrewer, które powinny być przyjemniejsze dla oka.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24
|
bwplot(cenam2~dataF, data=mieszkaniaKWW2011Warszawa3Zoliborz, scales=list(y=list(log=T, at=at), x=list(rot=90, cex=0.6)), ylim=c(7000, 14000), main="Zoliborz", panel = function(...) { kolory <- brewer.pal(4, "Set1") panel.abline(h=log(at,10), col="grey85") panel.abline(v=4.5 + c(0,12,24,36), col="grey85") panel.bwplot(..., col="grey20") args <- list(...) mod <- rlm(args$y~as.numeric(args$x)) panel.abline(mod, col=kolory[1], lwd=4, alpha=0.85) mediany <- sapply(miesiace, function(x) median(args$y[args$x == x], na.rm=T)) mod2 <- rlm(mediany~seq_along(miesiace)) panel.abline(mod2, col=kolory[2], lwd=4, alpha=0.85) indtmp <- c(1, length(mediany)) llines(indtmp, mediany[indtmp], col=kolory[4], lwd=4, alpha=0.85) panel.loess(..., col=kolory[3], lwd=4, alpha=0.85) } ) |

Solą w oku są już tylko te niebieskie wykresy pudełkowe, lepiej wyglądałyby one w kolorze szarym, mają być tłem dla linii trendu. Aby zmienić kolor tych punktów używamy funkcji trellis.par.set() i trellis.par.get().
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36
|
bwplot(cenam2~dataF, data=mieszkaniaKWW2011Warszawa3Zoliborz, scales=list(y=list(log=T, at=at), x=list(rot=90, cex=0.6)), ylim=c(7000, 14000), main="Zoliborz", panel = function(...) { tmp <- trellis.par.get("plot.symbol") tmp$pch=19 tmp$col="grey20" tmp$cex=1/2 trellis.par.set("plot.symbol",tmp) tmp <- trellis.par.get("box.rectangle") tmp$col="grey20" trellis.par.set("box.rectangle",tmp) tmp <- trellis.par.get("box.umbrella") tmp$col="grey20" trellis.par.set("box.umbrella",tmp) kolory <- brewer.pal(4, "Set1") panel.abline(h=log(at,10), col="grey85") panel.abline(v=4.5 + c(0,12,24,36), col="grey85") panel.bwplot(..., col="grey20") args <- list(...) mod <- rlm(args$y~as.numeric(args$x)) panel.abline(mod, col=kolory[1], lwd=4, alpha=0.85) mediany <- sapply(miesiace, function(x) median(args$y[args$x == x], na.rm=T)) mod2 <- rlm(mediany~seq_along(miesiace)) panel.abline(mod2, col=kolory[2], lwd=4, alpha=0.85) indtmp <- c(1, length(mediany)) llines(indtmp, mediany[indtmp], col=kolory[4], lwd=4, alpha=0.85) panel.loess(..., col=kolory[3], lwd=4, alpha=0.85) } ) |

Już jest nieźle. Ostatnia modyfikacja wykresu to dodanie opisu wykresu z liczbą procent o który zmieniła się cena mieszkania. Konstruujemy funkcję, która te procenty wyliczy i doklei do nazwy dzielnicy. Przy da się to do kolejnego przykładu, tutaj wyglądać może trochę sztucznie.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81
|
# # ta funkcja doda do nazw dzielnic informacje o zmianie wartosci mieszkan w danej dzielnicy # w badanym okresie czasu uzywajac czterech roznych metod # nazwyIProcenty <- function(mieszkaniaTmp) { # wektor z nazwami dzielnic ze zbioru danych i wektor nazw miesiecy ndzielnic <- levels(mieszkaniaTmp$dzielnica2) miesiace <- levels(mieszkaniaTmp$dataF) # anonimowa funkcja wyznacza cztery liczby # opisujace procentowa zmianę ceny ndzielnic <- paste(ndzielnic," ",sapply(ndzielnic, function(dzielnica) { tmp <- mieszkaniaTmp[mieszkaniaTmp$dzielnica2==dzielnica,] skala <- diff(range(as.numeric(tmp$dataF))) # odporna (robust) regresja liniowa # szukamy trendu liniowego biorąc pod uwagę wszystkie ceny mod <- rlm(tmp$cenam2 ~ I(as.numeric(tmp$dataF)/skala))$coeff proc1 <- round(1000*mod[2]/mod[1]) mediany <- sapply(miesiace, function(miesiac) median(tmp[tmp$dataF == miesiac,"cenam2"])) # ponownie odporna regresja liniowa, # tym razem liczymy wspóczynniki trendu tylko na podstawie median dla kadego miesiąca mod2 <- rlm(mediany~I(seq_along(miesiace)/length(mediany)))$coeff proc2 <- round(1000*mod2[2]/mod2[1]) # zmianę liczymy porównując medianę ceny z wrzenia 2007 z ceną z wrzesnia 2011 # ceny w okresie przejciowym są ignorowane mediany <- na.omit(mediany) proc4 <- round(1000*(mediany[length(mediany)]-mediany[1])/mediany[1]) # regresja lokalnie wygadzana wielomianami pierwszego stopnia # dodatkowe argumenty sa wymienione by zachowac zgodnosc z funkcj panel.loess tmp2 <- loess(tmp$cenam2 ~ I(as.numeric(tmp$dataF)/skala), span = 2/3, degree = 1, family = "symmetric" )$fitted # zmiane liczymy porówując oszacowany trend na wrzesien 2007 # z oszacowanym trendem na wrzesien 2011 tmp2 <- tmp2[order(I(as.numeric(tmp$dataF)/skala))] proc3 <- round(1000*(tmp2[length(tmp2)] - tmp2[1])/tmp2[1]) paste(proc1/10, "% / ", proc2/10, "% / ", proc3/10, "% / ", proc4/10, "% ", sep="") }), sep="") ndzielnic } levels(mieszkaniaKWW2011Warszawa3$dzielnica2) = nazwyIProcenty(mieszkaniaKWW2011Warszawa3) mieszkaniaKWW2011Warszawa3Zoliborz = mieszkaniaKWW2011Warszawa3[which(mieszkaniaKWW2011Warszawa3$dzielnica == "Zoliborz"),] bwplot(cenam2~dataF, data=mieszkaniaKWW2011Warszawa3Zoliborz, scales=list(y=list(log=T, at=at), x=list(rot=90, cex=0.6)), ylim=c(7000, 14000), main=as.character(mieszkaniaKWW2011Warszawa3Zoliborz$dzielnica2)[1], panel = function(...) { tmp <- trellis.par.get("plot.symbol") tmp$pch=19 tmp$col="grey20" tmp$cex=1/2 trellis.par.set("plot.symbol",tmp) tmp <- trellis.par.get("box.rectangle") tmp$col="grey20" trellis.par.set("box.rectangle",tmp) tmp <- trellis.par.get("box.umbrella") tmp$col="grey20" trellis.par.set("box.umbrella",tmp) kolory <- brewer.pal(4, "Set1") panel.abline(h=log(at,10), col="grey85") panel.abline(v=4.5 + c(0,12,24,36), col="grey85") panel.bwplot(..., col="grey20") args <- list(...) mod <- rlm(args$y~as.numeric(args$x)) panel.abline(mod, col=kolory[1], lwd=4, alpha=0.85) mediany <- sapply(miesiace, function(x) median(args$y[args$x == x], na.rm=T)) mod2 <- rlm(mediany~seq_along(miesiace)) panel.abline(mod2, col=kolory[2], lwd=4, alpha=0.85) indtmp <- c(1, length(mediany)) llines(indtmp, mediany[indtmp], col=kolory[4], lwd=4, alpha=0.85) panel.loess(..., col=kolory[3], lwd=4, alpha=0.85) } ) |

Wykres gotowy, w kolejnym odcinku pokażemy jak wygenerować taki wykres dla wszystkich dzielnic.