library(lattice) library(RColorBrewer) library(MASS) 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")) mieszkaniaKWW2011Warszawa <- mieszkaniaKWW2011[mieszkaniaKWW2011$miasto =="Warszawa", ] mieszkaniaKWW2011Warszawa$dzielnica <- factor(mieszkaniaKWW2011Warszawa$dzielnica) mieszkaniaKWW2011Warszawa <- mieszkaniaKWW2011Warszawa[mieszkaniaKWW2011Warszawa$powierzchnia<300, ] dataStart <- as.Date("13-09-2007","%d-%m-%Y") mieszkaniaKWW2011Warszawa2 <- mieszkaniaKWW2011Warszawa mieszkaniaKWW2011Warszawa2 <- mieszkaniaKWW2011Warszawa2[mieszkaniaKWW2011Warszawa2$data >= dataStart,] mieszkaniaKWW2011Warszawa2$dataF <- factor(as.character(mieszkaniaKWW2011Warszawa2$data, "%Y-%m")) # # mieszkaniaKWW2011Warszawa3 # tylko srednie 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"),] mieszkaniaKWW2011Warszawa3$dzielnica2 <- reorder(mieszkaniaKWW2011Warszawa3$dzielnica, mieszkaniaKWW2011Warszawa3$cenam2, median) mieszkaniaKWW2011Warszawa3Zoliborz <- mieszkaniaKWW2011Warszawa3[which(mieszkaniaKWW2011Warszawa3$dzielnica == "Zoliborz"),] at = seq(1000,24000,1000) png("mieszkaniaWarszawaZoliborzPudelka1.png",600,600) bwplot(cenam2~dataF, data=mieszkaniaKWW2011Warszawa3Zoliborz, main="Zoliborz") dev.off() png("mieszkaniaWarszawaZoliborzPudelka2.png",600,600) bwplot(cenam2~dataF, data=mieszkaniaKWW2011Warszawa3Zoliborz, scales=list(y=list(log=T, at=at), x=list(rot=90, cex=0.6)), main="Zoliborz") dev.off() png("mieszkaniaWarszawaZoliborzPudelka3.png",600,600) 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") dev.off() png("mieszkaniaWarszawaZoliborzPudelka4.png",600,600) 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) } ) dev.off() miesiace <- levels(factor(mieszkaniaKWW2011Warszawa3Zoliborz$dataF)) png("mieszkaniaWarszawaZoliborzPudelka5.png",600,600) 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(...) } ) dev.off() png("mieszkaniaWarszawaZoliborzPudelka6.png",600,600) 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) } ) dev.off() png("mieszkaniaWarszawaZoliborzPudelka7.png",600,600) 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) } ) dev.off() # # 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"),] png("mieszkaniaWarszawaZoliborzPudelka8.png",600,600) 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) } ) dev.off() png("mieszkaniaWarszawaZoliborzGrupa1.png",1400,1000) bwplot(cenam2~dataF|dzielnica2, data=mieszkaniaKWW2011Warszawa3, scales=list(y=list(log=T, at=at), x=list(rot=90, cex=0.6)), ylim=c(5000, 20000), 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) } ) dev.off() # # usun male dzielnice usun <- names(which(table(mieszkaniaKWW2011Warszawa3$dzielnica)<1000)) mieszkaniaKWW2011Warszawa3 <- mieszkaniaKWW2011Warszawa3[!(mieszkaniaKWW2011Warszawa3$dzielnica %in% usun),] mieszkaniaKWW2011Warszawa3$dzielnica <- factor(mieszkaniaKWW2011Warszawa3$dzielnica) # # wszystkie dzielnice png("mieszkaniaWarszawaZoliborzGrupa2.png",1400,1000) bwplot(cenam2~dataF|dzielnica2, data=mieszkaniaKWW2011Warszawa3, scales=list(y=list(log=T, at=at), x=list(rot=90, cex=0.6)), ylim=c(5000, 20000), 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) } ) dev.off()