wykres rozrzutu – SmarterPoland.pl

Ankiety studenckie a głaski dla pracowników dydaktycznych

W poprzednim wpisie dotyczącym wyników z ankiet studenckich zauważyliśmy, że w percepcji studenta (ocenianej przez korelację reszt) pytanie o ocenę przedmiotu jest czymś innym niż pozostałe pytania, dotyczące raczej sposobu prowadzenia zajęć.
Zamiast więc pokazywać wyniki ankiet we wszystkich ocenianych 11 wymiarach, skupię się na razie na dwóch, mianowicie ogólnej ocenie zajęć i ogólnej ocenie prowadzącego.

Jak juz pisałem, raporty opracowane przez studentów (patrz tutaj) są bardzo miłe dla oka, ale mam z nimi jeden szkopuł.
Mianowicie znaczna część tych raportów skupia się na pokazaniu par prowadzący-przedmiot, których oceniło przynajmniej 5 osób i który mieli ogólną ocenę w najwyższym kwartylu (takie top 25%).

Poniżej będę argumentował, że zarówno kryterium wyboru tych par do pokazania jak i sposób ich pokazania nie jest najlepszy.

Co jest złego w kryterium top 25% najlepszych z ocenionych przez przynajmniej 5 osób?

Po co przyjmować próg 5 osób? Jeżeli para prowadzący/przedmiot był oceniony przez mniej niż 5 osób to można argumentować, że taka ocena jest przypadkowa i nie jest reprezentatywna dla tego przedmiotu. Ten argument jednak ciężko obronić jeżeli mamy specjalistyczny kurs na który zapisało się tylko 4 osoby i wszystkie wypełniły ankiety.Można argumentować, że mała próbka spowoduje, że dobry prowadzący zostanie źle oceniony przypadkowo bo wypowie się tylko niewielka (niereprezentatywna) część studentów. Nie chcemy krzywdzić prowadzących przypadkowo złymi ocenami. Ale ten argument się też nie broni ponieważ pokazywanych jest tych 25% najlepszych. Nie znajduję uzasadnienia dla progu przynajmniej pięciu ankiet. Jestem może mało obiektywny, ponieważ prowadzę często przedmioty specjalistyczne (inna nazwa na niepopularne) na których liczba zapisanych osób waha się od 4 do 8 ;-).
Po co przyjmować próg top 25%? Jak rozumiem celem ankiety jest możliwość wskazania (choć nie jest jasne jeszcze komu, ale o tym później) dobrych dydaktyków. Ale i w grupie samych wybitnych dydaktyków i w grupie samych beznadziejnych zawsze da się wskazać top 25%. Lepszy wydaje się być bezwzględny próg. Jeżeli uznać średnią ocenę 5,5 albo 6 jako wynik ponadprzeciętny (w skali 1-7) to sensowniejsze było by pokazywanie wszystkich ponadprzeciętnych dydaktyków.
Dlaczego nie pokazywać wyników w rozbiciu na kursy. Wracamy do pytania, że nie jest jasne dla kogo przygotowane są te wyniki ankiet. Gdy byłem studentem, chciałem by ankiety były przeprowadzane po to by ułatwić studentom wybór ćwiczeniowca lub prowadzącego laboratorium. Jeżeli tak to dla każdego kursu chciałbym mieć wskazane kto jest najlepszym ćwiczeniowcem/laborantem i czy są pomiędzy ćwiczeniowcami duże różnice. Będąc prowadzącym uważam, że takie ankiety mogłyby mi pomóc ocenić, z drugiej strony barykad,y jakie są prowadzone przeze mnie zajęcia. Jak to jednak zrobić kiedy prowadząc ćwiczenia ze statystyki mogę moją średnia z ankiet mogę porównać ze średnią z wszystkich kursów. Wolałabym już porównanie ze średnią z wszystkich prowadzących statystykę.

Tak więc aktualny sposób prezentacji wyników wygląda trochę jak sposób w jaki firmy motywują sprzedawców, pokazując top 25% osób mających największe obroty w każdej kategorii sprzedawanych przedmiotów. I jest to raczej bonus łechczący ego osób, które znalazły się na tej liście, ale bonus z którego trudno wyciągnąć użyteczne informacje. (Oczywiście taki bonus jest bardzo ważny. Na niewiele więcej może liczyć starający się dydaktyk jak na wysokie miejsce w rankingach popularności ankietowanych studentów, ale to inna historia.)

Ok, ponarzekałem, czas na częśc konstruktywną.
Pierwszy pomysł to pokazanie każdej pary przedmiot-prowadzący na osi średnia ocena przedmiotu/średnia ocena prowadzącego. Dodatkowo zaznaczmy wielskością punktu jak duża jest grupa studentów daną parę oceniała. Po prawej stronie i u góry wykresu zaznaczyłem dodatkowo decyle (dziesiątki centyli) rozkładu tych ocen, aby łatwiej można było odczytać jaka część kurso/prowadzących miała wyniki lepsze niż k.

Ponieważ nie mam zgody na prezentowanie nazwisk prowadzących więc na poniższych wykresach jedyny punkt oznaczony czerwonym kolorem to wyniki ankiet z mojego kursu ,,Modele liniowe i mieszane”, który prowadziłem w semestrze zimowym 2010 (to jedyny kurs który prowadziłem w tamtym semestrze). Cytując za Maciejem Stuhrem ,,wstydu nie ma”, od siebie mogę dodać że szału też, ale czerwoną kropkę dorysowałem.

Nawet bez nazwisk prowadzących z takiego wykresu można odczytać kilka informacji. (Samorządzie, zdobądź zgode na pokazanie nazwisk dydaktyków o najlepszych ocenach)

Prowadzący są zazwyczaj lepiej oceniani niż przedmioty które prowadzą (cóż, dla mnie to zaskoczenie, właściwie jak to czytać, czy studenci chcieliby innych przedmiotów?).
Zdecydowana większość prowadzących/przedmiotów jest dobrze oceniana (na użytek tego postu uważam, że dobrze to średnia ocena 5 lub więcej).
Będąc prowadzącym i widząc oceny moich kursów moge sobie z tego wykresu odczytać gdzie się znajduję w rozkładzie wszystkich ocenianych przedmioto-prowadzących. Wciąż wolałbym taki wykres w rozbiciu na kursy, ale przynajmniej mogę więcej porównać niż tylko średnie moją i globalną.

Wszelkie uwagi co pokazać dodatkowo/inaczej mile widziane.

R, Kair, Cairo, wiek a liczba lat nauki

Dziś będzie o wielu tematach jednocześnie.
Po pierwsze przeglądając materiały z konferencji useR znalazłem informacje o pakiecie Cairo. Wstyd, że o tym pakiecie dowiedziałem się tak późno. Ale ciesze się że się dowiedziałem wystarczająco by napisać dlaczego.

Cairo to biblioteka do grafiki 2d umożliwiająca zapisywanie grafiki do różnych formatów, między innymi wektorowych formatów PS, PDF, SVG, ale również do rastrowych formatów PNG itp. Biblioteka pozwala na stosowanie takich miłych technik jak anty-aliasing itp. poprawiających wygląd grafiki, szczególnie rastrowej, szczególnie w porównaniu z tym co R produkuje domyślnie. Zobaczmy czy biblioteka ta odmieni rysunki na tym blogu. Zaczniemy od prostego przykładu, w którym porównamy liczbę lat nauki z liczbą przeżytych lat, oczywiście bazując na zbiorze danych o Diagnozie społecznej.

Aby nie zaciemniać tego co najważniejsze, nie będę tutaj wklejał kodu R, który generuje wykres. Osoby zainteresowane znajdą ten kod tutaj. W kodzie tym wykorzystujemy funkcję xyplot() z pakietu lattice do wygenerowania obiektu wyk opisującego wykres.

Porównamy wygląd grafiki zapisanej poleceniami png(grDevices), CairoPNG(Cairo) i CairoSVG(Cairo). W przypadku tej ostatniej funkcji rozmiary podajemy w calach, w poprzednich dwóch w pixelach.

# standardowa funkcja pakietu R
png("wiek_lnauki.png",600,600)
print(wyk)
dev.off()
# ten sam format png, ale z użyciem pakietu Cairo
CairoPNG("wiek_lnauki_cairo.png",600,600)
print(wyk)
dev.off()
# i tym razem grafika wektorowa SVG
CairoSVG("wiek_lnauki.svg",7,7)
print(wyk)
dev.off()

Poniżej wygenerowane pliki. Pierwszy w formacie png wygenerowany przez funkcję png.

Drugi w formacie png ale wygenerowany biblioteka Cairo.

I trzeci w formacie SVG, nie każda przeglądarka potrafi go wyświetlić, jeżeli poniżej nie widać obrazka to znaczy że trzeba zmienić przeglądarkę.

Czy widać różnice? Ogromne, szczególnie na krzywych które w pierwszym przypadku są niesamowicie spixelowane w drugim przypadku już znacznie gładsze. Widać to szczególnie w dużym powiększeniu. W dużym powiększeniu widać też zalety wektorowego formatu SVG, krzywe będą gładkie bez względu jak bardzo je powiększymy. Oczywiście za format wektorowy trzeba zapłacić. W zbiorze danych Diagnoza Społeczna znajdują się dane dla kilku tysięcy osób. Gdyby każdą z nich zaznaczyć punktem to wektorowy format w którym każdy z tych punktów byłby opisany, zajmowałby kilkanaście MB. Nie najlepiej jak na grafikę do umieszczenia w internecie. Dlatego też na trzecim z powyższych obrazków są tylko krzywe bez punktów.

Dyskusje o technikaliach mamy już za sobą, zobaczmy co w ogóle widać na tych wykresach. Porównujemy liczbę lat nauki versus wiek. Dane bierzemy zarówno dla osób ankietowanych w roku 2000 jak i 2011. Znaczna część osób uczestniczy w obu badaniach, dla tej grupy spodziewamy się, że lat przybędzie ale lat nauki niekoniecznie.
Ciągła linią zaznaczyłem wygładzona medianę, kropkowaną linią wygładzony kwantyl rzędu 90%.

Dla osób które obecnie są w wieku 40-60 lat, mediana liczby lat edukacji to 11, dla osób w wieku 25 lat ta mediana jest już o 2 lata wyższa, efekt coraz większej liczby osób studiujących i też zmian w systemie edukacji. Niższą medianę liczby lat edukacji u osób w wieku >70 lat można wytłumaczyć wojną.

Przyjrzyjmy się jeszcze kwantylowi rzędu 90%. Dla 30latków kwantyl ten wynosi 18 lat nauki, czyli ponad 10% dzisiejszych 30latków ma na koncie 18 lub więcej lat nauki. To sporo, ale dla niektórych uczenie się nigdy się nie kończy. W stosunku do stanu sprzed 11 lat coraz więcej lat spędzamy na edukacji. Czy jako społeczeństwo jesteśmy dzięki temu mądrzejsi?

Zmiana dochodów w ostatnich 8 latach w podziale na płeć

Analiza zróżnicowania dochodów bardzo mnie interesuje. To jednak większy temat i poświęcimy mu więcej czasu kiedy indziej. Dziś chciałbym podzielić się jednym wykresem, na którym zobaczymy jak wyglądał rozkład dochodu netto ankietowanego w zależności od wieku i w zależności od płci, tak w roku 2003 jak i w roku 2011. Oczywiście złotówka dzisiaj i wtedy to dwie różne złotówki, nie sposób przeliczyć siły nabywczej teraz i wtedy bo zależy ona od koszyka zakupów. Inaczej wyglądała zmiana cen artykułów luksusowych, inaczej zmiana cen jedzenia a inaczej zmiana cen mieszkań. W każdym razie na potrzeby tego rysunku wykorzystam wskaźnik inflacji dla tych 8% wynoszący 23.7% (na podstawie http://blog.opiekuninwestora.pl/index.php/inflacja/). Więc przedstawiane ceny to albo ceny podane w roku 2011, albo ceny podane w roku 2003 i skorygowane o inflacje.

W poniższych ilustracjach nie analizujemy osobno osób żyjących samotnie (tzw singli), par bez dzieci (DINKs) i par z dziećmi, choć pewnie pomiędzy tymi trzema grupami są znaczne różnice, ale nie wszystko na raz.

Ciągła linia to wygładzona ocena mediany, przerywana to kwantyl rzędu 90%.
Dużych niespodzianek nie ma, ale niektóre wyniki są ciekawe.
Po pierwsze, nawet po uwzględnieniu inflacji zarobki wzrosły (nie tylko ankietowanych, bo to jest oczekiwane, ale tez rozkład zarobków w grupie wiekowej np 30 latków). I to wzrosły znacznie. Oby związane to było z większa produktywnością a nie spadkiem wartości złotówki.
Po drugie kobiety zarabiają mniej. Dotyczy to i medianowych zarobków i kwantyla 90%.
Po trzecie, i chyba najciekawsze, w roku 2011 wyraźnie widać, że najwyższe dochody uzyskują osoby w wieku 30-40 lat.
Biorąc pod uwagę, że oś oY jest logarytmiczna, zarówno osoby młodsze jak i starsze zarabiają znacznie mniej. U kobiet ten okres wyższych dochodów kończy się wcześniej niż u mężczyzn, ale u obu płci kończy się czy to patrząc na medianę czy na kwantyl 90%.
Inaczej było w roku 2003. Nie było wtedy tak dużych dysproporcji związanych z wiekiem, szczególnie jeżeli porównywać mediany dochodów. Ciekawe czy to zróżnicowanie będzie się jeszcze pogłębiało.

To pierwsza ilustracja. Czas na dokładniejsze drążenie tematu. Kierunków badania czynników różnicujących wysokość dochodów jest dużo. Może macie propozycje od czego zacząć?

Rysujemy rozkład cen krok po kroku, część 2

Kontynuując temat z wczoraj, narysujemy rozkład cen mieszkania dla każdej z dzielnic Warszawy.

Punktem wyjścia jest przygotowanie danych, wykonajmy pierwsze 32 linie tak jak w poprzednim wpisie.

Aby wyświetlić na rożnych panelach dane dla kolejnych dzielnic, wystarczy zmodyfikowac formułę na cenam2~data|dzielnica, oraz za zabiór danych wskazać mieszkaniaKWW2011Warszawa2.

at = seq(1000,24000,1000)
 
xyplot(cenam2~data|dzielnica, group=rozmiar,
        data=mieszkaniaKWW2011Warszawa2,
        scales=list(y=list(log=T, at=at)), ylim=c(6000, 16000),
        type=c("smooth"), pch='.', lwd=3,
        auto.key=list(space="top", columns=3, pch=19),
        panel=function(...) {
  panel.abline(h=log(at,10), col="grey85")
  panel.abline(v=365.25 * (38:41), col="grey85")
  panel.xyplot(list(...)$x, list(...)$y, type="p",col="grey", cex=1/4)
  panel.xyplot(list(...)$x, list(...)$y, type="r",col="black", lwd=2,lty=2)
  panel.xyplot(...)
})

Kolejne panele przedstawiają kolejne dzielnice, ale ich kolejność jest alfabetyczna. Taka sama jak kolejność poziomów zmiennej czynnikowej dzielnica. Nie zawsze kolejność alfabetyczna będzie najlepsza. Użyjemy funkcji reorder by zmienić kolejność poziomów tak by odpowiadała medianie ceny metra kwadratowego w danej dzielnicy. Kod generujący obrazek będzie taki sam, zmieni się tylko kolejność dzielnic.

mieszkaniaKWW2011Warszawa2$dzielnica2 <- reorder(mieszkaniaKWW2011Warszawa2$dzielnica, mieszkaniaKWW2011Warszawa2$cenam2, median)
 
xyplot(cenam2~data|dzielnica2, group=rozmiar,
        data=mieszkaniaKWW2011Warszawa2,
        scales=list(y=list(log=T, at=at)), ylim=c(6000, 16000),
        type=c("smooth"), pch='.', lwd=3,
        auto.key=list(space="top", columns=3, pch=19),
        panel=function(...) {
  panel.abline(h=log(at,10), col="grey85")
  panel.abline(v=365.25 * (38:41), col="grey85")
  panel.xyplot(list(...)$x, list(...)$y, type="p",col="grey", cex=1/4)
  panel.xyplot(list(...)$x, list(...)$y, type="r",col="black", lwd=2,lty=2)
  panel.xyplot(...)
})

Dla niektórych dzielnic jest mało punktów, co powoduje, że trudno mieć zaufanie do wyznaczonego trendu. Tym razem usuniemy te dzielnice, dla których nie ma przynajmniej 2000 wierszy. Kod generujący wykres jest bez zmian, usuwamy tylko obserwacje z dzielnic w których obserwacji było mniej niż 2k.

usun <- names(which(table(mieszkaniaKWW2011Warszawa2$dzielnica)<2000))
mieszkaniaKWW2011Warszawa2 <- mieszkaniaKWW2011Warszawa2[!(mieszkaniaKWW2011Warszawa2$dzielnica %in% usun),]
mieszkaniaKWW2011Warszawa2$dzielnica2 <- factor(mieszkaniaKWW2011Warszawa2$dzielnica2)
 
 
xyplot(cenam2~data|dzielnica2, group=rozmiar,
        data=mieszkaniaKWW2011Warszawa2,
        scales=list(y=list(log=T, at=at)), ylim=c(6000, 16000),
        main="",
        type=c("smooth"), pch='.', lwd=3,
        auto.key=list(space="top", columns=3, pch=19),
        panel=function(...) {
  panel.abline(h=log(at,10), col="grey85")
  panel.abline(v=365.25 * (38:41), col="grey85")
  panel.xyplot(list(...)$x, list(...)$y, type="p",col="grey", cex=1/4)
  panel.xyplot(list(...)$x, list(...)$y, type="r",col="black", lwd=2,lty=2)
  panel.xyplot(...)
})

Rysujemy rozkład cen krok po kroku

Dostałem maila od Grzegorza L. z pytaniem jak zrobić krok po kroku wykresy z wpisów http://smarterpoland.pl/?p=169 i http://smarterpoland.pl/?p=172. W sumie cztery wykresy, każdy z nich postaram się opisać jak powstawał krok po kroku. Dziś zajmiemy się pierwszym z nich czyli zmiana cen mieszkań w dzielnicy Zoliborz w rożnych grupach wielkości mieszkań.

Zaczniemy od opisu jak dane były przygotowane a następnie pokażemy jak zrobić wykres krok po kroku. Skrypt R z wszystkimi komendami poniżej opisanymi znajduje się tutaj.

Część 1: Przygotowanie danych
Poniższy skrypt wczytuje do R dane odczytując je bezpośrednio z Internetu, następnie wybiera tylko dane dla Warszawy, dla dzielnicy Zoliborz. Wybiera tylko mieszkania o powierzchni do 300m2 i zgłoszone w ostatnich 3 latach. Następnie dodajemy do danych nowa zmienna opisującą wielkość mieszkania w jednej z trzech grup, małe do 49 m2, średnie od 49 do 68m2 i duże, o powierzchni ponad 68m2. Następnie dodatkowo wybieramy tylko podzbiór mieszkań z Zoliborza.

library(lattice)
#
# wczytujemy dane o mieszkaniach bezposrednio z internetu
# moze chwile potrwac
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 mieszkania tylko z Warszawy
# usuwamy informacje o dzielnicach z innych miast (funkcja factor to zalatwi)
# wybrane mieszkania zapisujemy do zmiennej mieszkaniaKWW2011Warszawa
#
mieszkaniaKWW2011Warszawa <- mieszkaniaKWW2011[mieszkaniaKWW2011$miasto =="Warszawa", ]
mieszkaniaKWW2011Warszawa$dzielnica <- factor(mieszkaniaKWW2011Warszawa$dzielnica)
 
#
# zawezimy horyzont do ofert sprzedazy z ostatnich 4 lat
# pominiemy oferty sprzed 13 wrzesnia 2007
# usuniemy tez mieszkania o powierzchni powyzej 300m2 (jest ich niewiele a psuja niektore wizualizacje)
#
dataStart <- as.Date("13-09-2007","%d-%m-%Y")
mieszkaniaKWW2011Warszawa2 <- mieszkaniaKWW2011Warszawa
mieszkaniaKWW2011Warszawa2 <- mieszkaniaKWW2011Warszawa2[mieszkaniaKWW2011Warszawa2$data >= dataStart,]
mieszkaniaKWW2011Warszawa2 <- mieszkaniaKWW2011Warszawa2[mieszkaniaKWW2011Warszawa2$powierzchnia<300, ]
 
#
# powierzchnie mieszkania zamienimy na zmienna jakosciowa
# wystepujaca na trzech poziomach, do 49m2, od 49 do 68 m2, powyzej 68m2
#
mieszkaniaKWW2011Warszawa2$rozmiar <- cut(mieszkaniaKWW2011Warszawa2$powierzchnia, c(0,49,68,200), labels=c("do 49 m2", "od 49 do 68 m2", "powyzej 68 m2"))
 
#
# wybierzemy tylko mieszkania z jednej dzielnicy, Zoliborz
#
mieszkaniaKWW2011Warszawa2Zoliborz <- mieszkaniaKWW2011Warszawa2[which(mieszkaniaKWW2011Warszawa2$dzielnica == "Zoliborz"),]

Część 2: Generowanie obrazków.

Ok, mamy przygotowane dane, teraz czas na zrobienie rysunku. Zaczniemy od prostego użycia funkcji xyplot{lattice}.

#
# pierwsze podejscie
# podajemy zmienne ktore maja byc narysowane na osi X i Y
# podajemy zmienna grupujaca i tytul wykresu
#
xyplot(cenam2~data, group=rozmiar,
        data=mieszkaniaKWW2011Warszawa2Zoliborz,
        main="Zoliborz")

Miło ze strony funkcji xyplot(), że jest w stanie narysować na osi OX zmienna typu Date, ale sam wykres nie powala na kolana.
Za bardzo w oczy rzucają się punkty, których jest tak dużo że nic ciekawego nie widać. Zmniejszymy więc wielkość punktów zaznaczając je jednopikselową kropka (pch=”.”), dodamy linie pomocnicze siatki i wygładzoną linię trendu (type=”g” i „smooth”).

#
# dodajemy pomocnicze linie siatki i wygladzona linie trendu
# zmieniamy tez sposob rysowania punktow na jednopunktowe kropki
#
xyplot(cenam2~data, group=rozmiar,
        data=mieszkaniaKWW2011Warszawa2Zoliborz,
        main="Zoliborz", type=c("p","g","smooth"), pch=".")

Już trochę lepiej, ale wciąż niewiele widać. Skala OY jest liniowa a dla prawoskośnych zmiennych dobrym pomysłem jest skala logarytmiczna. Dodajemy argument scales, który pozwala na zmianę osi OY na logarytmiczną.

#
# zmieniamy os OY na logarytmiczna
#
xyplot(cenam2~data, group=rozmiar,
        data=mieszkaniaKWW2011Warszawa2Zoliborz,
        scales=list(y=list(log=T)),
        main="Zoliborz", type=c("p","g","smooth"), pch=".")

Co nam się nie podoba tym razem? Wiele rzeczy. Po pierwsze przydała by się legenda. Dodamy ją argumentem auto.key, przy okazji określimy gdzie i jak ma być legenda rysowana, choć wcale nie musimy tego określać, wystarczyłoby auto.key=T. Druga rzecz to widoczność linii trendu w każdej z grup. Są mało widoczne z uwagi na grubość linii i z uwagi na duży rozstrzał na osi OY. Zmienimy i to i to, podamy grubość linii lwd=3 i ustalimy zakres zmienności na osi OY na 6-16k/m2. ostatnia poważna zmienna będzie dotyczyła znaczników na osi OY. Poprzednio wyglądały one słabo, jakieś potęgi 10 (ile to 10^4.2?), tym razem jawnie wskażemy w którym punktach mają pojawiać się znaczniki osi podając argument at.

#
# zmieniamy zakres osi OY
# jawnie wskazujemy gdzie maja byc rysowane znaczniki na osi OY
# dodajemy tez legende umieszczona na gorze wykresu rysowana poziomo (w trzech kolumnach, ale dla trzech wartosci wychodzi na jedno)
#
at = seq(1000,24000,1000)
xyplot(cenam2~data, group=rozmiar,
        data=mieszkaniaKWW2011Warszawa2Zoliborz,
        scales=list(y=list(log=T, at=at)), ylim=c(6000, 16000),
        main="Zoliborz",
        type=c("p", "g", "smooth"), pch='.', lwd=3,
        auto.key=list(space="top", columns=3, pch=19))

Wykres nabiera kształtów. Nadpiszemy teraz funkcję rysującą panel, w taki sposób by dodatkowo rysować trend wyznaczony metodą regresji liniowej. Rysujemy go na czarno przerywaną linią wykorzystując funkcję panel.xyplot().

#
# zmieniamy funkcje rysujaca panel, teraz ta funkcja dorysowuje krzywa regresji liniowej 
# przerywana linia
#
xyplot(cenam2~data, group=rozmiar,
        data=mieszkaniaKWW2011Warszawa2Zoliborz,
        scales=list(y=list(log=T, at=at)), ylim=c(6000, 16000),
        main="Zoliborz",
        type=c("p", "g", "smooth"), pch='.', lwd=3,
        auto.key=list(space="top", columns=3, pch=19),
        panel=function(...) {
  panel.xyplot(list(...)$x, list(...)$y, type="r",col="black", lwd=2,lty=2)
  panel.xyplot(...)
})

Czas na ostatnie modyfikacje. Co zmienimy? Po pierwsze punkty w legendzie są rysowane pustymi okręgami, a lepiej wyglądać będą wypełnione. Poniżej używam funkcji trellis.par.set() by o zmienić.
Po drugie, porzednie linie pomocnicze siatki nie były rysowane w tych samych miejscach co znaczniki na osi OY, usuwam więc argument type=”g” a zamiast tego w funkcji rysującej panel ręcznie dorysowuje linie pomocnicze w odpowiednich miejscach i pionowo i poziomo funkcją panel.abline(). Po trzecie zmieniam sposób rysowania punktów. Używam do tego funkcji panel.xyplot() z type=”p” i cex=1/4 dzięki czemu każda oferta sprzedaży jest zaznaczona jednym niewielkim szarym punktem.

#
# funkcja trellis.par.set() zmieniamy ksztalt io kolor kropek w legendzie
# z argumentu type usuwamy rysowanie punktow i lini pomocnicznych
# z tego tez powodu dodajemy do funkcji rysujacej panel kilka dodatkowych instrukcji
# rysujacych punkty, linie regresji liniowej i linie pomocnicze siatki
#
tt = trellis.par.get("superpose.symbol")
tt$pch=19
tt$col = c( "#0080ff", "#ff00ff", "darkgreen")
trellis.par.set("superpose.symbol",tt)
 
xyplot(cenam2~data, group=rozmiar,
        data=mieszkaniaKWW2011Warszawa2Zoliborz,
        scales=list(y=list(log=T, at=at)), ylim=c(6000, 16000),
        main="Zoliborz",
        type=c("smooth"), pch='.', lwd=3,
        auto.key=list(space="top", columns=3, pch=19),
        panel=function(...) {
  panel.abline(h=log(at,10), col="grey85")
  panel.abline(v=365.25 * (38:41), col="grey85")
  panel.xyplot(list(...)$x, list(...)$y, type="p",col="grey", cex=1/4)
  panel.xyplot(list(...)$x, list(...)$y, type="r",col="black", lwd=2,lty=2)
  panel.xyplot(...)
})

I koniec.

Jeszcze mała autoreklama. Więcej o pakiecie lattice można przeczytać w drugim wydaniu ,,Przewodnika po pakiecie R” wydanego przez wydawnictwo GiS w tym roku (jeden podrozdział dostępny w Internecie) lub w książce Lattice Multivariate Data Visualization with R, w całości poświęconej temu pakietowi (zobacz tutaj). Oczywiście można też znaleźć sporo informacji korzystając po prostu z googla.

Czym się różnią ceny mieszkań na Żoliborzu od cen mieszkań na Bemowie

Dzisiejszy wpis to przygotowanie gruntu do jutrzejszego, w ktorym rozliczymy sie ze zmianami cen mieszkan. Od jakiegos czasu na różnych stronach można znaleźć artykuły o nazwach ”realne ceny mieszkań spadły o X %’’ (dobrze że realne, co to by byly gdyby spadaly ceny nierealne). Celem tego i kolejnego wpisu jest pokazanie że zmianę ceny można liczyć na wiele różnych sposobów otrzymując wiele różnych wyników. Więc podanie zmiany nie mówiąc dokładnie jak była ona liczona to zwykłą propaganda.

Będziemy korzystać ze zbioru danych mieszkaniaKWW2011, przeanalizujemy tylko ceny z miasta Warszawa dla mieszkan o powierzchni do 300m2. Na początek przyjrzymy się dwóm dzielnicom. Jednej w której mieszkam (Bemowo) i jednej w której ceny się ciekawie zachowują (Żoliborz).

Przedstawimy zmianę trendu ceny mieszkań w czasie w rozbiciu na trzy grupy wielkości mieszkań. Kwantyle rzędu 1/3 i 2/3 z rozkładu wielkości mieszkań wynoszą 49m2 i 68m2, wiec będziemy analizować cenę m2 w grupie mieszkań w trzech grupach: do 49m2, w grupie mieszkań dużych o powierzchni powyżej 68m2 i w grupie pozostałych – średnich co do wielkości mieszkań.

Na poniższym wykresie każdy punkt to jedna oferta sprzedaży. Linia zielona, różowa i niebieska oznaczają lokalne wygładzenie wielomianowe dla ceny m2 dla różnych grup wielkości mieszkań, a czarna linia przerywana to globalny trend w cenach mieszkań bez uwzględnienia struktury mieszkań. Struktura jest istotna ponieważ hipotetycznie, jeżeli mieszkania duże maja niższa cene za m2, i przed czterema laty mieszkań dużych bylo sprzedawanych mniej niż teraz to zaobserwujemy zmiane w średniej cenie nie ze względu na rzeczywistą zmianę ceny ale ze wzgledu na zmianę struktury sprzedawanych mieszkań. Podobnie z analizą cen w Warszawie. Jeżeli w tym roku sprzedaje sie więcej mieszkań na obrzeżach Warszawy niż przed czterema laty (z roku na rok buduje się średnio coraz dalej bo tam jeszcze jest miejsce) to zmiana średniej ceny będzie zwiazana z tym ze w koszyku jest coraz więcej mieszkań dalekich od centrum a wiec tańszych.

No dobrze, to tyle tytułem teorii a teraz obrazki dla dwoch obiecanych dzielnic.

Dla Żoliborza jak widzimy ceny m2 dla mieszkań dużych sa wyzsze niz dla malych, moze byc to zwiazane ze duze mieszkania to juz apartamenty dlatego cenę ich m2 ciezko porównywać z cena m2 mieszkania o standardowym wykonczeniu. Pdobnie beda zachowywaly sie mieszkania w Śródmieściu. Te apartamenty tez najbardziej straciły na wartości. Na Bemowie jest inaczej. Nie buduje sie raczej apartamentow w dzielnicy sypialni wiec tutaj to male mieszkania maja wyższą cenę za m2. W obu przypadkach ceny maja tendencje do spadania ale w kazdej kategorii tempo spadku jest inne.

Poniżej prezentujemy wykresy dla wszystkich dzielnic Warszawy.

Mapa wartości samochodów

Wczoraj graficznie przedstawialiśmy jak wygląda rozkład cen w zależności od roku produkcji. Co ciekawe w skali logliniowej zależność pomiędzy ceną a rokiem produkcji jest bardzo bliska liniowej dla wielu marek. Zakładając zależność liniową możemy wyznaczyć dwie liczby dla każdej marki. Średnią procentową zmianę ceny w ciągu roku i orientacyjną średnią cenę auta wyprodukowanego w roku 2006. Rok 2006 został wybrany by móc porównywać ceny różnych marek ze sobą, oczekiwaną cenę w roku 2006 wyznaczmy bez względu na to czy dany samochód był produkowany w tym roku czy nie (wykorzystamy liniową zależność aby oszacować średnią cenę w roku 2006).
Na poniższej mapie zaznaczono każdą markę w układzie dwóch współrzędnych, utraty wartości w ciągu roku oraz orientacyjnej ceny w roku 2006. Należy zaznaczyć, że wartość procentowa na osi OX to różnica pomiędzy średnimi cenami dwóch kolejnych roczników, a nie różnica pomiędzy ceną danego rocznika w dwóch kolejnych latach. Tą drugą wartość będzie można ustalić za rok.

Dodatkowo na tej mapie wielkością punktu zaznaczono jak wiele aut było wykorzystanych do oszacowania obu wskaźników (cztery wielkości do 100, od 100 do 200, od 200 do 500, powyżej 500) oraz z jakiego kraju pochodzi dana marka. Jeżeli chodzi o pochodzenie to z uwagi na przejęcia oraz przenoszenie produkcji do innych krajów należy tę informację traktować z dystansem.

To co ciekawego widać na tej mapie to że są auta tanie, które nie tracą wiele z czasem na wartości (np Fiat Uno), drogie, które też dużo na wartości tracą (Porche Cayenne), tanie które dużo tracą (Dacia Logan), drogie które mało tracą (Volkswagen Bora) i całe spektrum innych możliwości.

Aby łatwiej ten wykres się analizowało możemy dokonać zbliżenia

Volkswagen Passat traci średni 16% na rok i kosztuje z roku 2006 około 40tys. W otoczeniu znaleźć można np Toyote Corolle Verso. Na prawo znajdziemy Forda S-Maxa o podobnej cenie ale za to znacznie większej utracie na wartości.

Zobaczmy co się stanie jeżeli zagregujemy modele z tej samej marki.

Otrzymujemy mapę marek, im bardziej na lewo tej mapy tym mniejsza utrata na wartości im niżej tym tańsze auta.

Zakup kontrolowany, czyli wybieramy auto dla rodziny

Od prawie czterech lat szukam rodzinnego auta. Ostatnio rozmawiałem ze znajomym, który był zdziwiony dlaczego to taki trudny wybór. Jak to stwierdził wystarczy wprowadzić dane do komputera, określić funkcje celu i wybrać najlepszą opcję. Hmmm, może i tak.

Mamy już dane zebrane z serwisu otomoto.pl, więc spróbujmy zobaczyć jak wyglądają rozkłady cen w czasie różnych marek. Informacje o zmianach cen znaleźć można na wielu serwisach, np. autocentrum.pl, ale tutaj przedstawimy te zmiany lepiej i czytelniej. Tak jest, konkurujemy z eutotax.

Aktualnie wysoko na liście rozważanych marek stoi Passat. Poniższe wykresy będą wykonane dla prawie każdej marki obecnej w zbiorze danych cenyAut2011. Ale dla ustalenia uwagi zobaczmy wygląda cena Passata i jak ta cena zależy od roku produkcji (oś OX), pojemności i typu silnika (poszczególne panele) oraz wersji marki (rożne kolory kropek).

Na panelu po prawej stronie podana jest nazwa marki, liczba ofert sprzedaży aut tej marki, oraz informacja o rozkładzie cen dla aut produkowanych w danym roczniku. Każda niebieska kropka to jedna oferta sprzedaży. Czerwona kropka oznacza medianę cen, czerwone kreski oznaczają odpowiednio kwantyle 5-25% i 75-95%. Ceny są w skali logarytmicznej. Na wykresie mamy tylko dane o autach używanych, w tej kategorii mediany cen układają się prawie idealnie w linii prostej, co odpowiada stałej procentowej zmianie w cenie. Tę informację jeszcze wykorzystamy.

Z wykresów widzimy, że passaty w dieslu z silnikiem 1900cm3 są dosyć popularne. 5 lat to jeszcze wersja B6. Zobaczmy teraz jak wyglądają podobne wykresy dla konkurencji.

Obrazki dla innych marek można znaleźć w katalogu tutaj.

Lub korzystając z poniższej listy marek, dla której zebrano takie wykresy.
Alfa Romeo 147, Alfa Romeo 156, Alfa Romeo 159, Alfa Romeo 166, Alfa Romeo GT,
Audi A2, Audi A3, Audi A4, Audi A5, Audi A6, Audi A6 Allroad, Audi A7, Audi A8, Audi Q5, Audi Q7, Audi TT,
BMW 118, BMW 316, BMW 318, BMW 320, BMW 325, BMW 330, BMW 520, BMW 525, BMW 530, BMW 535, BMW 730, BMW 740, BMW X3, BMW X5, BMW X6,
Chevrolet Aveo, Chevrolet Lacetti,
Chrysler 300C, Chrysler Grand Voyager, Chrysler Pacifica, Chrysler PT Cruiser, Chrysler Sebring, Chrysler Town & Country, Chrysler Voyager,
Citroen Berlingo, Citroen C1, Citroen C2, Citroen C3, Citroen C4, Citroen C4 Picasso, Citroen C5, Citroen C8, Citroen Saxo, Citroen Xsara, Citroen Xsara Picasso,
Dacia Logan,
Daewoo Lanos, Daewoo Matiz, Daewoo Nubira,
Dodge Grand Caravan, Dodge RAM,
Fiat 500, Fiat Brava, Fiat Bravo, Fiat Croma, Fiat Doblo, Fiat Ducato, Fiat Grande Punto, Fiat Marea, Fiat Multipla, Fiat Palio, Fiat Panda, Fiat Punto, Fiat Seicento, Fiat Stilo, Fiat Uno,
Ford C-MAX, Ford Fiesta, Ford Focus, Ford Focus C-Max, Ford Fusion, Ford Galaxy, Ford KA, Ford Mondeo, Ford Mustang, Ford S-Max, Ford Transit,
Honda Accord, Honda Civic, Honda CR-V, Honda Jazz,
Hyundai Accent, Hyundai Coupe, Hyundai Getz, Hyundai i30, Hyundai Matrix, Hyundai Santa Fe, Hyundai Tucson,
Jaguar S-Type, Jaguar X-Type,
Jeep Cherokee, Jeep Grand Cherokee,
Kia Carens, Kia Carnival, Kia Cee’d, Kia Picanto, Kia Rio, Kia Sorento, Kia Sportage,
Lancia Lybra,
Land Rover Discovery, Land Rover Freelander, Land Rover Range Rover,
Lincoln Town Car,
Mazda 2, Mazda 3, Mazda 323, Mazda 323F, Mazda 5, Mazda 6, Mazda 626, Mazda CX-7, Mazda MPV, Mazda Premacy, Mazda RX-8,
Mercedes-Benz A 140, Mercedes-Benz A 160, Mercedes-Benz A 170, Mercedes-Benz A 180, Mercedes-Benz B 180, Mercedes-Benz C 180, Mercedes-Benz C 200, Mercedes-Benz C 220, Mercedes-Benz CLK 200, Mercedes-Benz E 200, Mercedes-Benz E 220, Mercedes-Benz E 270, Mercedes-Benz E 280, Mercedes-Benz E 320, Mercedes-Benz E 350, Mercedes-Benz ML 270, Mercedes-Benz ML 320, Mercedes-Benz ML 350, Mercedes-Benz S 320, Mercedes-Benz S 350, Mercedes-Benz S 500, Mercedes-Benz SLK 200, Mercedes-Benz Sprinter, Mercedes-Benz Vito,
Mini Cooper,
Mitsubishi Carisma, Mitsubishi Colt, Mitsubishi Eclipse, Mitsubishi Galant, Mitsubishi L200, Mitsubishi Lancer, Mitsubishi Outlander, Mitsubishi Pajero, Mitsubishi Space Star,
Nissan Almera, Nissan Almera Tino, Nissan Micra, Nissan Navara, Nissan Note, Nissan Patrol, Nissan Primera, Nissan Qashqai, Nissan Terrano, Nissan X-Trail,
Opel Agila, Opel Astra, Opel Combo, Opel Corsa, Opel Frontera, Opel Insignia, Opel Meriva, Opel Omega, Opel Signum, Opel Tigra, Opel Vectra, Opel Vivaro, Opel Zafira,
Peugeot 106, Peugeot 107, Peugeot 206, Peugeot 206 CC, Peugeot 207, Peugeot 306, Peugeot 307, Peugeot 308, Peugeot 406, Peugeot 407, Peugeot 607, Peugeot 807,
Peugeot Partner,
Porsche 911, Porsche Cayenne,
Renault Clio, Renault Espace, Renault Grand Espace, Renault Grand Scenic, Renault Kangoo, Renault Laguna, Renault Megane, Renault Modus, Renault Scenic, Renault Thalia, Renault Trafic, Renault Twingo,
Renault Vel Satis,
Rover 25, Rover 45, Rover 75,
Saab 9-3, Saab 9-5,
Seat Alhambra, Seat Altea, Seat Arosa, Seat Cordoba, Seat Ibiza, Seat Leon, Seat Toledo,
Skoda Fabia, Skoda Felicia, Skoda Octavia, Skoda Roomster, Skoda Superb,
Smart Fortwo,
Subaru Forester, Subaru Impreza, Subaru Legacy,
Suzuki Grand Vitara, Suzuki Jimny, Suzuki Swift, Suzuki SX4,
Toyota Auris, Toyota Avensis, Toyota Aygo, Toyota Camry, Toyota Celica, Toyota Corolla, Toyota Corolla Verso, Toyota Land Cruiser, Toyota RAV-4, Toyota Yaris,
Volkswagen Bora, Volkswagen Caddy, Volkswagen Caravelle, Volkswagen Fox, Volkswagen Golf, Volkswagen Golf Plus, Volkswagen Jetta, Volkswagen Lupo, Volkswagen Multivan, Volkswagen New Beetle, Volkswagen Passat, Volkswagen Passat CC, Volkswagen Polo, Volkswagen Sharan, Volkswagen Tiguan, Volkswagen Touareg, Volkswagen Touran, Volkswagen Transporter,
Volvo C30, Volvo S40, Volvo S60, Volvo S80, Volvo V40, Volvo V50, Volvo V70, Volvo XC 70, Volvo XC 90

Ranking uczelni i paretooptymalność

Paretooptymalność oznacza, że nie ma innych obiektów lepszych w każdym rozważanym sensie. Poza tym to ładne słowo, jedno z niewielu których jeszcze nie zapomniałem ze wstępu do teorii gier. Front paretooptymalny to zbiór obiektów, takich że od żadnego z nich nie istnieje obiekt jednostajnie lepszy. Obiekty z frontu są nieporównywalne.

Tak tez pewnie jest z uczelniami, jedne są lepsze jeżeli chodzi o współpracę z przemysłem, inne mogą mieć więcej publikacji, jeszcze inne mogą mieć lepszy PR. Na kilku rysunkach będę chciał pokazać w jakich kategoriach które uczelnie są wysoko. W zbiorze danych o rankingu z Rzeczpospolitej mamy 32 atrybuty, to za dużo a techniki skalowania nie doprowadziły do niczego co byłoby łatwe w interpretowaniu. Dlatego poniżej te 32 atrybuty zamieniłem w 6 atrybutów odpowiadających 6 grupom atrybutów z oryginalnego rankingu. Nowy atrybut to suma punktów z atrybutów dla danej grupy, tz. atrybut prestiż to suma punktów za preferencje pracodawców, ocenę kadry akademickiej, uznanie międzynarodowe i wybór olimpijczyków. Reszta podobnie. Na każdym z poniższych wykresów pokazywane jest 90 uczelni w układzie opisanym przez dwa atrybuty. Każdy punkt odpowiada jednej uczelni. 10 uczelni które mają najwyższą sumę atrybutów jest oznaczonych przez nazwy a kropki odpowiadające tym uczelnią maja ciemniejsze kolory. Kolory oznaczają typ uczelni (wyróżniłem 6 typów, w tym uniwersytety, politechniki itp). Jasny pomarańcz to politechnika która nie jest jedną z 10 najlepszych uczelni w danym zestawieniu, ciemno pomarańczowy to politechnika która znalazła się w 10 najlepszych w danym zestawieniu. Obie osie są w skali pierwiastkowej bo tak lepiej wyglądają.

Jeżeli mamy 6 atrybutów to wszystkich par jest 30 (kolejność ma znaczenie bo mózg inaczej interpretuje składowa pozioma a inaczej pionową). Wykresy dla wszystkich par są w katalogu tutaj. Poniżej przedstawiam 6 wybranych najciekawszych wykresów. Skrypt użyty do wygenerowania tych wykresów jest tutaj. Informacja, które charakterystyki agregują poszczególne grupy można odczytać też z tego wykresu.

Efektywność naukowa a innowacyjność.

Żółte kropki na tym wykresie oznaczają uniwersytety medyczne. Są one oznaczone jako uczelnie o bardzo wysokiej efektywności naukowej ale bardzo niskiej innowacyjności. Liter rankingu Rzeczpospolitej (uniwersytet Warszawski) nie jest najlepszy ani pod względem efektywności naukowej ani pod względem innowacyjności. UJ bije go w obu przypadkach. Jeżeli chodzi o innowacyjność to najlepsze są politechniki (Wroclawska, Warszawska, Gdanska, Poznanska) stojące tuz za AGH. Ale do tego jeszcze wrócimy.

Umiędzynarodowienie a potencjał naukowy.

Umiędzynarodowienie nie jest jakoś szczególnie silnie skorelowane z potencjałem naukowym. Czołówka najbardziej umiędzynarodowionych uczelni zawiera i te o wysokim potencjale naukowym (UW, UJ, swoją drogą termin potencjał naukowy brzmi podejrzanie) i te o wynikach gorszych niż średnia w tej kategorii (Akademia Koźmińskiego).

Warunki studiowania a potencjał naukowy

Zależność pomiędzy warunkami studiowania a potencjałem naukowym jest ciekawa. Te dwie cechy wydają się ze sobą korelować, ale jest to wynikiem obecności dwóch skupisk (paradoks Simpsona). Uczelnie o przeciętnych lub niskich warunkach studiowania i przeciętnym lub niskim potencjale naukowym (w tej grupie nie widać korelacji) i grupie kilku uczelnie w którym oba te współczynniki są wysokie.

Btw: W kategorii warunki studiowania wygrywa UAM. Tam jeszcze nie pracowałem, ale proszę śmiało o przesyłanie ofert na gmaila.

Potencjał naukowy a prestiż

Słowo prestiż jest już tak nadużywane, że trudno nawet powiedzieć co to oznacza w kontekście uczelni wyższej. Jest to jeden z najbardziej skośnych charakterystyk, tylko 7 uczelni ma ten współczynnik powyżej 100 a zdecydowana większość ma ten współczynnik poniżej 50.

Ciekawym eksperymentem będzie poprowadzić linię regresji na powyższym rysunku i zobaczyć prestiż których uczelni przewyższa ,,rzeczywisty” potencjał naukowy (PR robi swoje) a dla których uczelni jest on zaniżony.

Innowacyjność a potencjał naukowy

Podobnie wracamy do innowacyjności. Moim zdaniem jest to jedna z ważniejszych wartości wartych promowania. Potrzebujemy jak ryba wody kreatywnych rozwiązań, pomysłów, odważnych działań, zamiast walenia pałkami po łapach lub głowach tych co się wychylają.

Na tym wykresie ciekawie rozkładają się trzy grupy punktów. Niebieskie kropki oznaczają akademie (wojskowe, rolnicze, humanistyczne), za wyjątkiem AGH wszystkie pozostałe akademie mają niską innowacyjność i niski potencjał naukowy. Pomarańczowe kropki to politechniki. Większość z nich cechuje się wysoką innowacyjnością i potencjałem naukowym średnio trochę lepszym niż średnia (Maturzyści, idźcie studiować na Politechniki, potrzebujemy w kraju więcej inżynierów). Fioletowe kropki to uniwersytety, te charakteryzują się wysokim potencjałem naukowym ale innowacyjność jest tylko trochę lepsza niż średnia.