Tag: kod w R

Kiedy dogonimy zachód… od podszewki

Dzisiaj pokażę jak powstawały wyniki przedstawione we wtorkowym wpisie Kiedy dogonimy zachód?. Będzie trochę o motywacji, dlaczego akurat taki temat, o tym jak ściągnąć dane automatycznie z internetu, o zaletach analizy ,,opóźnienia” i przede wszystkim będzie o robieniu wykresów w programie R.

Czytaj dalej Kiedy dogonimy zachód… od podszewki

Ważona średnia, analiza wariancji a wyniki PISA

Zgodnie z zapowiedziami dzisiaj będzie kilka technicznych komentarzy dotyczących wtorkowego wpisu. Pojawi się temat ważenia obserwacji, analizy wariancji, pojawi się bardzo przydatna funkcja by() z pakietu R i oczywiście opis krok po kroku jak zrobić wykres prezentowany w poprzednim odcinku.

Czytaj dalej Ważona średnia, analiza wariancji a wyniki PISA

Jak dojechać z Warszawy do Wrocławia zwiedzając po drodze Jasną Górę?

 

Przygotowując prezentację na Szkołę Infografiki (o której pisałem w poniedziałek)  jako pierwszą ,,infografikę” do pokazania wybrałem wizualizację rozkładu jazdy pociągów pomiędzy Paryżem a Lionem (patrz rysunek po prawej stronie) wykonaną przez E. J. Mareyego.

To dosyć kontrowersyjna grafika, dla jednych zupełnie nieczytelna, innych urzekająca elegancją i pomysłowością. Różne są gusta, w każdym razie mnie się bardzo podoba.

Na tym by się pewnie skończyło, gdyby nie to że z różnych powodów kursuję czasami pomiędzy Warszawą a Wrocławiem. Z różnych środków lokomocji najbardziej lubię jeździć pociągami więc pojawił się pomysł przygotowania podobnego wykresu dla trasy Warszawa-Wrocław. Pomysł zamienił się w chęć, chęć zamieniła się w obsesję której nie powstrzymało nawet dziwne API serwisu pkp.pl.

Niestety w przeciwieństwie do trasy Paryż – Lion trasa Wrocław – Warszawa nie jest trasą bliską najkrótszej możliwej. Trzeba jechać ,,dookoła” i można wybrać albo drogę przez Katowice albo Częstochowę albo Poznań. W wizualizacji którą zaraz pokażę zaznaczone będą połączenia przez wszystkie te trzy ścieżki, ale uniemożliwia to zachowanie takich samych proporcji w odległościach fizycznych pomiędzy miastami co w odległościach na osi OY wykresu. Dokładne mapowanie nie jest możliwe bo np. Częstochowa jest w 52% drogi z Warszawy do Wrocławia jeżeli nie jedzie się przez Katowice albo w 46% drogi jeżeli jedzie się jeszcze przez Katowice. Na trasie Paryż – Lion połączenia odpowiadają jednej możliwej drodze, dlatego też na wykresie Mareya nachylenia krzywych odpowiadają prędkościom, a u mnie nie. Ale postaram się trzymać proporcji możliwie blisko.

Dane zebrałem ze strony pkp.pl, niestety z powodu częstych zmian rozkładu należy te dane traktować jako jedynie ilustrację rozkładu jazdy pociągów na dzień 22 marca 2012 roku. Nie wiadomo jak długo ten rozkład będzie obowiązywał. Dostęp do danych wraz z wybranymi stacjami pośrednimi znajduje się w tym pliku csv.

Do narysowania wizualizacji rozkładu jazdy wykorzystałem oczywiście program R. Kod w R użyty do narysowania tego wykresu znajduje się tutaj.

Poniższa wizualizacja nie jest być może zbyt efektywna, ale można dzięki niej łatwo zaplanować sobie przejażdżkę do Wrocławia z kilkugodzinnym pobytem w Częstochowie, który to można wykorzystać na miłe spędzenie czasu w tymże mieście. Podobne wizualizacje powstają też np. do przedstawiania połączeń autobusowych pomiędzy różnymi punktami w mieście lub połączeń lotniczych pomiędzy miastami.

Różne kolory przedstawiają różne pociągi. Kropki odpowiadają stacjom. Pogrubione linie oznaczają pociągi jadące bezpośrednio z Warszawy do Wrocławia. Po kliknięciu powinna otworzyć się wersja svg tego rysunku, a jeżeli się nie otworzy to można znaleźć wersję png tutaj.


Ceny używanych aut po raz trzeci

 

W październiku poprzedniego roku i 1 stycznia tego roku zebrałem na bazie serwisu www.otomoto.pl zbiór opisujący cen i wiele dodatkowych parametrów z ofert sprzedaży samochodów używanych (zobacz między innymi te wpisy). Dzisiaj dodaję zbiór danych o cenach samochodów z połowy marca roku 2012. Na podstawie tych trzech pomiarów będziemy mogli się przyjrzeć temu jak się zmieniały ceny używanych aut w przeciągu ostatnich 6 miesięcy.

Do pierwszego zbioru danych dołączyłem kilka analiz statystycznych i wizualizacji. Do drugiego zbioru danych dodałem szczegółowy opis jak takie dane pobierać. Zabrakło mi jednak czasu by przedstawić wyniki analiz statystycznych danych zebranych w styczniu.

Okazało się jednak, że nawet surowe dane mogą być przydatne. W tak zwanym między czasie otrzymałem od trzech osób informację, że ten zbiór danych bardzo im się przydał w projekcie licencjackim. A od jednej osoby otrzymałem nawet napisaną już pracę licencjacką poświęconą różnym technikom analizy danych z przykładami na zbiorze danych o samochodach pobranym z tego bloga.
Super!!!
Więcej takich projektów!
Jeżeli autorzy się zgodzą to podlinkuję ich wyniki do odpowiednich wpisów.

Jeden z częściej pojawiających się komentarzy od osób pracujących na zebranych danych dotyczył braku informacji o typie nadwozia. W poprzednio zbieranych danych typ nadwozia nie był bezpośrednio identyfikowalny.

Mając to na uwadze tym razem listę zmiennych opisujących ofertę sprzedaży auta rozszerzyłem o dwie nowe zmienne. Pierwszą jest typ nadwozia a drugą jest kod pocztowy i nazwa miejscowości z której pochodzi dane ogłoszenie.

Ponieważ ostatnio interesowałem się analizą danych przestrzennych i funkcjami do przedstawiania danych na mapach, więc postaram się w najbliższych dniach pokazać co ciekawego można zrobić z cenami samochodów używanych mając też informacje o miejscu złożenia oferty sprzedaży.

Skrypt wczytujący do programu R dane zebrane w marcu 2012 znajduje się tutaj.

Dane w postaci tekstowej i binarnej znajdują się w tym katalogu.

Piątek, chmura słów, TextMining, morfologik i oczywiście R

Dzisiaj jest piątek, więc zamiast zaprzątać sobie głowę liczbami pooglądamy obrazki.

Dwa dni temu, w tym wpisie opisałem jak ściągnąć dane z Twittera i jako przykład ściągnąłem dane z kanału #debataACTA. Pokazałem też kilka podsumowań odkładając na później bardziej złożone analizy. Analizę nastawienia odłożę jeszcze na później, a dzisiaj pokażę jak danych tekstowych zrobić chmurę słów (ang. word cloud) używając R.

Samą chmurę słów można wykonać używając funkcji wordcloud() z pakietu wordcloud. Jako argumenty należy podać listę słów oraz współczynnik skalujący, odpowiadający wielkości danego słowa (najczęściej odpowiadający liczbie wystąpień danego słowa w tekście). Listę słów oraz częstości ich wystąpień można wygenerować używając funkcji str_split() i table(), ale aby było ciekawiej użyję w tym celu pakietu tm, który ma wiele przydatnych narzędzi do analizy tekstu (najwięcej przydatnych narzędzi ma do analizy tekstu angielskiego, ale z polskim też coś można zrobić).

Poniższy fragment kodu wczytuje dane, liczy tablicę częstości a następnie przedstawia ją graficznie z użyciem funkcji wordcloud(),

 

# Odczytaj dane
tweets = read.table(file="debataACTA4_13_luty.csv", sep=";", head=T)
 
# funkcje z pakietu tm i wordcloud
# do tworzenia korpusu wyrazów i ich wizualizacji
actaSlownik = Corpus(DataframeSource(data.frame(tweets[,2])))
actaSlownik = tm_map(actaSlownik, removePunctuation)
actaSlownik = tm_map(actaSlownik, tolower)
tdm = TermDocumentMatrix(actaSlownik)
m = as.matrix(tdm)
v = sort(rowSums(m),decreasing=TRUE)
 
# mając policzone częstości występowania możemy je zwizualizować
wordcloud(names(v), v^0.3, scale=c(5,0.5),random.order=F, colors="black")

 

Wykres otrzymany w ten sposób nie wygląda zbyt dobrze, ponieważ słowa występują w najróżniejszych odmianach a liczenie częstości osobno dla ,premierze’, 'premier’, 'premiera’, 'premierem’ nie wiele daje. Aby wykres wyglądał lepiej chcielibyśmy dla każdego ze słów znaleźć jego rdzeń i zliczać liczbę wystąpień rdzeni a nie różnych form. Jak to zrobić? Potrzebny będzie lematyzator, a dokładniej tzw. stemmer, czyli narzędzie, które dla każdego słowa wyodrębni tzw. 'stem’ (nie znam niestety polskiej terminologii a nie chcę tworzy kwiecistych tłumaczeń) czyli część słowa nie ulegającą odmianie.

Wykorzystam bezpłatny morfologik-stemmer dostępny na blogu http://morfologik.blogspot.com/ rozwijany przez Dawida Weissa i współpracowników. Narzędzie z którego skorzystałem to duży plik tekstowy mający w jednej kolumnie różne formy słów, w drugiej rdzenie a w trzeciej informacje o formie gramatycznej danego słowa. Oczywiście nie wszystkie słowa z twittera można znaleźć w tym zbiorze, nie ma tam np. wszystkich nazwisk. W każdym razie pierwsza transformacja polegała na przemapowaniu słowa na jego rdzeń, jeżeli słowo występuje w słowniku, lub pozostawienie słowa bez zmiany jeżeli w słowniku nie występuje. Dzięki temu z 13334 różnych słów zostaliśmy z 8597 słowami, z czego zdecydowana większość występuje tyko raz i nie znajdzie się na mapie tagów.

Zobaczmy jak wygląda rzeczona mapa.

[Rysunek 1. Chmura słów występujących przynajmniej dwa razy w zapisach z kanału debataACTA. Gdy było to możliwe słowa zostały przekształcone do swoich rdzeni. Wersja wektorowa tego rysunku (uwaga 7MB) znajduje się tutaj]

Dużo tych słów, mało widać, zróbmy więc jeszcze jedną iterację. Usuńmy wszystko co nie jest rzeczownikiem. Szczęśliwie słownik z pakietu morfologik ma informacje o tym czy dane słowo jest czy nie rzeczownikiem, więc zostawiamy tylko słowa, które znajdują się w słowniku i mają 'subst’ w trzeciej kolumnie.

[Rysunek 2. Chmura rzeczowników występujących przynajmniej dwa razy w zapisach z kanału debataACTA. Wersja wektorowa tego rysunku (uwaga 8MB) znajduje się tutaj]

Zróbmy jeszcze jeden eksperyment, mianowicie sprawdźmy jakie inne kanały występowały w wiadomościach z kanału #debataACTA. Zostawiamy więc tylko słowa zaczynające się od znaku #.

 

[Rysunek 3. Chmura nazw kanałów w wiadomościach na kanale debataACTA. Wersja wektorowa tego rysunku (uwaga 7MB) znajduje się tutaj]

#debataACTA, twitteR, Twitter i R

Zaczęło się niewinnie. Dwa tygodnie temu kolega Grzesiek P. powiedział, że analizę nastawienia (ang. Sentiment analysis) na dużych ilościach tekstu robi się banalnie. Tydzień temu kolega Paweł Ch. powiedział, że API twittera ma limit do 70 zapytań na minutę. W piątek przy okazji rozmowy o raporcie ,,Obiegi Kultury” Alek T. zapytał mnie czy mam doświadczenie w analizie danych z Twittera bo chętnie zobaczyłby się działo na kanale (moje autorskie tłumaczenie hashtaga) #debataACTA. Ponieważ nie interesowałem wcześniej się ani analizą nastawienia, ani API twittera ani hashtagami, stwierdziłem, że warto zobaczyć co w trawie ćwierka.

Tak się składa, że dzięki Jeffowi Gentry’emu API Twittera jest łatwo dostępne z poziomu R. Wystarczy załadować pakiet twitteR i kilka prostych funkcji pozwala na łatwą interakcję z ćwierkami z serwisu Twitter.com (ok, dalej będę używał tłumaczenia ,,wiadomościami”).

Zobaczmy prosty przykład.

> # interfejs do API twittera
> library(twitteR)
> # pobierz maksymalnie 1500 wiadomości zawierających etykietę #debataACTA
> # stworzonych 7 lutego 2012
> tweets = searchTwitter('#debataACTA', n=1500, until="2012-02-08", since="2012-02-06")
> # zamień wyniki na ramkę danych i wyświetl treść oraz autora pierwszej z wiadomości
> debataACTA =  twListToDF(tweets)
> debataACTA[1,c(1,4,10)]
                                                                                    text             created screenName
1 Mamy tyle Kultury i Edukacji ile na nia wydajemy 0,75% Budzetu #debataACTA #PremierRP 2012-02-07 15:37:46     TPHMAC

Ten pakiet ma te same ograniczenia co API, czyli 70 zapytań na minutę. Ale każde zapytanie to 25 wiadomości, a więc w sumie możemy ściągnąć do 1500 wiadomości na minutę.

Prostą pętlą, z pewną pomocą internetu, zebrałem wpisy z Twittera z okresu od 4 do 12 lutego 2012 zawierające etykietę #debataACTA. po oczyszczeniu zapisałem je w postaci pliku csv, który można pobrać z adresu tutaj.

Mamy więc dane, przyjrzyjmy się im. Do tematu analizy nastawienia podejdziemy w piątek, a dziś przyjrzyjmy się jak wyglądała aktywność ćwierkających. Na rysunku 1 można zobaczyć ile wpisów pojawiało się z tą etykietą. Temat jak widać skończył się z samą debatą. Na rysunku 2 można zobaczyć bliżej interesujący okres liczby wiadomości nadczas trwania debaty. Rzuca się w oczy szczególnie ta chwila ciszy przed burzą.

[Rysunek 1. Liczba wiadomości zawierających etykietę #debataACTA na godzinę. Kliknij by powiększyć.]

[Rysunek 2. Liczba wiadomości zawierających etykietę #debataACTA na 5 minut (debata rozpoczęła się po godzinie 14). Kliknij by powiększyć.]

Gdyby ktoś chciał sprawdzić o czym mówiono w której minucie to komentarze z znaczkami synchronizacyjnymi znaleźć można np. tutaj.

W sumie mamy ponad 6.5 tysiąca wiadomości. Zobaczmy kto je pisał. Zamiast wpisywać ksywy ćwierkających pokażę ilu z nich generuje odpowiedni procent wszystkich wiadomości.

[Rysunek 3. Na osi y przedstawiono liczbę wiadomości wygenerowanych przez X najaktywniejszych użytkowników. Okazuje się, że 7% piszących (dokładnie 72) generuje 50% wiadomości. A tylko 2% (dokładnie 21) generuje 25% wszystkich wiadomości. Kliknij by powiększyć.]

Ściągnęliśmy dane, zobaczyliśmy podstawowe podsumowania, czas zabrać się za analizę nastawienia. Więcej na ten temat w piątek (muszę ją jeszcze zrobić) a dzisiaj zakończę podsumowaniem z procentowym udziałem wiadomości zawierających ikonę buźki (;-), ;), :-), :))) itp).

[Rysunek 4. Procentowy udział wiadomości zawierających buźkę. Kliknij by powiększyć.]

Jak widać przez pierwsze dwie godziny debaty nie było ćwierkającym do śmiechu.

Informacje o wynagrodzeniach, jak je pokazywać i jak tego nie robić

Informacje o wynagrodzeniach z pewnością są interesujące dla bardzo wielu osób. Tym bardziej zaskakujące jest jak trudno dostać wiarygodne informacje na ten temat. A nawet jak już się uda komuś je zebrać to zaskakujące jest jak nieudolnie te cenne dane są prezentowane.

Jakiś czas temu na portalu gazetapraca znalazłem artykuł o zarobkach w dziale logistyki i transportu (link do artykułu).

Pierwszy wykres z artykułu prezentuje medianę i kwartyle i daje nadzieję, na to że czegoś będzie można się o wynagrodzeniach w tym sektorze dowiedzieć.

Niestety drugi wykres nadaje się jedynie jako przykład na zajęcia z wizualizacji pokazujący jak przedstawiać danych.

[Rys 1. Wykres z przytaczanego artykułu porównujący zarobki na różnych stanowiskach w firmach z kapitałem polskim vs. zagranicznym]

 

Pomijam sam fakt, że trudno odgadnąć czy wysokość słupka odpowiada medianie, czy średniej (btw: odpowiada medianie), bo to jest przesłonięte przez pytanie po co wykres został pokazany pod jakimś dziwnym kątem. Dlaczego zarobki 3100 specjalisty wglądają na mniejsze niż 2888 pracownika szeregowego? Czemu służy oś OY skoro w żaden sposób nie da się jej zastosować do odczytywania wysokości słupków? Jak zmieniłby się nasz odbiór tego wykresu gdyby zamienić kolejnością słupki niebieskie i czerwone?

Zobaczmy ten sam wykres w standardowych osiach współrzędnych.

[Rys21. Wykres słupkowy, dane jak powyżej, ale bez trzeciego wymiaru i obrotów.]

Pierwszą rzeczą, która rzuca się teraz w oczy to jak bardzo zarobki zarządu odstają od zarobków na innych stanowiskach. Drugą jest to, że można teraz odczytać proporcje nie patrząc na liczby a jedynie z wysokości słupków i lokalizacji pomocniczych linii.

Dla osób nieprzekonanych, proponuję mały eksperyment. Usuńmy liczby opisujące wysokości słupków i zobaczmy czy uda się odtworzyć zarobki na stanowisku kierowniczym lub starszego specjalisty (czy uda się zejść z błędem predykcji poniżej 20%). To ciekawy eksperyment do przeprowadzenia np. na wspomnianym  wykładzie z wizualizacji.

Poniżej wykres do eksperymentowania (proszę odgadnąć poziom słupków dla stanowiska kierowniczego lub starszego specjalisty).

zanuRkuj w R

Zostałem zaproszony do wygłoszenia referatu na seminarium poświęconym analizie danych. A ponieważ zapraszającym był sam Paweł Cichosz więc nie sposób było odmówić.

Seminarium odbędzie sie dziś o 12:15 (Gmach Elektroniki PW, sala 229 II piętro). Przy okazji powstała prezentacja, która może zaciekawi szersze grono odbiorców.

 

Slajdy są dostępne tutaj, kod w programie R jest dostępny tutaj.

Prezentacja miała nie być o statystyce, ale o języku R. Na omówienie języka jedno seminarium nie wystarczy, dlatego wybrałem siedem ciekawszych i bardziej zaawansowanych mechanizmów języka R. Dotyczą one takich aspektów jak waga atrybutów, funkcyjność, leniwa ewaluacja, przestrzenie nazw, profilowanie kod, automatyczne generowanie raportów z użyciem Sweave. Napisałem ,,zaawansowanych’; ponieważ nie są one omawiane w większości opracowań dotyczących R, a są moim zdaniem ciekawe i warto być ich świadomym.